控制系統(tǒng)仿真實(shí)驗(yàn)報(bào)告11456

1、控制系統(tǒng)仿真實(shí)驗(yàn)報(bào)告班級(jí)：測(cè)控1402班 姓名：王瑋 學(xué)號(hào)：1405040207 2018年01月實(shí)驗(yàn)一 經(jīng)典的連續(xù)系統(tǒng)仿真建模方法一 實(shí)驗(yàn)?zāi)康?1 了解和掌握利用仿真技術(shù)對(duì)控制系統(tǒng)進(jìn)行分析的原理和步驟。2 掌握機(jī)理分析建模方法。3 深入理解階常微分方程組數(shù)值積分解法的原理和程序結(jié)構(gòu)，學(xué)習(xí)用Matlab編寫(xiě) 數(shù)值積分法仿真程序。4 掌握和理解四階Runge-Kutta法，加深理解仿真步長(zhǎng)與算法穩(wěn)定性的關(guān)系。二 實(shí)驗(yàn)內(nèi)容:1. 編寫(xiě)四階 Runge_Kutta 公式的計(jì)算程序，對(duì)非線(xiàn)性模型（3）式進(jìn)行仿真。（1） 將閥位u 增大10和減小10，觀察響應(yīng)曲線(xiàn)的形狀;（2） 研究仿真步長(zhǎng)對(duì)穩(wěn)定性的影

2、響，仿真步長(zhǎng)取多大時(shí)RK4 算法變得不穩(wěn)定？（3） 利用 MATLAB 中的ode45()函數(shù)進(jìn)行求解，比較與（1）中的仿真結(jié)果有何區(qū)別。2. 編寫(xiě)四階 Runge_Kutta 公式的計(jì)算程序，對(duì)線(xiàn)性狀態(tài)方程（18）式進(jìn)行仿真（1） 將閥位增大10和減小10，觀察響應(yīng)曲線(xiàn)的形狀;（2） 研究仿真步長(zhǎng)對(duì)穩(wěn)定性的影響，仿真步長(zhǎng)取多大時(shí)RK4 算法變得不穩(wěn)定？（4） 閥位增大10和減小10，利用MATLAB 中的ode45()函數(shù)進(jìn)行求解階躍響應(yīng)，比較與（1）中的仿真結(jié)果有何區(qū)別。三 程序代碼:龍格庫(kù)塔:%RK4文件 clccloseH=1.2,1.4'u=0.55; h=1;TT=;XX=

