简介：本资源为信息光学课程配套的MATLAB编程实践压缩包，面向光学、光电、物理及工科高年级本科生与研究生，聚焦傅里叶光学、衍射建模与GUI交互式光学仿真等核心教学内容。压缩包共含多个.m主程序文件及可能的GUI界面脚本（如fig+callback结构），涵盖二维快速傅里叶变换（fft2）频谱分析、菲涅尔/夫琅禾费衍射数值模拟、光学图像处理与参数可视化等功能，助力理论理解与工程实现能力同步提升。资源大小22.71MB，文件总数未披露，但主体为可直接运行的MATLAB源码，无需额外依赖库，适合作为课程作业参考、实验复现或自主拓展学习。目前已有194人学习下载，提供完整可执行逻辑链：从光场建模→频域分析→界面交互→结果呈现，代码注释充分、模块划分清晰，便于逐层理解光学算法原理与MATLAB工程化表达方式。

1. 信息光学作业.zip：不是压缩包，是光场建模与衍射仿真的MATLAB实战沙盒

你手头这个标着“信息光学作业.zip”的压缩包，大概率不是某门课的PDF答案集，而是一套可直接运行、带完整注释的MATLAB光学仿真工程——它能让你在5分钟内复现菲涅尔衍射、夫琅禾费衍射、全息图生成、空间滤波器设计，甚至用纯数值方法“拍”出一张虚拟全息图。这不是教科书里的公式推演，而是把傅里叶光学、线性系统理论、复振幅传播这些抽象概念，全部落地成.m文件里可调试、可修改、可可视化的一行行代码。适合正在啃《信息光学》教材却卡在“怎么把课本公式变成图”的本科生，也适合需要快速验证光学系统响应、又不想从零搭COMSOL或Zemax的工程师。它不依赖硬件，不调用商业光学库，所有核心算法都基于MATLAB原生函数（fft2、ifft2、meshgrid、exp(1j*...)），意味着你改一个参数就能看到光场如何扭曲、聚焦、干涉——这才是信息光学该有的手感。

2. 光场建模基础：从复振幅定义到采样定理的硬约束

信息光学的本质，是把光看作一个二维复函数：$ U(x,y) = A(x,y)e^{j\phi(x,y)} $。这个复振幅承载了光强（$|U|^2$）和相位（$\arg(U)$）两个关键信息。而“信息光学作业.zip”里所有脚本的起点，就是对这个复振幅的精确离散化表达。别跳过这一步——很多人的仿真结果发虚、出现伪影、衍射环错位，根源都在这里。

2.1 空间域建模：网格尺寸与物理尺度的绑定逻辑

打开gen_aperture.m或create_object_plane.m，你会看到类似这样的初始化段：

lambda = 632.8e-9; % 波长：nm → m z = 0.5; % 传播距离：米 N = 512; % 图像尺寸（正方形） dx = 10e-6; % 像素物理尺寸：米（即采样间隔） x = dx * (-N/2 : N/2-1); y = x; [X, Y] = meshgrid(x, y);

注意：dx不是随便设的！它必须满足奈奎斯特采样定理对衍射场的约束。若你想观察距离z处的衍射图样，最大可观测空间频率为 $ f_{\max} = \frac{1}{\lambda z} \cdot \frac{D}{2} $（D为孔径直径）。对应到采样率，要求dx < lambda * z / D。实践中我一般取dx ≈ lambda * z / (2*D)作为安全余量。比如D=1mm,z=0.5m,lambda=633nm，则dx应小于约 15.8 μm —— 所以示例中10e-6是合理选择。若你盲目改成1e-3（1mm），那整个衍射图就是一堆马赛克，毫无物理意义。

2.2 频谱域建模：FFT归一化与物理量纲的还原

所有衍射计算最终都绕不开fft2和ifft2。但MATLAB的FFT默认不带物理量纲缩放，直接U_f = fft2(U)得到的是无量纲频谱，无法对应真实空间频率（单位：1/m）。必须手动补上归一化因子：

% 正向变换：空域 → 频域（需归一化以保持能量守恒） U_f = fftshift(fft2(ifftshift(U))) * dx^2; % 反向变换：频域 → 空域（同样需归一化） U_rec = ifftshift(ifft2(fftshift(U_f))) / dx^2;

