Matlab辅助激光光学分析与应用（ 第三版 word文件）

1、二零一四年 第三版Matlab辅助激光光学分析与应用作 者：刘良清 单 位：武汉凌云光电科技毕业院校：华中科技大学激光技术与工程学院学 历：硕士研究生研究方向：自适应光学、非线性光学、激光光学、固体激光器件、激光工业应用2014年五月 第三版前 言（序）每个人都会有自己的理想和梦想，只不过我们常常走着走着就把它们弄丢了。也许你和我一样，小时候也曾信誓旦旦地说过自己的将来的理想是要当科学家什么的。那时候不是因为我们理想很远大，其实我们根本就不知道理想是什么。我们不知道真正科学家到底是什么，也不知道除了科学家之外还有什么可以叫做理想！当然，我现在没能成为科学家，将来也没有太大的可能会成为科学家。我

2、既不是什么教授，也不是什么专家，顶多算是一个普通的技术员，普通得不能再普通了。记得研二的时候，导师还曾担心过我能不能正常毕业。整个大学过得似乎都有些浑浑噩噩，没什么明确的目标，只是盼望着大学毕业了能否找到一个温饱不愁的工作。现在我知道，这世界上根本就没有哪一样东西是真正会让你不愁的。很早，我就有想写一本书的想法。开始想写小说，结果构思来、构思去，发现好的情节似乎都被作家们写完了。突然有一个新奇的思路的时候，才发现自己文笔真的不怎地，开篇几段就写不下去了。罢了，作家之梦就这么完蛋了！研究生的最后一年，为了毕业论文和文章，慢慢自己学会了一点MATLAB编程，处理一些数据和图形。但是，后来发现之前做

3、过的一些研究，突然有一天再去看的时候，竟然怎么也看不懂。于是，我就想，干脆把自己研究过的或者是纯属娱乐研究过的东西都总结在一起，把整个思路都写下来，这样可以方便自己以后参考。既然要写，那就得像模像样地写，有多少就写多少。我不知道什么时候会写完，照这个想法，也许一辈子也写不完。因为，总会有新的东西值得我们去研究。那么就这样吧，这就是我的第一本书了。书中肯定会有各种问题和错误，很多结论只代表我个人观点，如有雷同请找专家。尽信书不是好书！目 录 TOC o 1-3 h z u HYPERLINK l _Toc389242565 Matlab辅助激光光学分析与应用 )我们可以用几行简单matlab程序

4、就可以画出具有高斯分布的电场强度，如图1.1所示，图形美观，方便对光强的分布有一个感性的视觉认识。程序代码为：clear;clc;w0=0.5;r=linspace(0,3*w0,200);eta=linspace(0,2*pi,200);rho,theta=meshgrid(r,eta);x,y=pol2cart(theta,rho);Iopt=exp(-2*rho.2/w02);surf(x,y,Iopt);shading interp;xlabel(位置 /mm);ylabel(位置 /mm);zlabel(相对强度 /a.u.);title(高斯强度分布);axis(-3*w0,3*w

5、0,-3*w0,3*w0,0,1);colorbar;colormap(hot);box on;grid off;图1.1 高斯光强分布另外，我们还可以画出高斯光束在自由传输过程中的强度变化，如图1.2所示，程序代码如下：clear;clc;lambda=1.064e-3;w0=0.5;ZR=pi*w02/lambda;z=linspace(-2*ZR,2*ZR,200);y=linspace(-4*w0,4*w0,200);py,pz=meshgrid(y,z);wz=w0*sqrt(1+(lambda*pz/pi/w02).2);Iopt=w02./wz.2.*exp(-2*py.2./w

6、z.2);surf(pz,py,Iopt);shading interp;xlabel(位置 /mm);ylabel(位置 /mm);zlabel(相对强度 /a.u.);title(高斯强度分布的传输);colorbar;colormap(hot);box on;grid off;图1.2 高斯光束自由传输强度变化以上我们以简单的例子展示了Matlab在可视化方面的强大功能，但本文不再对Matlab的基本功能和语法常识进行介绍，我们认为本书的读者已经具备了基本的Matlab编程技巧。或者说，我们所做的只是将我们的实际运用跟读者进行交流讨论，促进大家共同进步。当然，我们会在一些比较关键的地方指

