光学衍射模拟：衍射积分与傅立叶变换方法对比及MATLAB实现

简介这份资源面向学习大学物理光学、需要动手验证衍射理论的学生与教师提供一套用 MATLAB 编程模拟光学衍射的完整文档。内容围绕两条技术路线展开一是衍射积分方法通过计算接收屏上的光强分布模拟单缝衍射、杨氏双缝干涉与黑白光栅衍射二是傅立叶变换方法借助屏函数功率谱模拟简单孔径、黑白光栅及正弦光栅的夫琅和费衍射并附有可直接运行的示例代码与光强分布图。资源包共 1 个 doc 文档约 247KB以文字讲解配合程序清单的形式组织便于对照理解公式与实现。目前已有 793 人学习。读者可从中获得从物理公式到 MATLAB 代码的完整映射掌握用数值计算与图像显示复现衍射条纹的方法并借助傅立叶变换思路拓展到圆孔、矩孔、三角孔等孔径的衍射模拟适合作为课程实验的参考与自学材料。1. 光学衍射模拟为什么两套数学工具能画出同一张图把一束激光打向单缝屏幕上出现的不是一条亮线而是一串明暗相间的条纹。这个现象在物理课上被反复演示但真正让人困惑的地方在于同一张衍射图样既可以用「把缝拆成无数个点光源、逐个叠加复振幅」的衍射积分算出来也可以直接对孔径函数做一次二维傅立叶变换得到。两条路径的数学起点完全不同结果却严丝合缝地对上。这份 MATLAB 资源正是围绕这个「殊途同归」展开的它用衍射积分模拟单缝、杨氏双缝和黑白光栅又用傅立叶变换方法模拟圆孔、矩孔、三角孔、黑白光栅和正弦光栅的夫琅和费衍射。搞懂这两套方法为什么等价、各自在什么条件下成立比单纯跑通一段代码有价值得多——因为前者是底层机制后者只是调用习惯。2. 衍射积分把孔径拆成点光源再叠加2.1 单缝为什么被离散成 np 个相干光源衍射积分方法的核心思路非常朴素孔径上每一点都可以看作一个次级子波源接收屏上某点的复振幅就是所有这些子波源传到该点后的相干叠加。资源里的单缝程序把缝宽a切成np51个采样点每个点对应一个子波源然后对每个屏幕采样点累加所有子波的余弦和正弦分量lam500e-9; a1e-3; D1; % 波长、缝宽、屏距 ym3*lam*D/a; ny51; % 屏幕范围取到第三级暗纹附近 yslinspace(-ym,ym,ny); % 屏幕采样点 np51; yplinspace(0,a,np); % 缝内子波源采样点 for i1:ny sinphiys(i)/D; % 小角近似下的衍射角正弦 alpha2*pi*yp*sinphi/lam; % 每个子波源的相位 sumcossum(cos(alpha)); sumsinsum(sin(alpha)); B(i,:)(sumcos^2sumsin^2)/np^2; % 归一化光强 end这里有几个参数值得停下来看。ym3*lam*D/a不是随手取的它对应单缝衍射第一级暗纹位置lam*D/a的三倍保证画出来的图能覆盖到第三级暗纹条纹结构完整。np51是离散化精度子波源越多累加结果越接近连续积分但 51 个点对于单缝这种一维问题已经足够因为被积函数是光滑的三角函数采样密度远高于奈奎斯特要求。B(i,:)(sumcos^2sumsin^2)/np^2里的除以np^2是归一化让最大光强落在 1 附近方便后续映射到 0–255 的灰度。提示sinphiys(i)/D用了小角近似sinφ≈tanφ≈y/D。屏距D1米、屏幕范围只有毫米级时这个近似误差可以忽略但如果把D缩小到厘米级就该老老实实用ys(i)/sqrt(ys(i)^2D^2)。2.2 杨氏双缝相位差才是干涉的全部双缝程序换了一个角度不再离散缝内点而是直接算两条缝到屏幕点的光程差再用I4I0cos²(Φ/2)得到光强。这个公式本身就是两束等幅相干光叠加的解析结果代码只是把Φ2π(r2-r1)/λ代进去a2e-3; D1; % 缝间距、屏距 ym5*lam*D/a; n101; yslinspace(-ym,ym,n); for i1:n r1sqrt((ys(i)-a/2).^2D^2); % 到上缝距离 r2sqrt((ys(i)a/2).^2D^2); % 到下缝距离 phi2*pi*(r2-r1)./lam; % 相位差 B(i,:)sum(4*cos(phi/2).^2); % 双缝干涉光强 end注意这里r1、r2用的是精确的欧氏距离没有做小角近似因为双缝的缝间距a2e-3和屏距D1之比虽然小但干涉条纹对相位差极其敏感近似带来的误差会直接体现在条纹位置上。ym5*lam*D/a取到第五级条纹n101的采样密度保证每条条纹有足够多的点来描出余弦包络。这段代码里sum(4*cos(phi/2).^2)的sum其实作用在标量上是原程序的一个冗余写法不影响结果但读的时候别被它误导成在做什么累加。2.3 光栅衍射两个因子的乘积从哪来黑白光栅的公式II0(sinα/α)²(sin(Nβ)/sinβ)²是单缝衍射因子和多缝干涉因子的乘积。α(πa/λ)sinφ里的a是单缝缝宽负责调制包络β(πd/λ)sinφ里的d是光栅常数负责产生尖锐的主极大。资源里取N2、a2e-4、d5*a也就是缝宽是光栅常数的五分之一这个比例直接决定了缺级出现在第几级——第 5、10、15 级主极大正好落在单缝衍射的暗纹上被调制掉。N2; a2e-4; D5; d5*a; ym2*lam*D/a; n1001; yslinspace(-ym,ym,n); for i1:n sinphiys(i)/D; alphapi*a*sinphi/lam; betapi*d*sinphi/lam; B(i,:)(sin(alpha)./alpha).^2.*(sin(N*beta)./sin(beta)).^2; B1B/max(B); endn1001比前两个程序密得多因为光栅主极大非常窄采样点少了会直接漏掉峰值画出来的条纹忽宽忽窄这是数值模拟里典型的欠采样翻车现场。sin(beta)在beta接近整数倍 π 时分母趋零但分子sin(N*beta)同时趋零极限是NMATLAB 的浮点运算在这个邻域内通常能给出接近正确的结果但如果n取得太粗、采样点正好落在奇点上会出现NaN或异常尖峰这是后面避坑章节要展开的一个点。3. 傅立叶变换方法夫琅和费衍射就是孔径的功率谱3.1 为什么 fft2 之后还要 fftshift傅立叶光学里有一个结论夫琅和费衍射的复振幅分布正比于孔径函数的傅立叶变换光强就是它的模平方。资源里的做法是读入一张孔径图像做二维离散傅立叶变换取模再移位显示aimread(aperture.bmp); % 孔径图像二维矩阵 afftfft2(a); % 二维离散傅立叶变换 aabsabs(afft); % 取模得到振幅谱 aabssfftshift(aabs); % 直流分量移到中心 imshow(aabss,[]); % 显示频谱 colormap(gray); colorbarfft2输出的矩阵零频分量在四个角上这是因为离散傅立叶变换的求和从索引 1 开始正负频率分别排在矩阵两端。fftshift交换第一、三象限和第二、四象限把零频挪到中心视觉上才和实际衍射图样一致——人眼习惯中心最亮、向外衰减的分布。这一步不是数学上必需的纯粹是为了可视化对齐直觉但不做这一步初学者很容易把角上的亮点误认为某种异常。3.2 孔径形状如何决定衍射图样资源里用同一套流程跑了圆孔、矩孔、三角孔、单缝、黑白光栅和正弦光栅。圆孔的频谱是艾里斑中心亮斑加一圈圈暗环矩孔的频谱是两个 sinc 函数的乘积呈现十字形亮纹三角孔则出现三重对称的星状结构。这些差异全部来自孔径矩阵本身的形状fft2只是忠实地把空间域的边界条件翻译成频率域的分布。换句话说你不需要为每种孔径推导新的衍射公式只要把孔径画成黑白图像傅立叶变换就自动给出对应的夫琅和费衍射图样。这正是傅立叶方法相对于衍射积分的最大优势孔径任意复杂流程不变。3.3 正弦光栅只有三个斑说明了什么资源里先用余弦函数生成一张正弦光栅图像再做傅立叶变换结果只有三个衍射斑。这个结果不是巧合正弦光栅的透过率函数只含一个空间频率分量它的频谱自然只有零频和正负一级三个峰。对比黑白光栅——方波透过率包含无穷多奇次谐波频谱上就是一系列等间距的主极大。这个对比把「傅立叶变换把空间周期性翻译成频率离散峰」这件事变得非常直观也是理解光栅衍射缺级、谐波结构的入口。4. 避坑与排查数值模拟里那些不像 bug 的 bug4.1 条纹忽宽忽窄、峰值对不上现象光栅衍射程序跑出来主极大位置和理论值差了一截或者条纹宽度不均匀。原因通常是屏幕采样点数n不够窄主极大被欠采样峰值落在两个采样点之间被削平。解决把n从 1001 继续加大或者先用解析公式算出主极大位置再在附近局部加密采样确认峰值位置收敛。4.2 出现 NaN 或异常尖峰现象光栅公式里sin(N*beta)./sin(beta)在某些点算出NaN或极大值。原因是beta恰好等于整数倍 π 时分子分母同时为零浮点运算在奇点邻域内不稳定。解决在循环里判断abs(sin(beta))是否小于阈值小于时用极限值N替代或者直接对beta加一个极小偏移避开精确奇点。4.3 频谱图四个角亮、中心暗现象做完fft2直接imshow发现亮点在四个角。原因是没有做fftshift零频分量还在矩阵角落。解决显示前加一步fftshift或者接受这个布局但心里清楚角上才是低频。这个坑几乎每个第一次用fft2做光学模拟的人都踩过。4.4 衍射图样和理论对不上、整体偏移现象圆孔衍射的艾里斑不在画面中心。原因可能是读入的孔径图像本身不在矩阵中心或者图像有非零背景比如灰度不是纯黑白引入了额外的低频分量。解决确认孔径矩阵的中心对称性必要时用padarray把孔径补零到更大矩阵的中心背景统一归一化到 0 和 1。4.5 小角近似在近距离下失效现象把屏距D改小之后衍射积分算出的条纹位置和实验对不上。原因是sinφ≈y/D的近似在y与D可比时误差显著。解决改用精确的sinφy/sqrt(y^2D^2)或者干脆把屏距保持在远场条件Da^2/λ之下让夫琅和费近似本身成立。5. 两种方法的取舍什么时候该用积分什么时候该用 FFT衍射积分和傅立叶变换方法并不是互相替代的关系它们各自对应不同的物理条件和计算代价。衍射积分从子波叠加出发物理图像清晰适合菲涅耳区——也就是屏距不够远、夫琅和费近似不成立的情形因为你可以把传播因子完整地写进每个子波的相位里不需要假设远场。代价是计算量随孔径采样点和屏幕采样点的乘积增长二维问题很容易爆内存。傅立叶变换方法则天然对应夫琅和费远场一次fft2就把整个衍射场算完复杂度是O(N²logN)孔径多复杂都不增加流程步骤但它要求远场条件近场用它会得到错误的图样。维度衍射积分傅立叶变换适用区域菲涅耳区、近场夫琅和费区、远场孔径灵活性需逐点建模任意图像直接输入计算代价采样点乘积FFT 复杂度物理直觉子波叠加直观频谱映射抽象我一般会这么理解衍射积分是「从物理过程出发」傅立叶变换是「从数学结构出发」。前者帮你建立对相位、光程差、干涉的肌肉记忆后者帮你快速探索各种孔径的衍射花样。资源里两种方法都给了完整代码价值不在于哪段代码更短而在于你可以拿同一个孔径分别跑两条路径在远场条件下对比结果是否一致——这个交叉验证本身就是最好的理解检验。从那以后我每次做数值模拟都会先用解析解或另一种方法对一遍结果确认不是采样或近似在骗我。希望帮到你。本文还有配套的精品资源点击获取