3、;for i=1:h:200k1=f(H,u);k2=f(H+h*k1/2,u);k3=f(H+h*k2/2,u);k4=f(H+h*k3,u);H=H+h*(k1+2*k2+2*k3+k4)/6;TT=TT i;XX=XX H;end;hold onplot(TT,XX(1,:),'-',TT,XX(2,:);xlabel('time')ylabel('H')gtext('H1')gtext('H2')hold on水箱模型:function dH=f(H,u)k=0.2;u=0.5;Qd=0.15;A=2;a1

4、=0.20412;a2=0.21129;dH=zeros(2,1);dH(1)=1/A*(k*u+Qd-a1*sqrt(H(1);dH(2)=1/A*(a1*sqrt(H(1)-a2*sqrt(H(2);2編寫(xiě)四階 Runge_Kutta 公式的計(jì)算程序，對(duì)線(xiàn)性狀態(tài)方程（18）式進(jìn)行仿真：1 閥值u對(duì)仿真結(jié)果的影響U=0.45;h=1; U=0.5;h=1; U=0.55;h=1;2 步長(zhǎng)h對(duì)仿真結(jié)果的影響:U=0.5;h=5; U=0.5;h=20; U=0.5;h=39 U=0.5;h=50 由以上結(jié)果知,仿真步長(zhǎng)越大,仿真結(jié)果越不穩(wěn)定。采用ode45算法程序如下:function dH=

5、liu(t,H)k=0.2;u=0.45;Qd=0.15;A=2;a1=0.20412;a2=0.21129;dH=zeros(2,1);dH(1)=1/A*(k*u+Qd-a1*sqrt(H(1);dH(2)=1/A*(a1*sqrt(H(1)-a2*sqrt(H(2);在命令窗口運(yùn)行以下程序:t,H=ode45('liu',1 200,1.2 1.1);plot(t,H(:,1),'r','+',t,H(:,2),'g','*')u=0.45 u=0.5 u=0.55用ode45與用龍格庫(kù)塔法仿真結(jié)果基本一致。

6、2編寫(xiě)四階 Runge_Kutta 公式的計(jì)算程序，對(duì)線(xiàn)性狀態(tài)方程（18）式進(jìn)行仿真：%RK4文件 clcclearclosex=1.2,1.4'u=0.5; h=5;TT=;XX=;for i=1:h:200k1=f2(x,u);k2=f2(x+h*k1/2,u);k3=f2(x+h*k2/2,u);k4=f2(x+h*k3,u);x=x+h*(k1+2*k2+2*k3+k4)/6;TT=TT i;XX=XX x;end;hold onplot(TT,XX(1,:),'-',TT,XX(2,:);xlabel('time')ylabel('x&

7、#39;)gtext('x1')gtext('x2')hold on線(xiàn)性函數(shù)：function dx=f2(x,u)%線(xiàn)性Qd=0.1;a1=0.20412;a2=0.21129;A=2;k=0.2;dx=zeros(2,1);dx(1)=1/A*(k*u-x(1)/(2*sqrt(1.5)/a1)+Qd);dx(2)=1/A*(x(1)/(2*sqrt(1.5)/a1)-x(2)/(2*sqrt(1.4)/a2);1 閥值u對(duì)仿真結(jié)果的影響:U=0.45；h=1； U=0.5；h=1； U=0.55；h=1；2 步長(zhǎng)h對(duì)仿真結(jié)果的影響:U=0.5；h=5； U

8、=0.5；h=20； U=0.5；h=35； U=0.5；h=50； 當(dāng)步長(zhǎng)為50時(shí)仿真結(jié)果開(kāi)始不穩(wěn)定，仿真步長(zhǎng)越大,仿真結(jié)果越不穩(wěn)定。采用ode45算法程序如下:function dx=liu2(t,x)%ÏßÐÔQd=0.1;u=0.45;a1=0.20412;a2=0.21129;A=2;k=0.2;dx=zeros(2,1);dx(1)=1/A*(k*u-x(1)/(2*sqrt(1.5)/a1)+Qd);dx(2)=1/A*(x(1)/(2*sqrt(1.5)/a1)-x(2)/(2*sqrt(1.4)/a2);在命令窗口運(yùn)行以下程序:>

9、> t,x=ode45('liu2',1 200,1.2 1.4);>> plot(t,x(:,1),'r','+',t,x(:,2),'g','*')U=0.45； u=0.5； U=0.55用ode45與用龍格庫(kù)塔法仿真結(jié)果基本一致。四 思考題：1 討論仿真步長(zhǎng)對(duì)穩(wěn)定性和仿真精度的影響仿真步長(zhǎng)越短，穩(wěn)定性越高，精確度越高。2 你是怎樣實(shí)現(xiàn)閥位增大和減小的，對(duì)于非線(xiàn)性模型和線(xiàn)性模型方法一樣嗎？通過(guò)改變u來(lái)改變閥的開(kāi)度，線(xiàn)性系統(tǒng)和非線(xiàn)性系統(tǒng)方法一樣。實(shí)驗(yàn)二 面向結(jié)構(gòu)圖的仿真第一部分 線(xiàn)性系統(tǒng)仿真一

10、 實(shí)驗(yàn)?zāi)康? 掌握理解控制系統(tǒng)閉環(huán)仿真技術(shù)。2 掌握理解面向結(jié)構(gòu)圖的離散相似法的原理和程序結(jié)構(gòu)。3 掌握 MATLAB 中C2D 函數(shù)的用法，掌握雙線(xiàn)性變換的原理。二 實(shí)驗(yàn)內(nèi)容根據(jù)上面的各式，編寫(xiě)仿真程序，實(shí)現(xiàn)無(wú)擾動(dòng)時(shí)給定值階躍仿真實(shí)驗(yàn)1. 取K P = 1.78 ，T i = 85 s T = 10s，H2 S =H2set_ percent = 80, Q d = 0，tend = 700，進(jìn)行仿真實(shí)驗(yàn)，繪制響應(yīng)曲線(xiàn)。clcclear all A=2;ku=0.1/0.5;H10=1.5;H20=1.4;alpha12 = 0.25/sqrt(H10);alpha2 = 0.25/sqrt

11、(H20);R12=2*sqrt(H10)/alpha12;R2=2*sqrt(H20)/alpha2;H1SpanLo=0;H2SpanLo=0;H1SpanHi=2.52;H2SpanHi=2.52;Kp=1.78;Ti=85; R12*AR12 ad = 1/(A*R12);a1 = 1/(A*R12);a2 = 1/(A*R2);Kc=Kp/Ti;bc=Ti;Kd = 1/A;K1 = ku/A;K2 = 1/(A*R12); uc(1)=0;ud(1)=0;u1(1)=0;u2(1)=0;xc(1)=0;xd(1)=0;x1(1)=0;x2(1)=0;yd(1)=0;yc(1)=0;

12、y1(1)=0;y2(1)=0; nCounter = 70;T=10;k=1;deltaQd=0; H20_percent=(H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100;H2=80;tend = nCounter*T;for t=T:T:tend k=k+1; uc(k)= (H2 - (y2(k-1)+H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100)/100; ud(k)=deltaQd; u1(k)=yc(k-1); u2(k)=y1(k-1); xc(k) = xc(k-1) + Kc*T*uc(k-1); yc(k)=xc(

13、k)+bc*Kc*uc(k); xd(k) = exp(-ad*T)*xd(k-1) + Kd/ad*(1-exp(-ad*T)*ud(k);yd(k)=xd(k); x1(k) = exp(-a1*T)*x1(k-1) + K1/a1*(1-exp(-a1*T)*u1(k);y1(k)=x1(k); x2(k) = exp(-a2*T)*x2(k-1) + K2/a2*(1-exp(-a2*T)*u2(k);y2(k)=x2(k); endHlevel(:,1)=(y1+H10-H1SpanLo)/(H1SpanHi-H1SpanLo)*100;Hlevel(:,2)=(y2+H20-H2S

14、panLo)/(H2SpanHi-H2SpanLo)*100;yc=(yc+0.5)*100;y2sp=H2*ones(size(y1');yv=yc; textPositionH1=max(Hlevel(:,1); textPositionH2=max(Hlevel(:,2); H2Steady=Hlevel(size(Hlevel(:,1),1),1)*ones(size(y1'); xmax=max(0:T:tend); xmin=0; ymax=110; ymin=50; scrsz = get(0,'ScreenSize'); gca=figure(&

15、#39;Position',5 10 scrsz(3)-10 scrsz(4)-90); %gca=figure('Position',5 10 scrsz(3)/2 scrsz(4)/1.5) set(gca,'Color','w'); plot(0:T:tend,Hlevel(:,1),'r','LineWidth',2) hold on plot(0:T:tend,Hlevel(:,2),'b','LineWidth',2) hold on plot(0:T:tend,

16、yv,'k','LineWidth',2) hold on plot(0:T:tend,y2sp,'g','LineWidth',2) hold on plot(0:T:tend,H2Steady,'y','LineWidth',2) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/10 (ymax-ymin)/2+ymin-(ymax-ymin)/10,'Color','r','LineWidth&#

17、39;,6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/10,' 第一個(gè)水箱的液位H1','FontSize',16) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/6 (ymax-ymin)/2+ymin-(ymax-ymin)/6,'Color','b','LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/6,' 第二個(gè)水箱的

18、液位H2','FontSize',16) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/4.2 (ymax-ymin)/2+ymin-(ymax-ymin)/4.2,'Color','g','LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/4.2,' 第二個(gè)水箱的液位給定值','FontSize',16) line(tend/2 tend/2+27,(ymax-

19、ymin)/2+ymin-(ymax-ymin)/3.2 (ymax-ymin)/2+ymin-(ymax-ymin)/3.2,'Color','k','LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/3.2,'閥位變化情況','FontSize',16) axis(xmin xmax ymin ymax); text(tend/5,ymax+1.5,' 實(shí)驗(yàn)二 不考慮閥位飽和特性時(shí)的控制效果','FontSize

20、9;,22) grid2. 用 MATLAB 求出從輸入到輸出的傳遞函數(shù)，并將其用c2d 函數(shù)，利用雙線(xiàn)性變換法轉(zhuǎn)換為離散模型，再用dstep()函數(shù)求離散模型的階躍響應(yīng)，階躍幅值為3。clcclear allA=2;ku=0.1/0.5;H10=1.5;H20=1.4;alpha12 = 0.25/sqrt(H10);alpha2 = 0.25/sqrt(H20);R12=2*sqrt(H10)/alpha12;R2=2*sqrt(H20)/alpha2;H1SpanLo=0;H2SpanLo=0;H1SpanHi=2.52;H2SpanHi=2.52;Kp=1.78;Ti=85; R12*

21、AR12 ad = 1/(A*R12);a1 = 1/(A*R12);a2 = 1/(A*R2);Kc=Kp/Ti;bc=Ti;Kd = 1/A;K1 = ku/A;K2 = 1/(A*R12);numc=Kc*bc,Kc;% 用 MATLAB 求出從輸入到輸出的傳遞函數(shù)，denc=1;num1=K1;den1=1,a1;num2=K2;den2=1,a2;gc=tf(numc,denc);g1=tf(num1,den1);g2=tf(num2,den2);Sysq=gc*g1*g2;SysG=feedback(Sysq,1);gg=c2d(SysG,10,tustin);% 用c2d 函數(shù)，

22、利用雙線(xiàn)性變換法轉(zhuǎn)換為離散模型dstep(3*gg.num1,gg.den1);%用dstep()函數(shù)求離散模型的階躍響應(yīng)，階躍幅值為3 結(jié)果3 思考題在未考慮調(diào)節(jié)閥飽和特性時(shí)，討論一下兩個(gè)水箱液位的變化情況，工業(yè)上是否允許？討論閥位的變化情況，工業(yè)上是否能實(shí)現(xiàn)？在一開(kāi)始閥位大開(kāi)，H1，H2液位上升迅速，很快就達(dá)到預(yù)期值。但顯然不能在工業(yè)上實(shí)現(xiàn)。閥位有其本身的最大最小的限制，在仿真中出現(xiàn)的超過(guò)100%的情況在現(xiàn)實(shí)生活中不可能出現(xiàn)，因此這一部分對(duì)應(yīng)的控制效果也是無(wú)效的。第二部分 含有非線(xiàn)性環(huán)節(jié)的控制系統(tǒng)仿真一 實(shí)驗(yàn)?zāi)康? 掌握理解控制系統(tǒng)閉環(huán)仿真技術(shù)。5 掌握理解面向結(jié)構(gòu)圖的離散相似法的原理和程

23、序結(jié)構(gòu)。6 掌握理含有非線(xiàn)性環(huán)節(jié)的控制系統(tǒng)的仿真方法。二 實(shí)驗(yàn)內(nèi)容根據(jù)上面的各式，編寫(xiě)仿真程序，實(shí)現(xiàn)無(wú)擾動(dòng)時(shí)給定值階躍仿真實(shí)驗(yàn)1. 取K P = 1.78 ，T i = 85 s T = 10s，H2 S =H2set_ percent = 80, Q d = 0，tend = 700，進(jìn)行仿真實(shí)驗(yàn)，繪制響應(yīng)曲線(xiàn)。clcclear all A=2;ku=0.1/0.5;H10=1.5;H20=1.4;alpha12 = 0.25/sqrt(H10);alpha2 = 0.25/sqrt(H20);R12=2*sqrt(H10)/alpha12;R2=2*sqrt(H20)/alpha2;H1S

24、panLo=0;H2SpanLo=0;H1SpanHi=2.52;H2SpanHi=2.52; Kp=3.91/2.2;Ti=0.85*100;%Kp=3.21;%Ti=99999999999999; ad = 1/(A*R12);a1 = 1/(A*R12);a2 = 1/(A*R2);Kc=Kp/Ti;bc=Ti;Kd = 1/A;K1 = ku/A;K2 = 1/(A*R12); uc(1)=0;uv(1)=0;ud(1)=0;u1(1)=0;u2(1)=0;xc(1)=0;xv(1)=0;xd(1)=0;x1(1)=0;x2(1)=0;yc(1)=0;yv(1)=0;yd(1)=0;y

25、1(1)=0;y2(1)=0; nCounter = 70;T=10;k=1;deltaQd=0;c=0.5;H20_percent=(H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100;H2set_percent=80;tend = nCounter*T;for t=T:T:tend k=k+1; uc(k)= (H2set_percent - (y2(k-1)+H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100)/100; uv(k)=yc(k-1); ud(k)=deltaQd; if uv(k)>c yv(k)=c; end if

26、uv(k)<-c yv(k)=0; end if uv(k)<=c & uv(k)>=-c yv(k)=uv(k); end u1(k)=yv(k); u2(k)=y1(k-1); xc(k) = xc(k-1) + Kc*T*uc(k-1); yc(k)=xc(k)+bc*Kc*uc(k); xd(k) = exp(-ad*T)*xd(k-1) + Kd/ad*(1-exp(-ad*T)*ud(k);yd(k)=xd(k); x1(k) = exp(-a1*T)*x1(k-1) + K1/a1*(1-exp(-a1*T)*u1(k);y1(k)=x1(k); x2(

27、k) = exp(-a2*T)*x2(k-1) + K2/a2*(1-exp(-a2*T)*u2(k);y2(k)=x2(k);endHlevel(:,1)=(y1+H10-H1SpanLo)/(H1SpanHi-H1SpanLo)*100;Hlevel(:,2)=(y2+H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100;yv=(yv+0.5)*100;y2sp=H2set_percent*ones(size(y1'); textPositionH1=max(Hlevel(:,1); textPositionH2=max(Hlevel(:,2); H2Stea

28、dy=Hlevel(size(Hlevel(:,1),1),1)*ones(size(y1'); xmax=max(0:T:tend); xmin=0; ymax=110; ymin=50; scrsz = get(0,'ScreenSize'); gca=figure('Position',5 10 scrsz(3)-10 scrsz(4)-90) %gca=figure('Position',5 10 scrsz(3)/2 scrsz(4)/1.5) set(gca,'Color','w'); plot

29、(0:T:tend,Hlevel(:,1),'r','LineWidth',2) hold on plot(0:T:tend,Hlevel(:,2),'b','LineWidth',2) hold on plot(0:T:tend,yv,'k','LineWidth',2) hold on plot(0:T:tend,y2sp,'g','LineWidth',2) hold on plot(0:T:tend,H2Steady,'y','Line

30、Width',2) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/10 (ymax-ymin)/2+ymin-(ymax-ymin)/10,'Color','r','LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/10,' 第一個(gè)水箱的液位H1','FontSize',16) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin

31、)/6 (ymax-ymin)/2+ymin-(ymax-ymin)/6,'Color','b','LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/6,' 第二個(gè)水箱的液位H2','FontSize',16) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/4.2 (ymax-ymin)/2+ymin-(ymax-ymin)/4.2,'Color','g'

32、,'LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/4.2,' 第二個(gè)水箱的液位給定值','FontSize',16) line(tend/2 tend/2+27,(ymax-ymin)/2+ymin-(ymax-ymin)/3.2 (ymax-ymin)/2+ymin-(ymax-ymin)/3.2,'Color','k','LineWidth',6) text(tend/2+27,(ymax-ymin)/2+ymin-(ym

33、ax-ymin)/3.2,'閥位變化情況',16) axis(xmin xmax ymin ymax); text(tend/5,ymax+1.5,'實(shí)驗(yàn)三 考慮閥位飽和特性時(shí)的控制效果','FontSize',22) grid三 思考題與實(shí)驗(yàn)三相比，考慮調(diào)節(jié)閥飽和特性前后，響應(yīng)有何不同？H1 H2的液位在考慮飽和特性之后，響應(yīng)曲線(xiàn)比不考慮的時(shí)候略微平緩一些。實(shí)驗(yàn)三 采樣系統(tǒng)的仿真一 實(shí)驗(yàn)?zāi)康? 掌握理解數(shù)字控制系統(tǒng)的仿真技術(shù)。2 掌握理解增量式 PID 數(shù)字控制器的實(shí)現(xiàn)方法。二 實(shí)驗(yàn)內(nèi)容 根據(jù)上面的各式，編寫(xiě)仿真程序。取K P = 1.89 ，

34、 T i = 30 s, T d = 7.5s ，T = 10s，H2set _ percent = 80，tend = 700 ,進(jìn)行仿真實(shí)驗(yàn)，繪制仿真曲線(xiàn)。clcclear all A=2;ku=0.1/0.5;H10=1.5;H20=1.4;alpha12 = 0.25/sqrt(H10);alpha2 = 0.25/sqrt(H20);R12=2*sqrt(H10)/alpha12;R2=2*sqrt(H20)/alpha2;H1SpanLo=0;H2SpanLo=0;H1SpanHi=2.52;H2SpanHi=2.52;Kp=1.89;Ti=30;Td=10; ad = 1/(A*

35、R12);a1 = 1/(A*R12);a2 = 1/(A*R2);Kc=Kp/Ti;bc=Ti;Kd = 1/A;K1 = ku/A;K2 = 1/(A*R12); beta1=1/(a1*a2);beta3=beta1*a1/(a2-a1);beta2=-beta1-beta3; u(1)=0;u(2)=0;u(3)=0;y(1)=0;y(2)=0;y(3)=0;nCounter = 70;T=10;k=2;deltuU=0;e(1)=0;e(2)=0;e(3)=0;uc(1)=0; H20_percent=(H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100;H2

36、set_percent=80;tend = nCounter*T;for t=2*T:T:tend k=k+1; e(3)=(H2set_percent - (y(k-1)+H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100)/100; deltaU=Kp*(e(3)-e(2)+Kp*T/Ti*e(3)+Kp*Td/T*e(3)-2*e(2)+e(1); e(1)=e(2); e(2)=e(3); u(k)=u(k-1)+deltaU; y(k)=(exp(-a1*T)+exp(-a2*T)*y(k-1) - exp(-(a1+a2)*T)*y(k-2)+K1*K2*(

37、beta1+beta2+beta3)*u(k)- . K1*K2*(beta1*(exp(-a1*T)+exp(-a2*T)+beta2*(1+exp(-a2*T) . +beta3*(1+exp(-a1*T)*u(k-1)+K1*K2*(beta1*exp(-(a1+a2)*T)+beta2*exp(-a2*T)+ . beta3*exp(-a1*T)*u(k-2) ; y(k-2)=y(k-1); y(k-1)=y(k); u(k-2)=u(k-1); u(k-1)=u(k); endy=(y+H20-H2SpanLo)/(H2SpanHi-H2SpanLo)*100;y2sp=H2set_percent*ones(size(y');u=(u+0.5)*100; Hlevel(:,1)=y;Hlevel(:,2)=y;y1=y;y2=y; textPositionH1=max(Hlevel(:,1); textPositionH2=max(Hlevel(:,2); H2Steady=Hlevel(size(Hlevel(:,1),1),1)*ones(size(y1'); xmax=max(0:T:tend); xmin=0; ymax=90; ymin=50; scrsz = get(0,'ScreenSiz