7、出编程过程中需要注意的问题。1.2 波动方程当Maxwell统一了电磁理论以后，他马上意识到，波动可能是该方程组的解的形式。事实上，他希望找到一组满足波动形式的方程组，以辅助他完成找到真正的波动方程。既然已经知道了光是以波动方式传播的，基尔霍夫首先注意到了 SKIPIF 1 0 正好给出了精确的光速 SKIPIF 1 0 (之前就已经被测量过)，并且法拉第和克尔已经观测到强磁场和强电场会影响光在晶体中的传播。对初始接触Maxwell方程组的人来说，并不能一眼就看出它的解具有波动形式。但是经过适当的数学操作，我们就可以将它变为波动方程的形式。我们来推到电场E的波动方程，磁场B的波动方程的推到过程

8、是类似的。我们将方程 GOTOBUTTON ZEqnNum704203 * MERGEFORMAT REF ZEqnNum704203 * Charformat ! * MERGEFORMAT (1.3)进行卷积，可得： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 11)该方程可以由矢量微分恒等式简化： SKIPIF 1 0 MACROBUTT

9、ON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 12)卷积 SKIPIF 1 0 可由 GOTOBUTTON ZEqnNum413446 * MERGEFORMAT REF ZEqnNum413446 * Charformat ! * MERGEFORMAT (1.4)式代换，由此得到： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SE

10、Q MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 13)再由 GOTOBUTTON ZEqnNum378301 * MERGEFORMAT REF ZEqnNum378301 * Charformat ! * MERGEFORMAT (1.1)式代入上式，经过整理就可得到： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c

11、* Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 14)需要指出的是，上式中没有考虑到介质的极化。若考虑到介质的极化和实际一般光学问题中自由电荷为零的条件，上式修正为： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 15)式中， SKIPIF 1 0 为极化强度矢量。这样我们

12、得到了一般的电场传播方程，该方程在非线性光学中有很重要的地位。当光在真空中传播时，式 GOTOBUTTON ZEqnNum295979 * MERGEFORMAT REF ZEqnNum295979 * Charformat ! * MERGEFORMAT (1.15)中的右边所有项均为零，方程简化为： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMA

13、T 16)这样我们就得到了电场传播的波动方程形式。当然在有些实际问题中，式 GOTOBUTTON ZEqnNum295979 * MERGEFORMAT REF ZEqnNum295979 * Charformat ! * MERGEFORMAT (1.15)中右边的项并不是都为零，至少会有一项不为零，这与介质的性质有关。1.3 衍射考虑一个振动频率为 SKIPIF 1 0 的光场，其复振幅可以表述为 SKIPIF 1 0 ，则它也必须满足波动方程： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMA

14、T ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 17)由于(是假设)电场振幅的含时部分是显式给出的，则方程 GOTOBUTTON ZEqnNum182174 * MERGEFORMAT REF ZEqnNum182174 * Charformat ! * MERGEFORMAT (1.17)可以简化为： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arab

15、ic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 18)式中 SKIPIF 1 0 是波矢量的大小。 GOTOBUTTON ZEqnNum104057 * MERGEFORMAT REF ZEqnNum104057 * Charformat ! * MERGEFORMAT (1.18)式就是所谓的赫姆霍兹方程。如果我们忽略波动的矢量特性，而只考虑它的振幅(这里不再详细讨论其过程)，那么在标量近似下，就得到了标量赫姆霍兹方程： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTE

16、qn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 19)然后，我们考虑一束沿z轴传播的光束，它的电场复振幅写成 SKIPIF 1 0 的形式。我们将它代入标量赫姆霍兹方程 GOTOBUTTON ZEqnNum479321 * MERGEFORMAT REF ZEqnNum479321 * Charformat ! * MERGEFORMAT (1.19)式，得到： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT

17、 SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 20)在旁轴近似下，有 SKIPIF 1 0 。即是说，我们假设了电场的复振幅沿z轴传播方向是缓慢变化的，与平面波类似。但是我们允许振幅沿z轴在远大于波长量级的范围上有明显的变化。这样就得到了旁轴波方程： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabi

18、c * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 21)求解方程 GOTOBUTTON ZEqnNum233508 * MERGEFORMAT REF ZEqnNum233508 * Charformat ! * MERGEFORMAT (1.21)式，得到： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MER

