矩量法matlab程序设计实例
MATLAB程序设计及应用实例
c=input('请输入一个字符','s'); if c>='A' & c<='Z'
disp(setstr(abs(c)+1)); elseif c>='a'& c<='z'
disp(setstr(abs(c)-1)); elseif c>='0'& c<='9'
disp(abs(c)-abs('0')); else
disp(c); end
5.1.4 选择结构- switch语句
switch语句
其语句格式为： switch 表达式 case 值1 语句组1 case 值2 语句组2 …… case 值m 语句组m otherwise 语句组m+1 end
5.1.4 选择结构- switch语句
例 某商场对顾客所购买的商品实行打折销售，已知打折标 准，求所售商品的实际销售价格
例 矩阵乘法运算要求两矩阵的维数相容，否则会出错。 先求两矩阵的乘积，若出错，则自动转去求两矩阵的点乘
资料：23矩量法分析屏蔽带状线
作业23 用矩量法分析屏蔽带状线冯容辉一．题目：设有一长直屏蔽带状线，外导体为矩形金属柱面，矩形。
内导体为一宽度等于W的薄条的宽边和窄边分别为a和b，a b带，用矩量法分析该带状线。
屏蔽带状线的结构示意图:图1 屏蔽带状线的结构示意图二．使用MATLAB仿真可得：（1）单位电容、电感与特性阻抗与带状线宽度w的关系图2 单位电容与宽度w的关系图3 单位电感与宽度w的关系图4 单位特性阻抗与宽度w的关系三．源程序：function C=liu19(w,a,b)e0=8.854187818e-12;u0=pi*(4e-7);deltal=0.01;N1=w/deltal;N2=a/deltal;N3=b/deltal;N4=N2;N5=N3;N=N1+N2+N3+N4+N5;for m1=1:N1x1(m1)=-w/2+(m1-1/2)*deltal;y1(m1)=0;endfor m2=1:N2x2(m2)=-a/2+(m2-1/2)*deltal;y2(m2)=b/2;endfor m3=1:N3x3(m3)=a/2;y3(m3)=b/2-(m3-1/2)*deltal;endfor m4=1:N4x4(m4)=a/2-(m4-1/2)*deltal;y4(m4)=-b/2;endfor m5=1:N5x5(m5)=-a/2;y5(m5)=-b/2+(m5-1/2)*deltal;endx=[x1,x2,x3,x4,x5];y=[y1,y2,y3,y4,y5];dx=[deltal*ones(N1+N2,1);zeros(N3,1);deltal*ones(N4,1);zeros(N5,1)]; dy=[zeros(N1+N2,1);deltal*ones(N3,1);zeros(N4,1);deltal*ones(N5,1)]; c1=deltal^2;c2=c1;for m=1:Nfor n=1:Na1=(x(1)-x(n))^2+(y(1)-y(n))^2;b1=-2*(x(1)-x(n))*dx(n)+2*(y(1)-y(n))*dy(n);a2=(x(m)-x(n))^2+(y(m)-y(n))^2;b2=-2*(x(m)-x(n))*dx(n)-2*(y(m)-y(n))*dy(n);c1=deltal^2;c2=c1;p1=b1^2-4*a1*c1;if p1>=0s1(m,n)=(1/2+b1/2/c1)*log(a1+b1*1/2+c1*1/4)-2*1/2+sqrt(p1)/c1*atanh(( c1+b1)/sqrt(p1))-((-1/2+b1/2/c1)*log(a1-b1*1/2+c1*1/4)+2*1/2+sqrt(p1) /c1*atanh((-c1+b1)/sqrt(p1)));elses1(m,n)=(1/2+b1/2/c1)*log(a1+b1*1/2+c1*1/4)-2*1/2+sqrt(-p1)/c1*atan(( c1+b1)/sqrt(-p1))-((-1/2+b1/2/c1)*log(a1-b1*1/2+c1*1/4)+2*1/2+sqrt(-p 1)/c1*atan((-c1+b1)/sqrt(-p1)));endp2=b2^2-4*a2*c2;if p2>=0s2(m,n)=(1/2+b2/2/c2)*log(a2+b2*1/2+c2*1/4)-2*1/2+sqrt(p2)/c2*atanh(( c2+b2)/sqrt(p2))-((-1/2+b2/2/c2)*log(a2-b2*1/2+c2*1/4)+2*1/2+sqrt(p2) /c2*atanh((-c2+b2)/sqrt(p2)));elses2(m,n)=(1/2+b2/2/c2)*log(a2+b2*1/2+c2*1/4)-2*1/2+sqrt(-p2)/c2*atan(( c2+b2)/sqrt(-p2))-((-1/2+b2/2/c2)*log(a2-b2*1/2+c2*1/4)+2*1/2+sqrt(-p2)/c2*atan((-c2+b 2)/sqrt(-p2)));endS(m,n)=deltal*(s1(m,n)-s2(m,n));endendSS=S(2:N,:);SSS=[deltal*ones(1,N);SS];B=4*pi*e0*[ones(N1,1);zeros(N2+N3+N4+N5,1)];A=inv(SSS)*B;C=sum(A(1:N1))*deltal;clear all;clf;clc;d=1.2;e0=8.854*10^(-12);u0=4*pi*10^(-7);a=1;b=0.01;for k=1:6;w(k)=0.1+0.05*k;C(k)=liu19(w(k),a,b);L(k)=u0*e0/C(k);Z(k)=sqrt(L(k)/C(k));sumc=0;for n=1:2:1999CC(n)=4*a*sin(n*pi*w(k)/(2*a))*sinh(n*pi*b/(4*a))/((n*pi)^2*e0*cosh(n *pi*b/(2*a)));sumc=sumc+CC(n);endC1(k)=w(k)/sumc;L1(k)=u0*e0/C1(k);Z1(k)=sqrt(L1(k)/C1(k));%C1(k)=2*pi*e0/log(D(k)/d);%L1(k)=u0*e0/C1(k);%Z1(k)=sqrt(L1(k)/C1(k));%t(k)=C1(k)/C(k)endfigure(1)plot(w,C,'r-')hold onplot(w,C1,'b-')grid onlegend('½âÎö½â','¾ØÁ¿·¨')xlabel('´ø×´Ïß³¤¶Èw')ylabel('µ¥Î»³¤¶ÈµçÈÝ(F/m)')title('µ¥Î»µçÈÝËæ´ø×´Ïß³¤¶ÈwµÄ±ä»¯(a=1,b=0.01)') figure(2)plot(w,L,'r-')hold onplot(w,L1,'b-')grid onlegend('½âÎö½â','¾ØÁ¿·¨')xlabel('´ø×´Ïß³¤¶Èw')ylabel('µ¥Î»³¤¶Èµç¸Ð(H/m)')title('µ¥Î»µç¸ÐËæ´ø×´Ïß³¤¶ÈwµÄ±ä»¯(a=1,b=0.01)') figure(3)plot(w,Z,'r-')hold onplot(w,Z1,'b-')hold onlegend('½âÎö½â','¾ØÁ¿·¨')grid onxlabel('´ø×´Ïß³¤¶Èw')ylabel('ÌØÐÔ×迹(\Omega)')title('µ¥Î»ÌØÐÔ×迹Ëæ´ø×´Ïß³¤¶ÈwµÄ±ä»¯(a=1,b=0.01)')。
基于矩量法二维金属体散射(内含matlab程序)
基于矩量法的二维金属体散射计算1 问题的描述本题是用矩量法计算二维金属圆柱体的散射场，如图所示为一圆柱体和一个椭圆柱的截面，为了计算简单，选入射波为垂直z 轴入射的TM 或TE 平面波i z E i z E x22 矩量法求解过程电场积分方程问题的分析由麦克斯韦方程组H j E ωμ-=⨯∇ (1) J E j H +=⨯∇ωε (2)可得电场积分方程为''20')()(4)(ds K H J KZ E x z ρρρρ--=⎰⎰(3) 表示在圆柱表面的面电流在远处产生的总场。
设入射场为E i z ，散射场为E s z ，由金属表面的边界条件s z i z z E E E +==0 (4)得 ''20')()(4)(dl K H J KZE Cz iz ρρρρ-=⎰ (5) 离散化设入射波为)sin cos (φφy x jk iz eE +-=,将散射体截面C 分为N 份△C n ，用点匹配法对上述积分式子进行离散化， 即基函数可取{上在其它C n f ∆=10)(ρ (6)可得下列离散方程：[P]{J}={b} (7) 其中：dt y y x x K H KZ P m m C mn n))()((42220-+-=⎰∆ (8) )sin cos (i m i m y x ik m eb φφ+-= (9)当m ≠n 时，()))((42220m n m n n mn y y x x K H C KZ P -+-∆=(10) 当m=n 时 解析积分为⎥⎦⎤⎢⎣⎡⎪⎭⎫ ⎝⎛∆-∆=e C K j C ZK P n n nn 4lg 214γπ (11) 其中γ=,e=方程组的求解可用LU 分解求解方程组，即P=LU ，其中P 为可逆矩阵，L 为上三角矩阵，U 为下三角矩阵，则可利用这两个基本的三角矩阵进行求解J ，求出J 之后，就可求散射场dl e y x J KZF F y x jk z Cs ))sin()cos((),()(ϕϕφ+⎰= (12))4/3(81)(πρπρρ+-=K j e KF (13)与二维场中的散射截面2))sin()cos((2),(4)(dle y x J KZ y x jK Cz φφφσ+⎰=(14)输出结果的验证此散射问题也可用模式展开法进行求解，可用此结果对本问题进行验证。
matlab编程实例100例
1-32是：图形应用篇33-66是：界面设计篇67-84是：图形处理篇85-100是：数值分析篇实例1：三角函数曲线（1）function shili01h0=figure('toolbar','none',...'position',[198 56 350 300],...'name','实例01');h1=axes('parent',h0,...'visible','off');x=-pi:0.05:pi;y=sin(x);plot(x,y);xlabel('自变量X');ylabel('函数值Y');title('SIN( )函数曲线');grid on实例2：三角函数曲线（2）function shili02h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例02');x=-pi:0.05:pi;y=sin(x)+cos(x);plot(x,y,'-*r','linewidth',1);grid onxlabel('自变量X');ylabel('函数值Y');title('三角函数');实例3：图形的叠加function shili03h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例03');x=-pi:0.05:pi;y1=sin(x);y2=cos(x);plot(x,y1,...'-*r',...x,y2,...'--og');grid onxlabel('自变量X');ylabel('函数值Y');title('三角函数');实例4：双y轴图形的绘制function shili04h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例04');x=0:900;a=1000;b=0.005;y1=2*x;y2=cos(b*x);[haxes,hline1,hline2]=plotyy(x,y1,x,y2,'semilogy','plot'); axes(haxes(1))ylabel('semilog plot');axes(haxes(2))ylabel('linear plot');实例5：单个轴窗口显示多个图形function shili05h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例05');t=0:pi/10:2*pi;[x,y]=meshgrid(t);subplot(2,2,1)plot(sin(t),cos(t))axis equalsubplot(2,2,2)z=sin(x)-cos(y);plot(t,z)axis([0 2*pi -2 2])subplot(2,2,3)h=sin(x)+cos(y);plot(t,h)axis([0 2*pi -2 2])subplot(2,2,4)g=(sin(x).^2)-(cos(y).^2);plot(t,g)axis([0 2*pi -1 1])实例6：图形标注function shili06h0=figure('toolbar','none',...'position',[200 150 450 400],...'name','实例06');t=0:pi/10:2*pi;h=plot(t,sin(t));xlabel('t=0到2\pi','fontsize',16);ylabel('sin(t)','fontsize',16);title('\it{从0to2\pi 的正弦曲线}','fontsize',16) x=get(h,'xdata');y=get(h,'ydata');imin=find(min(y)==y);imax=find(max(y)==y);text(x(imin),y(imin),...['\leftarrow最小值=',num2str(y(imin))],...'fontsize',16)text(x(imax),y(imax),...['\leftarrow最大值=',num2str(y(imax))],...'fontsize',16)实例7：条形图形function shili07h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例07');tiao1=[562 548 224 545 41 445 745 512];tiao2=[47 48 57 58 54 52 65 48];t=0:7;bar(t,tiao1)xlabel('X轴');ylabel('TIAO1值');h1=gca;h2=axes('position',get(h1,'position'));plot(t,tiao2,'linewidth',3)set(h2,'yaxislocation','right','color','none','xticklabel',[]) 实例8：区域图形function shili08h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例08');x=91:95;profits1=[88 75 84 93 77];profits2=[51 64 54 56 68];profits3=[42 54 34 25 24];profits4=[26 38 18 15 4];area(x,profits1,'facecolor',[0.5 0.9 0.6],...'edgecolor','b',...'linewidth',3)hold onarea(x,profits2,'facecolor',[0.9 0.85 0.7],...'edgecolor','y',...'linewidth',3)hold onarea(x,profits3,'facecolor',[0.3 0.6 0.7],...'edgecolor','r',...'linewidth',3)hold onarea(x,profits4,'facecolor',[0.6 0.5 0.9],...'edgecolor','m',...'linewidth',3)hold offset(gca,'xtick',[91:95])set(gca,'layer','top')gtext('\leftarrow第一季度销量')gtext('\leftarrow第二季度销量')gtext('\leftarrow第三季度销量')gtext('\leftarrow第四季度销量')xlabel('年','fontsize',16);ylabel('销售量','fontsize',16);实例9：饼图的绘制function shili09h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例09');t=[54 21 35;68 54 35;45 25 12;48 68 45;68 54 69];x=sum(t);h=pie(x);textobjs=findobj(h,'type','text');str1=get(textobjs,{'string'});val1=get(textobjs,{'extent'});oldext=cat(1,val1{:});names={'商品一：';'商品二：';'商品三：'};str2=strcat(names,str1);set(textobjs,{'string'},str2)val2=get(textobjs,{'extent'});newext=cat(1,val2{:});offset=sign(oldext(:,1)).*(newext(:,3)-oldext(:,3))/2; pos=get(textobjs,{'position'});textpos=cat(1,pos{:});textpos(:,1)=textpos(:,1)+offset;set(textobjs,{'position'},num2cell(textpos,[3,2]))实例10：阶梯图function shili10h0=figure('toolbar','none',...'position',[200 150 450 400],...'name','实例10');a=0.01;b=0.5;t=0:10;f=exp(-a*t).*sin(b*t);stairs(t,f)hold onplot(t,f,':*')hold offglabel='函数e^{-(\alpha*t)}sin\beta*t的阶梯图'; gtext(glabel,'fontsize',16)xlabel('t=0:10','fontsize',16)axis([0 10 -1.2 1.2])实例11：枝干图function shili11h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例11');x=0:pi/20:2*pi;y1=sin(x);y2=cos(x);h1=stem(x,y1+y2);hold onh2=plot(x,y1,'^r',x,y2,'*g');hold offh3=[h1(1);h2];legend(h3,'y1+y2','y1=sin(x)','y2=cos(x)') xlabel('自变量X');ylabel('函数值Y');title('正弦函数与余弦函数的线性组合');实例12：罗盘图function shili12h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例12');winddirection=[54 24 65 84256 12 235 62125 324 34 254];windpower=[2 5 5 36 8 12 76 14 10 8];rdirection=winddirection*pi/180;[x,y]=pol2cart(rdirection,windpower); compass(x,y);desc={'风向和风力','气象台','10月1日0：00到','10月1日12：00'};gtext(desc)实例13：轮廓图function shili13h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例13');[th,r]=meshgrid((0:10:360)*pi/180,0:0.05:1); [x,y]=pol2cart(th,r);z=x+i*y;f=(z.^4-1).^(0.25);contour(x,y,abs(f),20)axis equalxlabel('实部','fontsize',16);ylabel('虚部','fontsize',16);h=polar([0 2*pi],[0 1]);delete(h)hold oncontour(x,y,abs(f),20)实例14：交互式图形function shili14h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例14');axis([0 10 0 10]);hold onx=[];y=[];n=0;disp('单击鼠标左键点取需要的点'); disp('单击鼠标右键点取最后一个点'); but=1;while but==1[xi,yi,but]=ginput(1);plot(xi,yi,'bo')n=n+1;disp('单击鼠标左键点取下一个点');x(n,1)=xi;y(n,1)=yi;endt=1:n;ts=1:0.1:n;xs=spline(t,x,ts);ys=spline(t,y,ts);plot(xs,ys,'r-');hold off实例14：交互式图形function shili14h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例14');axis([0 10 0 10]);hold onx=[];y=[];n=0;disp('单击鼠标左键点取需要的点'); disp('单击鼠标右键点取最后一个点'); but=1;while but==1[xi,yi,but]=ginput(1);plot(xi,yi,'bo')n=n+1;disp('单击鼠标左键点取下一个点');x(n,1)=xi;y(n,1)=yi;endt=1:n;ts=1:0.1:n;xs=spline(t,x,ts);ys=spline(t,y,ts);plot(xs,ys,'r-');hold off实例15：变换的傅立叶函数曲线function shili15h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例15');axis equalm=moviein(20,gcf);set(gca,'nextplot','replacechildren') h=uicontrol('style','slider','position',...[100 10 500 20],'min',1,'max',20) for j=1:20plot(fft(eye(j+16)))set(h,'value',j)m(:,j)=getframe(gcf);endclf;axes('position',[0 0 1 1]);movie(m,30)实例16：劳伦兹非线形方程的无序活动function shili15h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例15');axis equalm=moviein(20,gcf);set(gca,'nextplot','replacechildren') h=uicontrol('style','slider','position',...[100 10 500 20],'min',1,'max',20) for j=1:20plot(fft(eye(j+16)))set(h,'value',j)m(:,j)=getframe(gcf);endclf;axes('position',[0 0 1 1]);movie(m,30)实例17：填充图function shili17h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例17');t=(1:2:15)*pi/8;x=sin(t);y=cos(t);fill(x,y,'r')axis square offtext(0,0,'STOP',...'color',[1 1 1],...'fontsize',50,...'horizontalalignment','center') 例18：条形图和阶梯形图function shili18h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例18');subplot(2,2,1)x=-3:0.2:3;y=exp(-x.*x);bar(x,y)title('2-D Bar Chart')subplot(2,2,2)x=-3:0.2:3;y=exp(-x.*x);bar3(x,y,'r')title('3-D Bar Chart')subplot(2,2,3)x=-3:0.2:3;y=exp(-x.*x);stairs(x,y)title('Stair Chart')subplot(2,2,4)x=-3:0.2:3;y=exp(-x.*x);barh(x,y)title('Horizontal Bar Chart')实例19：三维曲线图function shili19h0=figure('toolbar','none',...'position',[200 150 450 400],...'name','实例19');subplot(2,1,1)x=linspace(0,2*pi);y1=sin(x);y2=cos(x);y3=sin(x)+cos(x);z1=zeros(size(x));z2=0.5*z1;z3=z1;plot3(x,y1,z1,x,y2,z2,x,y3,z3)grid onxlabel('X轴');ylabel('Y轴');zlabel('Z轴');title('Figure1:3-D Plot')subplot(2,1,2)x=linspace(0,2*pi);y1=sin(x);y2=cos(x);y3=sin(x)+cos(x);z1=zeros(size(x));z2=0.5*z1;z3=z1;plot3(x,z1,y1,x,z2,y2,x,z3,y3)grid onxlabel('X轴');ylabel('Y轴');zlabel('Z轴');title('Figure2:3-D Plot')实例20：图形的隐藏属性function shili20h0=figure('toolbar','none',...'position',[200 150 450 300],...'name','实例20');subplot(1,2,1)[x,y,z]=sphere(10);mesh(x,y,z)axis offtitle('Figure1:Opaque')hidden onsubplot(1,2,2)[x,y,z]=sphere(10);mesh(x,y,z)axis offtitle('Figure2:Transparent') hidden off实例21PEAKS函数曲线function shili21h0=figure('toolbar','none',...'position',[200 100 450 450],...'name','实例21');[x,y,z]=peaks(30);subplot(2,1,1)x=x(1,:);y=y(:,1);i=find(y>0.8&y<1.2);j=find(x>-0.6&x<0.5);z(i,j)=nan*z(i,j);surfc(x,y,z)xlabel('X轴');ylabel('Y轴');zlabel('Z轴');title('Figure1:surfc函数形成的曲面') subplot(2,1,2)x=x(1,:);y=y(:,1);i=find(y>0.8&y<1.2);j=find(x>-0.6&x<0.5);z(i,j)=nan*z(i,j);surfl(x,y,z)xlabel('X轴');ylabel('Y轴');zlabel('Z轴');title('Figure2:surfl函数形成的曲面')实例22：片状图function shili22h0=figure('toolbar','none',...'position',[200 150 550 350],...'name','实例22');subplot(1,2,1)x=rand(1,20);y=rand(1,20);z=peaks(x,y*pi);t=delaunay(x,y);trimesh(t,x,y,z)hidden offtitle('Figure1:Triangular Surface Plot'); subplot(1,2,2)x=rand(1,20);y=rand(1,20);z=peaks(x,y*pi);t=delaunay(x,y);trisurf(t,x,y,z)title('Figure1:Triangular Surface Plot'); 实例23：视角的调整function shili23h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例23');x=-5:0.5:5;[x,y]=meshgrid(x);r=sqrt(x.^2+y.^2)+eps;z=sin(r)./r;subplot(2,2,1)surf(x,y,z)xlabel('X-axis')ylabel('Y-axis')zlabel('Z-axis')title('Figure1')view(-37.5,30)subplot(2,2,2)surf(x,y,z)xlabel('X-axis')ylabel('Y-axis')zlabel('Z-axis')title('Figure2')view(-37.5+90,30)subplot(2,2,3)surf(x,y,z)xlabel('X-axis')ylabel('Y-axis')zlabel('Z-axis')title('Figure3')view(-37.5,60)subplot(2,2,4)surf(x,y,z)xlabel('X-axis')ylabel('Y-axis')zlabel('Z-axis')title('Figure4')view(180,0)实例24：向量场的绘制function shili24h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例24');subplot(2,2,1)z=peaks;ribbon(z)title('Figure1')subplot(2,2,2)[x,y,z]=peaks(15);[dx,dy]=gradient(z,0.5,0.5); contour(x,y,z,10)hold onquiver(x,y,dx,dy)hold offtitle('Figure2')subplot(2,2,3)[x,y,z]=peaks(15);[nx,ny,nz]=surfnorm(x,y,z);surf(x,y,z)hold onquiver3(x,y,z,nx,ny,nz)hold offtitle('Figure3')subplot(2,2,4)x=rand(3,5);y=rand(3,5);z=rand(3,5);c=rand(3,5);fill3(x,y,z,c)grid ontitle('Figure4')实例25：灯光定位function shili25h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例25');vert=[1 1 1;1 2 1;2 2 1;2 1 1;1 1 2;12 2;2 2 2;2 1 2];fac=[1 2 3 4;2 6 7 3;4 3 7 8;15 8 4;1 2 6 5;5 6 7 8];grid offsphere(36)h=findobj('type','surface');set(h,'facelighting','phong',...'facecolor',...'interp',...'edgecolor',[0.4 0.4 0.4],...'backfacelighting',...'lit')hold onpatch('faces',fac,'vertices',vert,...'facecolor','y');light('position',[1 3 2]);light('position',[-3 -1 3]); material shinyaxis vis3d offhold off实例26：柱状图function shili26h0=figure('toolbar','none',...'position',[200 50 450 450],...'name','实例26');subplot(2,1,1)x=[5 2 18 7 39 8 65 5 54 3 2];bar(x)xlabel('X轴');ylabel('Y轴');title('第一子图');subplot(2,1,2)y=[5 2 18 7 39 8 65 5 54 3 2];barh(y)xlabel('X轴');ylabel('Y轴');title('第二子图');实例27：设置照明方式function shili27h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例27');subplot(2,2,1)sphereshading flatcamlight leftcamlight rightlighting flatcolorbaraxis offtitle('Figure1')subplot(2,2,2)sphereshading flatcamlight leftcamlight rightlighting gouraudcolorbaraxis offtitle('Figure2')subplot(2,2,3)sphereshading interpcamlight rightcamlight leftlighting phongcolorbaraxis offtitle('Figure3')subplot(2,2,4)sphereshading flatcamlight leftcamlight rightlighting nonecolorbaraxis offtitle('Figure4')实例28：羽状图function shili28h0=figure('toolbar','none',...'position',[200 150 450 350],...'name','实例28');subplot(2,1,1)alpha=90:-10:0;r=ones(size(alpha));m=alpha*pi/180;n=r*10;[u,v]=pol2cart(m,n);feather(u,v)title('羽状图')axis([0 20 0 10])subplot(2,1,2)t=0:0.5:10;x=0.05+i;y=exp(-x*t);feather(y)title('复数矩阵的羽状图')实例29：立体透视（1）function shili29h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例29');[x,y,z]=meshgrid(-2:0.1:2,...-2:0.1:2,...-2:0.1:2);v=x.*exp(-x.^2-y.^2-z.^2);grid onfor i=-2:0.5:2;h1=surf(linspace(-2,2,20),...linspace(-2,2,20),...zeros(20)+i);rotate(h1,[1 -1 1],30)dx=get(h1,'xdata');dy=get(h1,'ydata');dz=get(h1,'zdata');delete(h1)slice(x,y,z,v,[-2 2],2,-2)hold onslice(x,y,z,v,dx,dy,dz)hold offaxis tightview(-5,10)drawnowend实例30：立体透视（2）function shili30h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例30');[x,y,z]=meshgrid(-2:0.1:2,...-2:0.1:2,...-2:0.1:2);v=x.*exp(-x.^2-y.^2-z.^2); [dx,dy,dz]=cylinder;slice(x,y,z,v,[-2 2],2,-2)for i=-2:0.2:2h=surface(dx+i,dy,dz);rotate(h,[1 0 0],90)xp=get(h,'xdata');yp=get(h,'ydata');zp=get(h,'zdata');delete(h)hold onhs=slice(x,y,z,v,xp,yp,zp);axis tightxlim([-3 3])view(-10,35)drawnowdelete(hs)hold offend实例31：表面图形function shili31h0=figure('toolbar','none',...'position',[200 150 550 250],...'name','实例31');subplot(1,2,1)x=rand(100,1)*16-8;y=rand(100,1)*16-8;r=sqrt(x.^2+y.^2)+eps;z=sin(r)./r;xlin=linspace(min(x),max(x),33); ylin=linspace(min(y),max(y),33); [X,Y]=meshgrid(xlin,ylin);Z=griddata(x,y,z,X,Y,'cubic'); mesh(X,Y,Z)axis tighthold onplot3(x,y,z,'.','Markersize',20)subplot(1,2,2)k=5;n=2^k-1;theta=pi*(-n:2:n)/n;phi=(pi/2)*(-n:2:n)'/n;X=cos(phi)*cos(theta);Y=cos(phi)*sin(theta);Z=sin(phi)*ones(size(theta));colormap([0 0 0;1 1 1])C=hadamard(2^k);surf(X,Y,Z,C)axis square实例32：沿曲线移动的小球h0=figure('toolbar','none',...'position',[198 56 408 468],...'name','实例32');h1=axes('parent',h0,...'position',[0.15 0.45 0.7 0.5],...'visible','on');t=0:pi/24:4*pi;y=sin(t);plot(t,y,'b')n=length(t);h=line('color',[0 0.5 0.5],...'linestyle','.',...'markersize',25,...'erasemode','xor');k1=uicontrol('parent',h0,...'style','pushbutton',...'position',[80 100 50 30],...'string','开始',...'callback',[...'i=1;',...'k=1;,',...'m=0;,',...'while 1,',...'if k==0,',...'break,',...'end,',...'if k~=0,',...'set(h,''xdata'',t(i),''ydata'',y(i)),',...'drawnow;,',...'i=i+1;,',...'if i>n,',...'m=m+1;,',...'i=1;,',...'end,',...'end,',...'end']);k2=uicontrol('parent',h0,...'style','pushbutton',...'position',[180 100 50 30],...'string','停止',...'callback',[...'k=0;,',...'set(e1,''string'',m),',...'p=get(h,''xdata'');,',...'q=get(h,''ydata'');,',...'set(e2,''string'',p);,',...'set(e3,''string'',q)']); k3=uicontrol('parent',h0,...'style','pushbutton',...'position',[280 100 50 30],...'string','关闭',...'callback','close');e1=uicontrol('parent',h0,...'style','edit',...'position',[60 30 60 20]);t1=uicontrol('parent',h0,...'style','text',...'string','循环次数',...'position',[60 50 60 20]);e2=uicontrol('parent',h0,...'style','edit',...'position',[180 30 50 20]); t2=uicontrol('parent',h0,...'style','text',...'string','终点的X坐标值',...'position',[155 50 100 20]); e3=uicontrol('parent',h0,...'style','edit',...'position',[300 30 50 20]); t3=uicontrol('parent',h0,...'style','text',...'string','终点的Y坐标值',...'position',[275 50 100 20]);实例33：曲线转换按钮h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例33');x=0:0.5:2*pi;y=sin(x);h=plot(x,y);grid onhuidiao=[...'if i==1,',...'i=0;,',...'y=cos(x);,',...'delete(h),',...'set(hm,''string'',''正弦函数''),',...'h=plot(x,y);,',...'grid on,',...'else if i==0,',...'i=1;,',...'y=sin(x);,',...'set(hm,''string'',''余弦函数''),',...'delete(h),',...'h=plot(x,y);,',...'grid on,',...'end,',...'end'];hm=uicontrol(gcf,'style','pushbutton',...'string','余弦函数',...'callback',huidiao);i=1;set(hm,'position',[250 20 60 20]);set(gca,'position',[0.2 0.2 0.6 0.6])title('按钮的使用')hold on实例34：栅格控制按钮h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例34');x=0:0.5:2*pi;y=sin(x);plot(x,y)huidiao1=[...'set(h_toggle2,''value'',0),',...'grid on,',...];huidiao2=[...'set(h_toggle1,''value'',0),',...'grid off,',...];h_toggle1=uicontrol(gcf,'style','togglebutton',...'string','grid on',...'value',0,...'position',[20 45 50 20],...'callback',huidiao1);h_toggle2=uicontrol(gcf,'style','togglebutton',...'string','grid off',...'value',0,...'position',[20 20 50 20],...'callback',huidiao2);set(gca,'position',[0.2 0.2 0.6 0.6])title('开关按钮的使用')实例35：编辑框的使用h0=figure('toolbar','none',...'position',[200 150 350 250],...'name','实例35');f='Please input the letter';huidiao1=[...'g=upper(f);,',...'set(h2_edit,''string'',g),',...];huidiao2=[...'g=lower(f);,',...'set(h2_edit,''string'',g),',...];h1_edit=uicontrol(gcf,'style','edit',...'position',[100 200 100 50],...'HorizontalAlignment','left',...'string','Please input the letter',...'callback','f=get(h1_edit,''string'');',...'background','w',...'max',5,...'min',1);h2_edit=uicontrol(gcf,'style','edit',...'HorizontalAlignment','left',...'position',[100 100 100 50],...'max',5,...'min',1);h1_button=uicontrol(gcf,'style','pushbutton',...'string','小写变大写',...'position',[100 45 100 20],...'callback',huidiao1);h2_button=uicontrol(gcf,'style','pushbutton',...'string','大写变小写',...'position',[100 20 100 20],...'callback',huidiao2);实例36：弹出式菜单h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例36');x=0:0.5:2*pi;y=sin(x);h=plot(x,y);grid onhm=uicontrol(gcf,'style','popupmenu',...'string',...'sin(x)|cos(x)|sin(x)+cos(x)|exp(-sin(x))',...'position',[250 20 50 20]);set(hm,'value',1)huidiao=[...'v=get(hm,''value'');,',...'switch v,',...'case 1,',...'delete(h),',...'y=sin(x);,',...'h=plot(x,y);,',...'grid on,',...'case 2,',...'delete(h),',...'y=cos(x);,',...'h=plot(x,y);,',...'grid on,',...'case 3,',...'delete(h),',...'y=sin(x)+cos(x);,',...'h=plot(x,y);,',...'grid on,',...'case 4,',...'y=exp(-sin(x));,',...'h=plot(x,y);,',...'grid on,',...'end'];set(hm,'callback',huidiao)set(gca,'position',[0.2 0.2 0.6 0.6])title('弹出式菜单的使用')实例37：滑标的使用h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例37');[x,y]=meshgrid(-8:0.5:8);r=sqrt(x.^2+y.^2)+eps;z=sin(r)./r;h0=mesh(x,y,z);h1=axes('position',...[0.2 0.2 0.5 0.5],...'visible','off');htext=uicontrol(gcf,...'units','points',...'position',[20 30 45 15],...'string','brightness',...'style','text');hslider=uicontrol(gcf,...'units','points',...'position',[10 10 300 15],...'min',-1,...'max',1,...'style','slider',...'callback',...'brighten(get(hslider,''value''))');实例38：多选菜单h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例38');[x,y]=meshgrid(-8:0.5:8);r=sqrt(x.^2+y.^2)+eps;z=sin(r)./r;h0=mesh(x,y,z);hlist=uicontrol(gcf,'style','listbox',...'string','default|spring|summer|autumn|winter',...'max',5,...'min',1,...'position',[20 20 80 100],...'callback',[...'k=get(hlist,''value'');,',...'switch k,',...'case 1,',...'colormap default,',...'case 2,',...'colormap spring,',...'case 3,',...'colormap summer,',...'case 4,',...'colormap autumn,',...'case 5,',...'colormap winter,',...'end']);实例39：菜单控制的使用h0=figure('toolbar','none',...'position',[200 150 450 250],...'name','实例39');x=0:0.5:2*pi;y=cos(x);h=plot(x,y);grid onset(gcf,'toolbar','none')hm=uimenu('label','example');huidiao1=[...'set(hm_gridon,''checked'',''on''),',...'set(hm_gridoff,''checked'',''off''),',...'grid on'];huidiao2=[...'set(hm_gridoff,''checked'',''on''),',...'set(hm_gridon,''checked'',''off''),',...'grid off'];hm_gridon=uimenu(hm,'label','grid on',...'checked','on',...'callback',huidiao1);hm_gridoff=uimenu(hm,'label','grid off',...'checked','off',...'callback',huidiao2);实例40：UIMENU菜单的应用。
Matlab编程技巧与实例分析
Matlab编程技巧与实例分析引言Matlab是一种强大的数学软件工具，广泛应用于科学、工程、金融等领域。
本文将分享一些Matlab编程的技巧和实例分析，帮助读者更好地使用和理解这个工具。
一、向量化操作在Matlab中，向量化操作是一种重要的编程技巧。
通过将循环操作转换为矢量操作，可以大大提高代码的效率。
例如，假设我们有一个包含10000个元素的向量x，我们想将其中大于0的元素设置为1，小于0的元素设置为-1，可以使用以下代码实现：```Matlabx(x > 0) = 1;x(x < 0) = -1;```通过使用向量化操作，我们避免了使用循环，并且代码更简洁。
二、矩阵运算Matlab中矩阵运算相比其他编程语言更加方便和高效。
例如，矩阵乘法可以使用`*`操作符实现，而不需要使用循环。
另外，Matlab还提供了一些内置的矩阵函数，例如`transpose`（转置矩阵）和`inv`（求逆矩阵）等。
在实例分析中，我们可以考虑一个线性回归的问题。
假设我们有一组数据集X和对应的目标值y，我们想通过线性回归拟合一个模型。
使用Matlab的矩阵运算，可以简洁地表示为：```Matlabtheta = inv(X' * X) * X' * y;```这里，`X'`表示矩阵X的转置，`inv`表示求逆矩阵。
通过这种方式，我们可以快速地得到线性回归的参数。
三、绘图和数据可视化Matlab提供了丰富的绘图和数据可视化功能，可以帮助我们更直观地理解数据和结果。
例如，我们可以使用`plot`函数绘制函数的曲线，使用`scatter`函数绘制散点图，并使用`histogram`函数绘制直方图。
在实例分析中，我们可以考虑一个简单的例子：绘制正弦函数的曲线。
可以使用以下代码实现：```Matlabx = linspace(0, 2*pi, 100);y = sin(x);plot(x, y);```通过这种方式，我们可以看到正弦函数的周期性和变化趋势。
学习matlab程序-简单示例
Matlab 编程示例.程序结构及函数作用在软件Matlab 中实现主成分分析可以采取两种方式实现：一是通过编程来实现；二是直接调用Matlab 种自带程序实现。
下面主要主要介绍利用Matlab 的矩阵计算功能编程实现主成分分析。
1程序结构主函数子函数2函数作用Cwstd.m ——用总和标准化法标准化矩阵Cwfac.m ——计算相关系数矩阵；计算特征值和特征向量；对主成分进行排序；计算各特征值贡献率；挑选主成分（累计贡献率大于85%），输出主成分个数；计算主成分载荷Cwscore.m ——计算各主成分得分、综合得分并排序Cwprint.m ——读入数据文件；调用以上三个函数并输出结果3.源程序3.1 cwstd.m 总和标准化法标准化矩阵%cwstd.m,用总和标准化法标准化矩阵function std=cwstd(vector)cwsum=sum(vector,1); %对列求和[a,b]=size(vector); %矩阵大小,a 为行数,b 为列数for i=1:afor j=1:bCwprint.m Cwstd.m Cwfac.m Cwscore.mstd(i,j)= vector(i,j)/cwsum(j);endend3.2 cwfac.m计算相关系数矩阵%cwfac.mfunction result=cwfac(vector);fprintf('相关系数矩阵:\n')std=CORRCOEF(vector) %计算相关系数矩阵fprintf('特征向量(vec)及特征值(val)：\n')[vec,val]=eig(std) %求特征值(val)及特征向量(vec)newval=diag(val) ;[y,i]=sort(newval) ; %对特征根进行排序，y为排序结果，i为索引fprintf('特征根排序：\n')for z=1:length(y)newy(z)=y(length(y)+1-z);endfprintf('%g\n',newy)rate=y/sum(y);fprintf('\n贡献率：\n')newrate=newy/sum(newy)sumrate=0;newi=[];for k=length(y):-1:1sumrate=sumrate+rate(k);newi(length(y)+1-k)=i(k);if sumrate>0.85 break;endend %记下累积贡献率大85%的特征值的序号放入newi中fprintf('主成分数：%g\n\n',length(newi));fprintf('主成分载荷：\n')for p=1:length(newi)for q=1:length(y)result(q,p)=sqrt(newval(newi(p)))*vec(q,newi(p));endend %计算载荷disp(result)3.3 cwscore.m%cwscore.m,计算得分function score=cwscore(vector1,vector2);sco=vector1*vector2;csum=sum(sco,2);[newcsum,i]=sort(-1*csum);[newi,j]=sort(i);fprintf('计算得分：\n')score=[sco,csum,j]%得分矩阵：sco为各主成分得分；csum为综合得分；j为排序结果3.4 cwprint.m%cwprint.mfunction print=cwprint(filename,a,b);%filename为文本文件文件名，a为矩阵行数(样本数)，b为矩阵列数(变量指标数)fid=fopen(filename,'r')vector=fscanf(fid,'%g',[a b]);fprintf('标准化结果如下：\n')v1=cwstd(vector)result=cwfac(v1);cwscore(v1,result);4.程序测试例题4.1原始数据中国大陆35个大城市某年的10项社会经济统计指标数据见下表。
matlab编程悬臂梁例题
matlab编程悬臂梁例题
悬臂梁是结构工程中常见的一个问题，可以用来研究梁的挠曲
和应力分布。
在MATLAB中，我们可以通过编程来模拟悬臂梁的行为。
首先，我们需要确定悬臂梁的几何形状、材料特性和受力情况。
然后，我们可以使用有限元分析或者基本的梁理论来建立模型。
在MATLAB中，我们可以使用不同的工具箱或者自己编写代码来
解决悬臂梁的例题。
首先，我们需要建立悬臂梁的几何形状和材料
特性，包括梁的长度、截面形状、材料的弹性模量和截面惯性矩等。
然后，我们可以通过应力分析公式或者有限元分析来计算悬臂梁上
的应力分布。
接着，我们可以编写MATLAB代码来绘制悬臂梁上的应
力分布图，并进一步分析梁的挠曲情况。
另外，我们还可以使用MATLAB的符号计算工具箱来求解悬臂梁
的挠曲方程，并进行参数化分析。
这样可以帮助我们更好地理解悬
臂梁的行为，比如在不同的受力情况下梁的挠曲情况会如何变化。
此外，我们还可以通过MATLAB的优化工具箱来进行悬臂梁的优化设计，比如在保证强度的前提下最小化梁的重量。
总之，MATLAB提供了丰富的工具和功能来解决悬臂梁的例题，
我们可以通过编程来模拟悬臂梁的行为，并进行多方面的分析和优化设计。
希望这些信息能够帮助你更好地理解如何在MATLAB中进行悬臂梁的编程例题。
基于matlab和矩量法的电磁仿真开发
算的函数写成基函数之间的线性组合，接着匹配方程，
再由离散的线性方程组求解出展开系数，便能由此近似
的解出未知函数。
矩量法电磁计算时常用的基函数为 RWG 基函数，
关于 RWG 基函数已有大量文献进行了详细描述，本文
采用的也是此基函数。
3 仿真应用 3.1 散射
图 2 平面波参数输入框
完成入射电磁波的设定后就可以开始计算，依据模 型网格剖分后的网格数据采用矩量法进行阻抗矩阵的 计算，得到网格的电流大小，由此计算出散射场的结果。
Lu Ju-cheng1,Yang Shuo2(1.China Airborne Missile Academy,Henan Luoyang 471009;2.First military representative office of the air force in Luoyang,Henan Luoyang 471009)
图 3 球体
图 4 方体
模型进行网格剖分，网格大小与仿真波长有关，网 格越小越好，但计算机计算能力有限制，应折中考虑。网 格剖分后模型如图 5、图 6 所示。
矩量法[5-9]是将算子方程化为矩阵方程，然后求解该
矩阵方程的方法。
对于实际的电磁散射或辐射问题，数学上可以用一
个一般的算子方程来描述，即：
L*f= g
(1)
一般要获得式(1)的精确解是非常困难的，除非 L 为
非常简单的线性算子。在计算电磁学中它们可以是空间
和时间的函数。
使用矩量法求解问题的一般步骤是，先要将所需计
摘 要：利用数值计算软件 MATLAB 开发基于矩量法的电磁仿真软件，具备建模、仿真设置、电磁散射和辐
射的计算、辐射场的绘制等功能，经验证计算结果较为可靠。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
矩量法m ａtla ｂ程序设计实例：
Ha ｌlen 方程求对称振子天线
一、条件与计算目标 已知:
对称振子天线长为Ｌ,半径为a ，且天线长度与波长得关系为,，设,半径a=0、０００0001，因此波数为。
目标：
用H ａll ｅn 方程算出半波振子、全波振子以及不同值得对应参数值。
求:（1）电流分布
（2）E 面方向图 (二维）,H 面方向图（二维）,半波振子空间方向性图(三维）
二、对称振子放置图
图1 半波振子得电流
分布
半波振子天线平行于z 轴放置，在ｘ轴与ｙ轴上得分量都为零，坐标选取方式有两种形式，一般选取图1得空间放置方
式。
图１给出了天线得电流分布情况，由图可知，当天线很细时,电流分布近似正弦分布。
三、Ha ｌｌen 方程
得解题思路
()()()()2
1
'
'
'
'
12,cos sin sin 'z z
i
z z
z
z i z k
z G z z dz c kz c kz E k z z dz j ωμ'++=-⎰⎰
对于中心馈电得偶极子，Hallen 方程为
()22'1222
('),'cos sin sin ,2L L i
L L V i z G z z dz c kz c kz k z z j η
+
--
++=
<<+⎰
脉冲函数展开与点选配，得到
()1121
,''cos sin sin ,1,2,,2n
n
N
z i
n m m m m z n V I G z z dz c kz c kz k z m N j η
+''=++=
=⋅⋅⋅∑⎰
上式可以写成 矩阵形式为
四、结果与分析 （1）电流分布
图２
不同电流
分布图
分
析：由图
2可知半
波振子
天线＝
0、5得
电流分
布最大，
馈点电流最大，时辐射电阻近似等于输入电阻,因为半波振子得输入电流正好就
是波腹电流。
(2)E面方向图（二维)
图５不同得E 面方向图（1)
分析：
(a)θ＝0时，辐射场为0。
（b)当（短振子)时,方向
函数与方向图与电流元
得近似相同。
（c）时，最大辐射方向为,主瓣随增大变窄。
后开始出现副瓣。
由图6可以瞧出。
（d）时,随增大,主瓣变窄变小,副瓣逐渐变大;继续增大,主瓣转为副瓣，而原副瓣
变为主瓣。
(如图6所示）
图6不同得E面方向图（2)
H面方向图（二维）
图７未归一化得不同
得H面方向图
图8 归一化得不同得H面方向图空间方向性图（三维)
图９半波振子得空间方向图
图１0 半波
振子得空间
剖面图
附程序：
ｃlｃ;
ｃｌｅar
all
ｃlｆ;
tｉc; %计时
lamｂda=1；
N=31；ａ＝0、000000１;％已知天线与半径
ｉi=１；
for h＝０、2：0、1：0、9
L=h＊lambda;
len=Ｌ/N;％将线分成奇数段,注意首末两端得电流为0
ｅ０=8、８54ｅ-０12；u０=4＊ｐi＊10^（—7);ｋ=２＊pi/lambda;
c=3e＋０08;ｗ＝2＊pi*c；％光速，角频率
ａta=ｓqrｔ（u０／e0）;
z(1）=—L／2+len/2；
fｏr n=2：Ｎ
z（n)=z（n—1）＋len；
end
fｏr m=１:N
ｆｏｒn＝１：Ｎ
if （m＝＝ｎ)
p（m，n）＝lｏg（lｅｎ/a）／(２＊pi）-j*k＊len/4/pi；
else
ｒ(m,ｎ）=sｑｒt((z(ｍ）－z（ｎ））＾2+a^2）；
p(ｍ，n）=leｎ＊ｅｘp（—ｊ＊ｋ*r（m，ｎ）)/（4＊pi＊ｒ（m，n)）；
enｄ
end
ｅnｄ
fｏr m=1：N
ｑ(m）=ｃos（k*z（ｍ）)；
ｓ（m)=sin(k＊z（m）);
ｔ(m）=ｓin（ｋ＊ａbs(ｚ(m)）)/(ｊ*2＊ａta）；
ｅｎd
pp=p（Ｎ+1:N^2—N）；
ｐp=reshａｐｅ（ｐp，N,N-2)；
maｔ=［pｐ，q',s']；%构造矩阵
I=mａt＼t＇；
II＝［０；I(1：N-２）；0];％加上两端零电流
Cuｒreｎt=abs（II）;
x＝liｎspace（－Ｌ/2,Ｌ／２,N）;
fｉgｕre(1)；
striｎg=［'b'，＇g',＇r'，’y'，＇c’，'ｋ’,’m＇,＇r']；
sｔrｉng1=[＇ｋo’，'ｂo’,’yo’，'co＇,＇mo’，’rｏ’，'gｏ’,'bｏ']；
plot(x,Cuｒrｅnt，sｔring(ii）,'lｉnｅwｉｄtｈ＇,1、3)；
xlabel（’L／＼lambda’),ylabeｌ(’电流分布’)；
ｇｒid on
hｏld ｏn
%legenｄ(’L=0、1\lambｄa＇,’L=0、2\lambda'，’L=0、3\lａmbｄa＇,'Ｌ=0、4\ｌambｄa'，’L＝0、5\lamｂda'，＇L＝0、6\ｌambda’，'L=０、7＼lambｄa’，’Ｌ=0、8\laｍbda＇,'L=0、9＼lamｂｄa＇,'L＝1\ｌambda')
ｌegend('L=0、1\laｍbｄa’,’L=0、3\laｍbda’，'L=0、5＼ｌambda＇，＇L＝0、7\l ａｍbda'，'L＝０、9\lambdａ'，＇L＝１、１\ｌamｂｄa’，'L＝1、3\lamｂda'，’L=１、5\lambｄa'）
Zmn=１/I((N+1）/2）;％%%％％％V=1v
theta=linｓpace（0，2＊pi，360)；
ｆor m＝１：３60
for n=1：N
F1（ｍ,n)＝II（n）、*exｐ(j*k*z（n)＊cos（m＊pi／180)）*len*sin（m ＊pｉ/180);
end
enｄ
F２=—sum（Ｆ1’)；
F=F2/max（F2）;％％％归一化
ｆｉgure（２)；
polar(ｔheta，abs（F），string（ii))；
ｔitle('E面归一化方向图’）
ｖiew（90，-90）
%lｅｇend(’L＝h\lａmbｄａ’，’Ｌ=０、3＼laｍｂdａ’，＇Ｌ=0、3\lamｂda’,'L=０、4\lambda＇，’Ｌ=0、5＼lambｄａ'，’L=0、6\lambｄa','L＝0、7\lambdａ',’Ｌ=0、8\lａmbda＇，’L＝0、9＼lamｂda'，＇Ｌ=1＼ｌaｍｂda'）
legend（'L=０、1＼lａmbdａ',’Ｌ=０、3＼ｌａmbda＇,’Ｌ=0、5\lambda'，’L=0、7\lambdａ＇，’L=0、９＼lamｂdａ＇，’L＝1、1＼ｌａmbdａ’，'L=1、3\lambdａ’,’L=1、5\lambda')
hold oｎ
fiｇure(3）
ｋk＝１;
for phi=0：pi/180:２＊ｐi
foｒn=1：Ｎ
FF(ｎ)=ＩＩ(n)＊len＊exp(i*k＊ｌen＊n＊cos(pi/２）)*sin(ｐｉ/２）；
eｎd;
FＦF(kk）=ｓｕｍ（FＦ);
kk=kk+1；
end；
phi=0：pi／180:2＊ｐｉ；
ｐolaｒ（phi，ＦFF/maｘ(abs(ＦFＦ）)，sｔrｉｎg(ii））；tｉtle（＇不同Ｌ/\lambd ａH－pｌane pattern，F（{＼theta｝，｛＼pｈi}）,\thｅta=90＇）；
legenｄ（'L＝0、1＼ｌａmbda’,'Ｌ=0、３\lａmbda',＇L＝0、5\lａmbda'，'Ｌ=0、7\ｌaｍbda’,'L=0、9\laｍbda',’Ｌ=1、1\lamｂda’，'L=１、3\lamｂda＇,'Ｌ=1、5＼lambｄa＇）
ｈolｄoｎ
ｆigurｅ(4）
ｐolaｒ(phi,FＦF／max((FＦＦ)),string(ｉi））；titｌe（'归一化H－ｐｌane ｐａttern,F（｛\theｔa｝，{\phi})，\theｔa=9０'）;
ｈｏld on
figｕｒe（５）
mｍ=1;
foｒtheｔa=0:０、0１*pｉ：pi；
for n=１：N
E(1，n)=2*pi*ｃ＊u0＊ｌｅn/（４＊pi*1)＊（exp(-ｉ*ｋ*１）*exp(i＊ｋ＊len*n ＊cos（thｅta)）*siｎ(thｅｔa））；
end
ＥE＝E＊II;
G（mm)=（4＊pi*1＾2)/atａ/ａbｓ（IＩ(（N—１)/2+1)）^2/（—real(Zmn）)＊abs(EE）^2;
ｍm＝mm+１;
end。