%   t=[0:0.01:20];
  if stepabcd==1
   windstep=figure('color',[0.1,0.7,0.8],'position',[0.1*winw,0.1*winh,0.8*winw,0.8*winh],...
       'name','離散系統W域輸出Y(t)和狀態X(t)之階躍響應','numbertitle','off','menubar','figure',...
       'resize','off');
 end
  if stepabcd==0;
    figure(windstep)
 end

 stepabcd=0;
 
%以上完成系統的階躍響應曲線作圖
rr1=size(bnum);
rr2=size(bden);
finderr=0;
%判斷是否為比例環節
if ((rr1(1,2)==1)&(rr2(1,2)==1)),
   finderr=1;
   dcg=dcgain(bnum,bden);
   [yc x]=step(bnum,bden);
   nsn=length(yc);t=1:nsn;
%   [yc x]=step(bnum,bden,t);
   yc1=yc/dcg;
   subplot(1,1,1);
   plot(t,yc1,'r',t,yc,'b');
   grid on;hold on
   xlabel('s/t');
   ylabel('Y(t)/dcg');
   legend('紅線:Y(t)/dcg階躍響應曲線','藍線:Y(t)階躍響應曲線',0);
   strtext=['這是比例環節，系統指標: 穩態值=' num2str(dcg) ' ; tr,tp,Mp,ts 無意義 '];
   title(strtext);
else,
%判斷是否為積分環節
findjf1=0;
findjf2=0;
for i=1:1:rr1(1,2)-1
   if bnum(1,i)==0,
   else,findjf1=1;break;
   end
end
if bden(1,rr2(1,2))==0,
   for i=1:1:rr2(1,2)-2
     if bden(1,i)==0,
     else,findjf2=1;break;
     end
   end
else,
   findjf2=1;
end
if (findjf1==0&findjf2==0),
   finderr=1;
   dcg=dcgain(a,b,c,d);
   [yc x]=step(a,b,c,d,1);
 %  [yc x]=step(a,b,c,d,1,t);
   nsn=length(yc);t=1:nsn;
   
   yc1=yc/dcg;
   subplot(2,1,1);
   plot(t,x,'r');
   grid on;hold on
   xlabel('s/t');
   ylabel('X(t)');
   legend('紅線:X(t)狀態曲線',0);
   title('X(t)狀態曲線');
   subplot(2,1,2);
   grid on;hold on
   step(a,b,c,d,1);
 %  plot(t,yc1,'r',t,yc,'b');
   grid on;hold on
   xlabel('s/t');
   ylabel('Y(t)和Y(t)/dcg');
   legend('紅線:Y(t)/dcg階躍響應曲線','藍線:Y(t)階躍響應曲線',0);
   strtext=[' 這是積分環節，無穩態值且tr,tp,Mp,ts 無意義 '];
   title(strtext);
end
end
if finderr==0,
dcg=dcgain(a,b,c,d);
%[yc x]=step(a,b,c,d,1,t);

[yc x]=step(a,b,c,d,1);
nsn=length(yc);t=1:nsn;

yc1=yc/dcg;
subplot(2,1,1);
plot(t,x,'r');
grid on;hold on
xlabel('s/t');
ylabel('X(t)');
legend('紅線:X(t)狀態曲線',0);
title('X(t)狀態曲線');
subplot(2,1,2);
%plot(t,yc1,'r',t,yc,'b');
grid on;hold on
step(a,b,c,d,1)
%plot(t,yc,'b')
grid on;hold on
ylabel('Y(t)和Y(t)/dcg');
%legend('紅線:Y(t)/dcg階躍響應曲線','藍線:Y(t)階躍響應曲線',0);
%********************
[z p]=tf2zp(bnum,bden);
ii=find(real(p)>0);n1=length(ii);
if (n1<=0),%系統穩定
  syms s t;
  [numlength1 numlength2]=size(bnum);
  [denlength1 denlength2]=size(bden);
  l1=['('];
  for i=1:numlength2-1
     l1=[l1 num2str(bnum(1,i)) '*s^' num2str(numlength2-i) '+'];
  end;
  l1=[l1 num2str(bnum(1,numlength2)) ')'];
  l2=['('];
  for i=1:denlength2-1
     l2=[l2 num2str(bden(1,i)) '*s^' num2str(denlength2-i) '+'];
  end;
  l2=[l2 num2str(bden(1,denlength2)) ')'];
  l3=[l1 '/(s*' l2 ')'];
  fun=ilaplace(eval(l3));
  %以上進行逆拉氏變換得出系統在階躍輸入下的時間響應函數
  ycn=step(bnum,bden);
  tn=length(ycn);
  ty=0:0.01:tn;cin=length(ty);
  ycnn=step(bnum,bden,ty);
  [z,p]=tf2zp(bnum,bden);
  
  
uicontrol(windstep,'style','push','string','穩態值',...
       'position',[0.06,0.01,0.11,0.05],'units','normalized',...
      'callback','mysteppre_dcg');
uicontrol(windstep,'style','push','string','上升時間tr',...
       'position',[0.22,0.01,0.11,0.05],'units','normalized',...
      'callback','mysteppre_tr');
uicontrol(windstep,'style','push','string','峰值時間tp和最大超調量Mp',...
       'position',[0.37,0.01,0.25,0.05],'units','normalized',...
      'callback','mysteppre_tpmp');
uicontrol(windstep,'style','push','string','調整時間ts',...
       'position',[0.66,0.01,0.11,0.05],'units','normalized',...
      'callback','mysteppre_ts');
uicontrol(windstep,'style','push','string','結束',...
       'position',[0.83,0.01,0.11,0.05],'units','normalized',...
      'callback','closestep');
else,
  strtext=['系統不穩定，系統指標dcg,tr,tp,mp,ts無法計算！'];title(strtext);
end
end
%else
%   abcdefg=0;
%   ksys=tf(bnum,bden);
%   ltiview('step',ksys)
%end