19、GEFORMAT 22)于是电场的表达式为： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 23)值得一提的是，基尔霍夫早在1887年就提出了著名的菲涅耳-基尔霍夫衍射公式： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec

20、c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 24)式 GOTOBUTTON ZEqnNum599911 * MERGEFORMAT REF ZEqnNum599911 * Charformat ! * MERGEFORMAT (1.23)和 GOTOBUTTON ZEqnNum309431 * MERGEFORMAT REF ZEqnNum309431 * Charformat ! * MERGEFORMAT (1.24)在分母 SKIPIF 1 0 时具有一致性，并在指数上： SKIPIF 1 0 MACROB

21、UTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 25)同时，式 GOTOBUTTON ZEqnNum599911 * MERGEFORMAT REF ZEqnNum599911 * Charformat ! * MERGEFORMAT (1.23)是 GOTOBUTTON ZEqnNum309431 * MERGEFORMAT REF ZEqnNum309431 * Charfor

22、mat ! * MERGEFORMAT (1.24)式在满足 SKIPIF 1 0 条件下的菲尼尔旁轴近似。另外，如果进一步满足远场条件 SKIPIF 1 0 ，就得到夫琅禾费衍射近似： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 26) 小孔衍射假设光场透过一个圆柱对称的小孔，这时，孔径上的场分布可以写为： SKIPIF 1 0 MACROB

23、UTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 27)这样，二维衍射积分可以简化为一维衍射积分。将 GOTOBUTTON ZEqnNum677149 * MERGEFORMAT REF ZEqnNum677149 * Charformat ! * MERGEFORMAT (1.27)式代入到菲尼尔衍射积分公式中，得到简化衍射积分式： SKIPIF 1 0 MACROBUTTON M

24、TPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 28)对角度的积分项，我们可以借助下面的公式完成： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT

25、29)式中， SKIPIF 1 0 称为零阶Bessel函数。这样， GOTOBUTTON ZEqnNum798311 * MERGEFORMAT REF ZEqnNum798311 * Charformat ! * MERGEFORMAT (1.28)式可以简化为： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 30)式 GOTOBUTTON

26、ZEqnNum730949 * MERGEFORMAT REF ZEqnNum730949 * Charformat ! * MERGEFORMAT (1.30)中的积分项也称为 SKIPIF 1 0 的汉克尔变换。在夫琅禾费衍射近似下， SKIPIF 1 0 项等于1，积分项变为 SKIPIF 1 0 的汉克尔变换。于是夫琅禾费柱对称圆孔衍射方程为： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEq

27、n c * Arabic * MERGEFORMAT 31)虽然经过了一系列简化，然而 SKIPIF 1 0 是复振幅通常都是不确定的，即使知道了强度分布，相位分布也可能是比较难预测的。当然，也可以通过辅助手段测量强度分布和相位分布。这里，我们以平面波入射为例，讨论圆孔夫琅禾费衍射问题。这时 SKIPIF 1 0 可用常数代替，不妨设为1。利用Bessel函数的递推关系，我们可以得到解析的圆孔夫琅禾费衍射公式： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Ar