为什么是dx^2？因为连续傅里叶变换定义为 $ \mathcal{F}{U(x,y)} = \iint U(x,y) e^{-j2\pi(f_x x + f_y y)} dx dy $，离散化后dx dy就是积分步长。漏掉它，你的衍射强度会随分辨率变化剧烈漂移——比如N=256时强度是100，N=1024时突然变成1.6，这就是没归一化的典型翻车现场。

2.3 衍射核构造：菲涅尔近似 vs 夫琅禾费近似的切换开关

propagate_fresnel.m和propagate_fraunhofer.m是两个核心传播函数。它们的区别不在算法复杂度，而在物理适用条件：

• 夫琅禾费衍射：适用于z >> D²/λ（远场），此时传播核简化为纯相位二次项，等价于对孔径场做一次FFT。代码极简：

  H_fraunhofer = exp(-1j * pi / (lambda * z) * (FX.^2 + FY.^2)); U_out = ifft2(fft2(U_in) .* H_fraunhofer);

• 菲涅尔衍射：适用于中间距离，传播核含三次项（实际常省略），核心是二次相位调制：

  H_fresnel = exp(1j * k * z) .* exp(1j * k / (2*z) * (X.^2 + Y.^2)); U_out = ifft2(fft2(U_in .* exp(1j * k / (2*z) * (X.^2 + Y.^2))) .* H_fresnel);

关键参数k = 2*pi/lambda必须用国际单位（m⁻¹），否则相位爆炸。我见过太多人用632.8直接代入，结果exp(1j*k*z)的指数部分动辄1e8，cos/sin完全失真。务必确认所有长度单位统一为米。

3. 全息图生成实战：从物光参考光干涉到数字再现

“信息光学作业.zip”里最惊艳的部分，是generate_hologram.m—— 它不调用任何工具箱，仅用基础MATLAB就完成离轴全息图的端到端生成与再现。这不是图像叠加，而是严格遵循 $ I(x,y) = |U_o + U_r|^2 = |U_o|^2 + |U_r|^2 + U_o U_r^* + U_o^* U_r $ 的干涉原理。

3.1 物光场构建：多目标叠加与相位编码

物光场U_object并非单个点源。脚本中通常这样构造：

% 创建三个不同位置、不同相位的点源（模拟三维物体） U1 = point_source(X, Y, 0.1e-3, 0.15e-3, lambda, z_obj); % (x0,y0) U2 = point_source(X, Y, -0.05e-3, 0.2e-3, lambda, z_obj) .* exp(1j*pi/3); U3 = point_source(X, Y, 0.12e-3, -0.08e-3, lambda, z_obj) .* exp(1j*pi/2); U_object = U1 + U2 + U3;

point_source函数内部实现的是球面波：U = (1./R) .* exp(1j*k*R)，其中R = sqrt((X-x0).^2 + (Y-y0).^2 + z_obj^2)。注意z_obj是物距（非传播距离），它决定球面波曲率——这是实现“三维感”的关键。若全设为同一z_obj，再现图就是扁平的；故意拉开z_obj差值（如 0.1m / 0.15m / 0.2m），再现时不同深度的像会出现在不同焦平面。

3.2 参考光设计：离轴角与载频控制

离轴全息的核心是让参考光与物光成夹角θ，在频域中将孪生像、零级像、物像分离开。脚本中参考光常设为倾斜平面波：

theta_x = 0.02; % 弧度，约1.15度 theta_y = 0.01; U_ref = exp(1j * k * (sin(theta_x)*X + sin(theta_y)*Y)); % 近似 sinθ≈θ

载频位置计算：该参考光在频域的峰值位于(fx, fy) = (sin(theta_x)/lambda, sin(theta_y)/lambda)。以lambda=633nm、theta_x=0.02计，fx ≈ 31600 1/m。对应FFT后的像素索引为idx_x = round(fx * dx * N) + N/2。若你的N=512,dx=10e-6，则idx_x ≈ 372—— 这就是你后续用fftshift后应在右上象限看到的亮斑位置。如果没看到，要么theta太小（载频淹没在零级），要么dx太大（频域分辨率不足）。

3.3 干涉图生成与数字再现流程

全息图H是干涉强度图，但数字再现时不能直接fft2(H)—— 因为H是实数，其频谱含共轭对称，直接FFT会混叠。标准做法是先转为复全息图（记录相位信息），再滤波：

