function  y  = exam_3_26( zeta )
%EXAM_3_26 二阶系统的阶跃时间响应
%   zeta阻尼系数
%   y阶跃响应
%   徐国保于2017年2月26日编写
t=0:0.1:30;
if(zeta>=0)&(zeta<1)
    y=y1(zeta,t);
elseif zeta==1
    y=y2(zeta,t);
else
    y=y3(zeta,t);
end
plot(t,y)
title(['The second order response of zeta=',num2str(zeta)])
xlabel('time(t)')
ylabel('response(y)')
    function y=y1(zeta,t)
        %阻尼系数0<zeta<1的二阶系统阶跃响应
        y=1-1/sqrt(1-zeta^2)*exp(-zeta*t).*sin(sqrt(1-zeta^2)*t+atan(sqrt(1-zeta^2)/zeta));
    end
    function y=y2(zeta,t)
        %阻尼系数zeta=1的二阶系统阶跃响应
        y=1-exp(-zeta*t).*(1+t);
    end
    function y=y3(zeta,t)
        %阻尼系数zeta>1的二阶系统阶跃响应
        sq=sqrt(zeta^2-1);
        y=1-1/(2*(1-zeta^2+zeta*sq))*exp(-(zeta-sq)*t)-1/(2*(1-zeta^2-zeta*sq))*exp(-(zeta+sq)*t);
    end
end