28、abic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 32)图1.3 夫琅禾费圆孔衍射 (孔径，距离1m远)即使是解析式，我们还是不能直观地感受到衍射斑的样式。下面我们利用Matlab给出夫琅禾费圆孔衍射的强度分布。程序代码如下：R=0.1;lambda=1.064e-3;k=2*pi/lambda;z=1.0e3;r=linspace(0,2*1.22*lambda/2/R*z,201);eta=linspace(0,2*pi,201);rho,theta=meshgrid(r,eta);x,y=pol2cart(theta,rho

29、);Bess=besselj(1,rho*R*k/z);Ie=4*pi2*R2*Bess.2./(rho*k).2/lambda2;surf(x,y,Ie);axis(-max(r),max(r),-max(r),max(r),0,max(Ie(:);shading interp;box on; grid off;figure;plot(x(1,:),Ie(1,:),k,x(101,:),Ie(101,:),k);程序中使用了这样一个参数 SKIPIF 1 0 ，因为平面波假设具有最小的衍射角，这个参数就称为衍射极限角，所以在衍射区域里，只取二倍衍射极限范围就可以大致看出小孔的衍射特性。并且，

30、对于解析表达，Matlab提供了Bessel函数工具箱，可以直接调用。当然，为了能够从图上看到跟实际观测相近的衍射环，可以将相对强度分布取四次开根号，如图1.3下图所示。图1.4 数值积分法计算的夫琅禾费圆孔衍射上面，我们使用的解析表达式绘图。接下来，我们采用数值积分方法直接求解 GOTOBUTTON ZEqnNum189768 * MERGEFORMAT REF ZEqnNum189768 * Charformat ! * MERGEFORMAT (1.26)式夫琅禾费圆孔衍射分布。这需要将目标平面划分网格，然后用网格格点上的复振幅代替领域内的平均振幅分布求相面上的振幅分布。为了简化计算量，

31、我们不再计算二维分布，只计算一维分布，计算结果如图1.4所示。程序代码如下：R=0.1;lambda=1.064e-3;k=2*pi/lambda;z=1.0e3;r=linspace(0,2*1.22*lambda/2/R*z,201);eta=linspace(0,2*pi,201);rho,theta=meshgrid(r,eta);x,y=pol2cart(theta,rho);r0=linspace(0,R,201);eta0=linspace(0,2*pi,201);rho0,theta0=meshgrid(r0,eta0);x0,y0=pol2cart(theta0,rho0);

32、deta=R/200*2*pi/200;E2=zeros(201,1);for gk=1:201 for m=1:200 for n=1:201 E2(gk)=E2(gk)-j/lambda/z*exp(x(1,gk)2+y(1,gk)2)/z/2+z). *j*k)*exp(j*k*(x(1,gk)*x0(m,n)+y(1,gk)*y0(m,n)/z). *deta*rho0(m,n); end endendIe=conj(E2).*E2;plot(r,Ie,k,-r,Ie,k);比较图1.3和图1.4，可以看出，数值积分法计算出的衍射斑与理论解析结果是一致的。另外，值得注意的是，从图上看到

33、，一级衍射峰值比零级衍射峰值要低很多，其他级次就更低了。那为什么我们在实验中仍然可以清晰地看到许多级次的衍射环呢？那是因为人眼对光子的感光度比较高，而实验所用的激光都具有很高的光子数密度，因此还是可以看到很多的衍射级次。对于旁轴远场条件不够满足的情况，夫琅禾费近似不成立，即使是菲涅耳近似也可能不成立，这时必须使用更为严格的菲涅耳-基尔霍夫标量衍射公式： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEq

34、n c * Arabic * MERGEFORMAT 33)式中， SKIPIF 1 0 。我们仍然以平面波的圆孔衍射为例，对 GOTOBUTTON ZEqnNum918285 * MERGEFORMAT REF ZEqnNum918285 * Charformat ! * MERGEFORMAT (1.33)式数值计算不同相面上的衍射图样，如图1.5所示。程序代码如下：%exmp1_3_1R=0.1;lambda=1.064e-3;k=2*pi/lambda;for z=0.5:5:25; r=linspace(0,2*1.22*lambda/2/R*z,201); eta=linspace