% 1. 生成实干涉图（可保存为uint16图像） H_real = real(U_object + U_ref) .^ 2 + imag(U_object + U_ref) .^ 2; H_real = mat2gray(H_real); % 归一化到[0,1] % 2. 数字再现：提取+1级衍射（物像） H_complex = ifft2(fft2(H_real) .* filter_mask); % filter_mask是手动画的圆盘 U_recon = ifft2(fft2(H_complex) .* exp(1j * k * z_rec / (2*f) * (FX.^2 + FY.^2)));

filter_mask通常是一个半径R_mask的圆形区域，中心对准+1级载频点。R_mask要足够大以包住整个物频谱，又不能太大以免引入噪声。经验值：R_mask ≈ 0.3 * N（N为图像边长）。这是全息再现成败的玄学边界——太小切掉高频细节（边缘模糊），太大带进零级杂光（背景发灰）。

4. 空间滤波与频谱操控：阿贝成像原理的MATLAB验证

“信息光学作业.zip”里spatial_filtering.m是理解阿贝成像理论最直观的入口。它不讲“频谱面放小孔”，而是让你亲手在频域删掉某些频率分量，再看空域图像如何坍缩——这才是真正的所见即所得。

4.1 频谱面定位与可视化：为什么你总找不到正确的频谱图？

新手常犯的错误：对原始图像I直接fft2(I)就以为得到频谱。错！正确流程必须包含三步：

I_centered = imsubtract(I, mean2(I)); % 去直流（否则零频炸屏） I_fft = fftshift(fft2(ifftshift(I_centered))); % 标准化FFT流 I_spectrum = log(1 + abs(I_fft)); % 对数压缩，否则看不清

为什么必须ifftshift再fft2？因为MATLAB的fft2默认把(1,1)当作原点，但图像坐标系原点在左上角，而物理频谱原点在中心。ifftshift把图像中心移到(1,1)，fft2后再fftshift把频谱中心移回图像中心。漏掉任一环节，频谱图就是旋转/错位的。我曾为这一步调试两小时，最后发现只是少了一个ifftshift。

4.2 滤波器设计：低通/高通/带通的MATLAB实现

脚本中提供三种模板：

% 低通：保留中心圆盘 LPF = (FX.^2 + FY.^2) <= (0.1*N)^2; % 高通：切除中心 HPF = (FX.^2 + FY.^2) > (0.05*N)^2; % 方向滤波（提取45度线条）：扇形掩模 theta = atan2(FY, FX); bandpass_45 = (abs(theta - pi/4) < 0.1) | (abs(theta + 3*pi/4) < 0.1);

FX,FY由freqspace(N,'meshgrid')生成，单位是“归一化频率”（-0.5 ~ 0.5）。所以0.1*N表示半径为图像宽度10%的圆——这是经验安全值。若设为0.5*N，等于没滤波；若设为0.01*N，只剩一个白点。

4.3 阿贝实验复现：从光栅到莫尔条纹的逐级崩塌

最硬核的验证，是用grating_generator.m生成正弦光栅（周期p=50px），然后做三级滤波：

运行时你会亲眼看到：当滤掉±1级，条纹立刻从锐利变模糊；当只留0级，整张图坍缩为一片死灰。这不是理论推导，是你键盘敲出来的阿贝成像铁证。

5. 避坑指南：那些让信息光学作业跑不通的5个致命细节

信息光学仿真的脆弱性远超一般信号处理——微小的单位错、符号错、归一化漏，都会导致结果完全失真。以下是我在帮A同学调试时踩过的血泪坑，按出现频率排序：

5.1 现象：衍射图样严重偏心，主极大不在图像中心

原因：孔径函数U_aperture的定义未居中。例如用U = zeros(N); U(1:100,1:100)=1;创建方孔，其几何中心在(50,50)，但FFT要求物理中心在(N/2,N/2)。
解决：所有孔径必须用meshgrid构造，并显式偏移：

x = dx*(-N/2:N/2-1); [X,Y] = meshgrid(x,x); U_ap = (abs(X) < 0.5e-3) & (abs(Y) < 0.5e-3); % 中心在(0,0)

5.2 现象：全息再现图一片漆黑，或只有中心一个亮点

原因：U_object和U_ref的振幅量级差异过大（如U_ref是U_object的1000倍），导致干涉项 $U_o U_r^*$ 被淹没。
解决：在干涉前强制归一化：