35、(0,2*pi,201); rho,theta=meshgrid(r,eta); x,y=pol2cart(theta,rho); r0=linspace(0,R,201); eta0=linspace(0,2*pi,201); rho0,theta0=meshgrid(r0,eta0); x0,y0=pol2cart(theta0,rho0); deta=R/200*2*pi/200; E2=zeros(201,1); for gk=1:201 for m=1:200 for n=1:201 Rrho=sqrt(x(1,gk)-x0(m,n)2+(y(1,gk)-y0(m,n)2+z2);

36、Rtheta=z/Rrho; E2(gk)=E2(gk)-j/lambda/2*exp(Rrho*j*k). *(1+Rtheta)/Rrho*deta*rho0(m,n); end end end Ie=conj(E2).*E2; Ie=Ie/max(Ie); plot3(z*ones(size(r),r,Ie,k,z*ones(size(r),-r,Ie,k); hold on;end图1.5 不同传输距离处的圆孔衍射样式(孔径)可见，在远场条件不满的情况下，衍射样式有较大的差异。远场条件由菲涅耳数决定： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEF

37、ORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 34)式中，a为对光束起实际限制作用元件的横向限度(半径)，L为传输距离。当F满足约为1的数量级时，菲涅耳衍射近似成立，当F满足远小于1的数量级时，夫琅禾费衍射近似成立。(a)(b)图1.6 矩孔衍射图样D= (a)z=20,(b)z=500当衍射孔为矩形孔时，同样利用 GOTOBUTTON ZEqnNum918285 * MERGEFORMAT REF ZEqnNum918285

38、* Charformat ! * MERGEFORMAT (1.33)式，可以数值求解不同相面上的衍射样式。如图1.6所示。其中左图为衍射样式的相对强度分布，右图为了突出与实际观测的近似，对强度分布取四次开根号。程序代码如下：%exmp1_3_2 矩孔衍射clear;tic;R=0.1;lambda=1.064e-3;k=2*pi/lambda;z=500;xmax=2*1.22*lambda/2/R*z;x=linspace(-xmax,xmax,61);y=x;x,y=meshgrid(x,y);x0=linspace(-R,R,61);y0=x0;x0,y0=meshgrid(x0,y0

39、);deta=(2*R/60)2;E2=zeros(61,61);for bk=1:61 for gk=1:61 for m=1:60 for n=1:60 Rrho=sqrt(x(gk,bk)-x0(m,n)2+(y(gk,bk)-y0(m,n)2+z2); Rtheta=z/Rrho; E2(gk,bk)=E2(gk,bk)-j/lambda/2*exp(Rrho*j*k)*(1+Rtheta). /Rrho*deta; end end endendIe=conj(E2).*E2;surf(x,y,Ie);shading interp;axis(-xmax,xmax,-xmax,xmax)

40、;figure;Ie=conj(E2).*E2;surf(x,y,Ie.(1/4);shading interp;axis(-xmax,xmax,-xmax,xmax);toc;另外，矩孔衍射也可以求得平面波衍射夫琅禾费近似解析表达式： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 35)则远场衍射的强度分布为： SKIPIF 1 0 MACRO

41、BUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 36)式中，a和b分别为矩孔的边长。它的衍射图样如图1.7所示，单缝衍射与此有相类似的分布，只不过去除一个自由度，这里就不再赘述。程序代码如下：%矩孔衍射的解析计算 exmp1_3_3R=0.1;lambda=1.064e-3;k=2*pi/lambda;z=1.0e3;xmax=8*1.22*lambda/2/R*z;x=lins

42、pace(-xmax,xmax,200);y=x;x,y=meshgrid(x,y);IF=sin(k*x*R/z).2.*sin(k*y*R/z).2./(k*x/2/z).2./(k*y/2/z).2/lambda2/z2;surf(x,y,IF.(1/2);colormap(hot)axis equalshading interp;图1.7 夫琅禾费矩孔衍射振幅绝对值图样 双缝衍射对于双缝夫琅禾费衍射，根据相干叠加原理，只需要将 GOTOBUTTON ZEqnNum704272 * MERGEFORMAT REF ZEqnNum704272 * Charformat ! * MERGEF

43、ORMAT (1.35)稍作调整并叠加： SKIPIF 1 0 MACROBUTTON MTPlaceRef * MERGEFORMAT SEQ MTEqn h * MERGEFORMAT ( SEQ MTSec c * Arabic * MERGEFORMAT 1. SEQ MTEqn c * Arabic * MERGEFORMAT 37)式中，a为缝宽，d为缝的中心间距。计算结果如图1.8所示。程序代码如下：%双缝衍射解析计算exmp1_3_4clear;R=0.1;d=4*R;lambda=1.064e-3;k=2*pi/lambda;z=1.0e3;xmax=3*1.22*lambda/2/R*z;x=linspace(-xmax,xmax,200);y=x;x,y=meshgrid(x,y);IF1=j*exp(j*k*(x-d/2).2/2/z+z).*sin(k*(x-d/2)*R/z)./(k*(x-d/2)/2/z)/lambda/z;IF2=j*exp(j*k*(x+d/2).2/2/z+z).*sin(k*(x+d/2)*R/z)./(k*(x+d/2)/2/z)/lambda/z;IF=IF1+IF2;surf(y,x,abs(IF);colo