U_object = U_object / max(abs(U_object(:))); U_ref = U_ref / max(abs(U_ref(:))); H = abs(U_object + U_ref).^2;

5.3 现象：菲涅尔传播后出现诡异环状伪影，且随z增大而加剧

原因：传播距离z与dx不匹配，违反采样定理。当z很大时，球面波曲率小，但dx过大会导致空间频率混叠。
解决：动态调整dx或N。保守策略是固定N=1024，按z增大而增大dx（降低分辨率保物理正确）：

dx = 5e-6 * (1 + z/0.5); % z每增0.5m，dx增5μm

5.4 现象：空间滤波后图像出现明显棋盘格噪声

原因：滤波器掩模filter_mask是二值图（0/1），在频域造成陡峭截断，等效于空域卷积一个sinc函数，引发吉布斯现象。
解决：改用软阈值或高斯窗：

% 替代硬截断 sigma = 0.02 * N; LPF_soft = exp(-(FX.^2 + FY.^2)/(2*sigma^2));

5.5 现象：fft2结果全是Inf或NaN

原因：复振幅中存在1/0或exp(1j*Inf)。常见于点源模型1./R中R=0（即(x0,y0)落在某个像素上）。
解决：给R加一个小量避免除零：

R = sqrt((X-x0).^2 + (Y-y0).^2 + z_obj^2) + eps('single');

6. 进阶技巧：用相位恢复反演物体三维形貌

“信息光学作业.zip”里最被低估的功能，是phase_retrieval.m—— 它实现了Gerchberg-Saxton（GS）算法，仅凭一幅衍射强度图（即|U_out|^2），就能迭代反演出丢失的相位信息arg(U_in)。这不仅是作业需求，更是数字全息显微、无透镜成像的底层引擎。

6.1 GS算法核心循环：四步迭代的物理意义

GS算法本质是在空域（已知强度|U_in|^2）和频域（已知强度|U_out|^2）之间来回投影。脚本中的主循环如下：

U = sqrt(I_in) .* exp(1j*rand(size(I_in))); % 初始化随机相位 for iter = 1:100 % 步骤1：空域约束——保持振幅，随机相位→已知振幅+当前相位 U = sqrt(I_in) .* exp(1j*angle(U)); % 步骤2：传播到频域（FFT） U_f = fftshift(fft2(ifftshift(U))); % 步骤3：频域约束——保持振幅，替换为已知衍射强度 U_f = sqrt(I_out) .* exp(1j*angle(U_f)); % 步骤4：反向传播回空域（IFFT） U = ifftshift(ifft2(fftshift(U_f))); end

关键洞察：步骤1的sqrt(I_in)是物面强度（如透射光阑的明暗分布），步骤3的sqrt(I_out)是你实测/仿真的衍射图强度。算法不关心相位，只用强度约束“逼出”唯一解。收敛标志是mean(abs(angle(U)-angle(U_prev))) < 1e-3。

6.2 收敛性优化：为什么你的GS算法永远不收敛？

GS算法 notorious 地慢且易陷局部极小。脚本中预置了三个加速技巧：

support_mask是物面非零区域的二值图，必须手动提供（如support_mask = I_in > 0.1*max(I_in(:))）。没有它，算法会在背景噪声上浪费大量迭代。

6.3 三维形貌重建：从相位到高度的标定链

GS算法输出的是物面复振幅U_in = A_in * exp(j*phi_in)。若物体是反射式（如微机电结构），phi_in直接正比于表面高度h：
$$ \phi_{in}(x,y) = \frac{4\pi}{\lambda} h(x,y) $$
因此高度图可直接计算：

h = (lambda / (4*pi)) * unwrap(angle(U_in)); % unwrap消除2π跳变 surf(X*1e3, Y*1e3, h*1e6); % 单位：mm × mm × μm

标定陷阱：unwrap必须沿正确方向（通常是dim=1行方向）。若angle(U_in)有陡峭斜坡，unwrap可能误判跳变点。我的习惯是：先imshow(angle(U_in))查看相位图，若发现大面积“条纹断裂”，就在unwrap前加smooth3(angle(U_in),'gaussian', [3 3 1])预平滑。从那以后我每次做相位恢复，都强制先imshow检查原始相位图，再决定是否平滑——这步省不得，否则重建高度全是台阶状伪影。
希望帮到你

本文还有配套的精品资源，点击获取