尧图网页设计 
简介本资源是一套面向光学工程、物理光学及研究生科研实践的MATLAB涡旋光束仿真工具聚焦涡旋相位生成、全息衍射建模与拓扑荷调控三大核心问题适用于光学信息处理、量子通信与微纳光操控等前沿方向的教学与实验设计。压缩包为1KB的RAR格式仅含1个主程序文件.m完整实现了基于贝塞尔函数的涡旋光束相位构建、全息光栅周期性调制建模、二维傅里叶变换fft2/ifft2驱动的一阶衍射场计算以及不同拓扑荷ℓ±1,±2,…下LG模式光束的对比可视化。已有1936人学习下载程序结构简洁、注释清晰可直接运行观察相位螺旋性、衍射光斑分布与拓扑荷对衍射级次能量分配的影响是理解涡旋光束动力学特性与全息调控机制的轻量级实操入口。1. 从一束“拧巴”的光说起涡旋光束的独特魅力如果你玩过手电筒或者观察过激光笔的光斑你会发现它们投射出的光斑通常是一个均匀的亮圆点。但自然界和实验室里存在一种“不走寻常路”的光——涡旋光束。想象一下你面前有一束光它的波前不是平坦的而是像螺旋楼梯或者开瓶器的螺纹一样绕着传播轴旋转前进。这种光在垂直于传播方向的横截面上中心点是一个暗点因为所有相位的波在这里相互抵消了而光强则呈一个明亮的环状分布。这个螺旋的“拧巴”程度在物理学上用一个整数l来量化它就是拓扑荷。拓扑荷为1意味着波前像右旋螺纹为-1则是左旋螺纹为2则拧得更紧相位变化更快。这束“拧巴”的光可不仅仅是物理学家实验室里的玩具。它的核心特性是携带轨道角动量。普通的光如平面波、高斯光束只有自旋角动量对应左旋或右旋圆偏振而涡旋光束多了一份“绕轴旋转”的角动量。这个特性让它成为了一个无比强大的信息载体。在光通信领域不同拓扑荷的涡旋光束可以作为相互正交的独立信道理论上可以无限叠加从而极大地提升通信容量这就是轨道角动量复用技术。在微观操控领域比如光学镊子利用涡旋光束的环状光强分布和角动量可以稳定地捕获并旋转微米甚至纳米尺度的粒子仿佛一双无形的手在操控微观世界。此外在显微成像、量子信息处理、激光加工等领域涡旋光束都展现出了独特的优势。那么如何“制造”并研究这种神奇的光束呢在实验上通常使用空间光调制器加载计算全息图来生成而在理论研究和教学验证中MATLAB成为了绝佳的模拟工具。它强大的矩阵运算和可视化能力让我们可以在电脑上“无中生有”清晰地看到涡旋光束的相位分布、光强图案并计算其拓扑荷为深入理解其物理本质和设计实际应用方案提供了低成本、高效率的途径。接下来我将带你一步步用 MATLAB 构建一个完整的涡旋光束模拟与全息生成程序并深入探讨其中的关键细节和容易踩坑的地方。2. 构建基石模拟程序的核心原理与数学模型在动手写代码之前我们必须搞清楚要在计算机里“构造”什么。模拟一束涡旋光束本质上是生成一个二维复数矩阵这个矩阵代表了光场的复振幅分布U(x, y) A(x, y) * exp(i * φ(x, y))。其中A(x, y)是振幅分布φ(x, y)是相位分布。对于最典型的拉盖尔-高斯涡旋光束其振幅就是高斯函数而相位则包含了涡旋的核心——螺旋相位项。2.1 螺旋相位与拓扑荷的数学表达涡旋光束的螺旋相位项是模拟的关键它在极坐标(r, θ)下表达最为简洁φ_spiral l * θ这里l就是拓扑荷可以是任意整数正、负或零。θ是方位角。当l ≠ 0时相位随角度θ从 0 到 2π 线性变化l * 2π。这意味着沿着一个围绕中心的闭合路径走一圈累积的相位变化是2πl。这正是拓扑荷的物理定义相位绕奇点中心暗点一周的缠绕数。在 MATLAB 的笛卡尔网格(X, Y)上我们需要计算每个像素点对应的方位角θ atan2(Y, X)。atan2函数返回的是(-π, π]区间内的角度值这直接生成了一个从-π跳变到π的相位分布即所谓的相位包裹现象。当我们用imagesc或imshow显示这个phase矩阵时就会看到熟悉的螺旋图案。l的绝对值决定了螺旋条纹的密度符号决定了旋转方向。2.2 全息图生成原理如何用计算全息“编码”涡旋光束模拟出涡旋光束的复振幅场后如何生成能用于空间光调制器加载的全息图呢最常用的是计算全息法特别是相位型全息图。因为大多数SLM只调制相位而不调制振幅。我们的目标是将目标光场涡旋光束的相位信息φ_target编码到一个可供SLM显示的灰度图中。一种经典方法是GS算法Gerchberg-Saxton算法它是一种迭代傅里叶变换算法可以在已知目标光强通常设为均匀和约束条件下如SLM仅能调制相位求解出最优的纯相位分布。但对于像涡旋光束这样相位形式明确的场我们常用更简单的闪耀光栅法或相位叠加法。相位叠加法非常直观全息图的相位φ_hologram由两部分组成一部分是承载涡旋信息的相位φ_vortex另一部分是起到偏折作用的闪耀光栅相位φ_gratingφ_hologram mod(φ_vortex φ_grating, 2π)φ_grating 2π * (fx * X fy * Y)其中(fx, fy)是光栅的频率决定了衍射光的方向。mod(..., 2π)运算确保了相位值在[0, 2π)之间这正好对应SLM上一个周期如0到2π对应灰度0到255。将这个相位分布φ_hologram线性映射到0-255的灰度区间就得到了8位灰度全息图。注意这里有一个关键细节。atan2给出的φ_vortex范围是(-π, π]而mod(..., 2π)运算期望输入是[0, 2π)或(-∞, ∞)。直接使用atan2的结果进行mod运算在θ接近-π的跳变边界处可能会产生不连续的人工条纹。更稳健的做法是使用angle(exp(1i * l * θ))函数它直接返回(-π, π]的相位且兼容复数运算或者先计算l*θ再使用mod(l*theta, 2*pi)将其转换到[0, 2π)区间。2.3 模拟观测傅里叶变换与衍射传播生成了全息图相位分布后我们如何在模拟中“看到”它产生的涡旋光束呢这需要模拟光通过一个光学系统后的传播。最常用的方法是角谱衍射理论。在MATLAB中我们可以利用快速傅里叶变换来高效实现。基本思路是将全息图表示的相位板视为一个复振幅透射率函数t(x,y) exp(i * φ_hologram)。当一束平面波模拟为均匀振幅垂直入射后出射光场U0 t(x,y)。为了观察在某个距离z处的衍射场我们对U0做二维傅里叶变换到频域乘以角谱传递函数H(fx, fy) exp(i * 2π * z * sqrt(1/λ² - fx² - fy²))然后再做逆傅里叶变换回空域得到Uz。Uz的模平方就是观察面上的光强分布I(x,y) |Uz|²。在实际编写用于教学的简化程序时我们常常使用更直接的傅里叶变换法来模拟一个4f滤波系统将全息图放在输入面经过第一个透镜进行傅里叶变换到频谱面在频谱面上我们通常只取其中一个衍射级通过滤波再经过第二个透镜进行逆傅里叶变换得到输出面的像。在这个像面上我们应该能看到清晰的涡旋光束光强图一个亮环和相位图螺旋条纹。这种方法的代码更简洁物理图像清晰非常适合原理演示。3. 手把手实现MATLAB核心代码逐行解析理论清晰后我们进入实战环节。下面我将分模块构建一个完整的模拟程序并解释每一行代码的意图和注意事项。3.1 环境初始化与参数设置任何模拟的第一步都是定义“世界”的尺度。我们需要设定模拟的分辨率、物理尺寸和波长。%% 1. 参数设置 clear; close all; clc; % 清空环境好习惯 % 物理参数 lambda 632.8e-9; % 波长单位米这里用He-Ne激光的典型波长 k 2 * pi / lambda; % 波数 % 模拟窗口参数 N 1024; % 像素数推荐使用2的幂次便于FFT且速度快 L 0.01; % 模拟窗口的物理尺寸单位米 (10mm x 10mm) dx L / N; % 空间采样间隔 dy dx; x linspace(-L/2, L/2, N); % 生成坐标轴 y x; [X, Y] meshgrid(x, y); % 生成二维网格 [Theta, R] cart2pol(X, Y); % 转换为极坐标Theta范围[-pi, pi] % 涡旋光束参数 topological_charge 3; % 拓扑荷 l可以改为任意整数测试 beam_waist 0.001; % 高斯光束的束腰半径单位米 (1mm)关键点解析N1024选择2的幂次如51210242048能最大化FFT算法的效率。太小如128会导致采样不足图像粗糙太大如4096会显著增加计算时间对于教学演示1024是个平衡点。L0.01这个窗口尺寸需要根据光束的束腰和预期观察的衍射范围来调整。如果L太小光束可能被截断太大则采样点过于稀疏。通常L应比光束主要部分如3倍束腰大一些。cart2pol这个函数非常关键它一次性计算了网格上每一点的半径R和角度Theta。注意Theta的范围是[-π, π]这直接对应atan2(Y,X)的结果。3.2 生成理想涡旋光束场接下来我们构造一个理想的拉盖尔-高斯模式的涡旋光束。%% 2. 生成理想涡旋光束拉盖尔-高斯近似 % 高斯振幅部分 A_gaussian exp(-R.^2 / beam_waist^2); % 螺旋相位部分 % 方法1直接使用 angle(exp(1i * l * theta))得到范围在[-pi, pi]的相位 phase_vortex angle(exp(1i * topological_charge * Theta)); % 方法2使用 mod 运算得到范围在[0, 2*pi)的相位 % phase_vortex mod(topological_charge * Theta, 2*pi); % 注意当 Theta 从 -pi 跳变到 pi 时方法2在边界可能产生微小跳变方法1更稳健。 % 构建复振幅场 U_vortex A_gaussian .* exp(1i * phase_vortex); % 计算光强 I_vortex abs(U_vortex).^2;关键点解析A_gaussian这里用了最简单的高斯分布忽略了拉盖尔多项式项这模拟的是拓扑荷为l、径向指数p0的拉盖尔-高斯光束。对于p0的模式光强会是多环结构代码需要引入拉盖尔多项式laguerreL(p, abs(l), 2*R.^2/beam_waist^2)。phase_vortex的计算我强烈推荐使用angle(exp(1i * l * Theta))。虽然exp(1i*l*Theta)再取angle看起来多了一步但它完美处理了相位缠绕得到数学上严格且连续的(-π, π]相位分布避免了mod函数在-π边界可能引入的数值不连续。这种不连续在后续计算中可能被放大产生伪影。点乘.*注意是点乘这是MATLAB中对矩阵元素进行逐元素运算的运算符。*是矩阵乘法在这里用会报错。3.3 生成用于SLM的相位全息图现在我们将这个涡旋相位编码到一张全息图中。%% 3. 生成计算全息图相位型 % 添加闪耀光栅将涡旋光束衍射到特定方向与零级光分离 grating_frequency 50; % 光栅频率单位周期/米 % 生成线性相位光栅 phase_grating 2 * pi * grating_frequency * X; % 只沿X方向偏折 % 将涡旋相位与光栅相位叠加并取模2π % 使用 atan2 生成的 phase_vortex 范围是[-pi, pi]先将其平移至[0, 2pi) phase_vortex_0_to_2pi mod(phase_vortex pi, 2*pi); % 先加pi使范围变为[0, 2pi) hologram_phase mod(phase_vortex_0_to_2pi phase_grating, 2*pi); % 将相位映射到8位灰度图 (0-255) hologram_gray uint8(255 * hologram_phase / (2*pi)); % 另一种常见编码方式余弦型全息图振幅型可用于理论理解 % hologram_cosine 0.5 0.5 * cos(phase_vortex phase_grating); % hologram_gray_cos uint8(255 * hologram_cosine);关键点解析grating_frequency这个参数控制衍射角。频率越大衍射角越大。需要根据SLM的像素尺寸和光学系统的设计来选择合适的值确保衍射后的光斑能完整进入后续透镜或探测器。在模拟中可以调整它来观察衍射光斑在输出面的位置变化。phase_vortex_0_to_2pi这是一个重要的预处理步骤。因为phase_vortex来自angle()函数范围是(-π, π]。而mod(..., 2π)函数对负数取模的结果在[0, 2π)例如mod(-0.1, 2*pi)会得到一个接近2π的正数。如果直接将phase_vortex与phase_grating相加再取模在phase_vortex为负值的区域会与正值区域产生一个2π的跳变编码到灰度图上就是一条明显的黑线。因此先通过mod(phase_vortex pi, 2*pi)将其整体平移至[0, 2π)区间可以避免这个跳变生成更平滑的全息图。uint8将[0, 2π)的相位线性映射到[0, 255]的整数灰度值。这是大多数SLM驱动软件所期望的输入格式8位BMP或PNG图像。3.4 模拟光学系统与衍射重建我们使用一个简化的4f系统模型来模拟通过全息图后的光场重建。%% 4. 模拟4f滤波系统重建涡旋光束 % 假设平面波垂直入射全息图的复振幅透射率为 exp(i * hologram_phase) U_input exp(1i * hologram_phase); % 入射光为均匀振幅1 % 第一步傅里叶变换到频谱面 U_freq fft2(fftshift(U_input)); % fftshift 将零频移到中心便于观察 I_freq log(1 abs(fftshift(U_freq)).^2); % 对数显示增强对比度 % 在频谱面进行滤波选择1级衍射 % 首先找到零频中心位置 center_freq floor(N/2) 1; % 创建一个圆形滤波窗口选择1级衍射斑 filter_radius 30; % 滤波孔径半径像素 [FX, FY] meshgrid(1:N, 1:N); filter_mask ((FX - center_freq).^2 (FY - center_freq).^2) filter_radius^2; % 将滤波窗口移动到1级衍射斑的大致位置由光栅频率决定 shift_pixels round(grating_frequency * L / (1/dx)); % 根据光栅频率估算偏移像素数 filter_mask_shifted circshift(filter_mask, [0, shift_pixels]); % 应用滤波 U_freq_filtered U_freq .* fftshift(filter_mask_shifted); % 注意滤波在未fftshift的频域进行 % 第二步逆傅里叶变换回空域输出面 U_output ifft2(U_freq_filtered); I_output abs(fftshift(U_output)).^2; % 输出光强 phase_output angle(fftshift(U_output)); % 输出相位关键点解析fftshift与ifftshift这是FFT模拟光学衍射中最容易混淆的点。fft2默认输出的零频在矩阵角落。fftshift将其移到中心方便我们观察和操作频谱。但在进行滤波时滤波掩模必须与fft2输出后的未移位频谱对齐。所以流程是U_freq fft2(U_input)- 对U_freq应用移位后的滤波掩模filter_mask_shifted-U_output ifft2(U_freq_filtered)- 对U_output使用fftshift将结果移到中心显示。简单记忆fftshift只用于显示和生成与显示坐标对齐的掩模实际的滤波乘法在“原始”频域进行。filter_radius和shift_pixels滤波是关键步骤。filter_radius决定了我们选择多少高频信息太小会损失细节使光斑失真太大会引入其他级次的噪声。shift_pixels需要根据光栅频率精确计算shift N * dx * fx其中fx grating_frequency是光栅空间频率。这里N*dx L所以shift_pixels round(grating_frequency * L)。由于矩阵索引是整数需要四舍五入。如果偏移不准滤波窗口可能套不住1级衍射斑导致重建失败。对数显示log(1...)频谱的强度动态范围很大直接显示abs(U_freq).^2可能只能看到最亮的零级看不到较弱的涡旋光束级次。取对数可以压缩动态范围让所有级次都可见。3.5 结果可视化与拓扑荷验证最后我们将所有结果绘制出来并验证生成的涡旋光束是否具有正确的拓扑荷。%% 5. 可视化 figure(Position, [100, 100, 1400, 800]); % 子图1理想涡旋光束的相位和光强 subplot(2,3,1); imagesc(x*1e3, y*1e3, phase_vortex); % 坐标转换为毫米显示 axis image; colorbar; colormap hsv; % HSV色彩映射最适合显示相位周期 title([理想涡旋相位 (l, num2str(topological_charge), )]); xlabel(x (mm)); ylabel(y (mm)); subplot(2,3,4); imagesc(x*1e3, y*1e3, I_vortex); axis image; colorbar; colormap gray; title(理想涡旋光强 (高斯包络)); xlabel(x (mm)); ylabel(y (mm)); % 子图2计算全息图 subplot(2,3,2); imshow(hologram_gray); title(相位型计算全息图 (8位灰度)); % 子图3频谱面光强对数 subplot(2,3,3); imagesc(log(1 abs(fftshift(U_freq)).^2)); axis image; colorbar; colormap hot; title(频谱面光强 (对数刻度)); hold on; % 绘制滤波窗口位置 viscircles([center_freqshift_pixels, center_freq], filter_radius, Color, g, LineWidth, 1.5); hold off; % 子图4重建输出光强 subplot(2,3,5); imagesc(I_output); axis image; colorbar; colormap hot; title(重建输出光强); xlabel(像素); ylabel(像素); % 子图5重建输出相位 subplot(2,3,6); imagesc(phase_output); axis image; colorbar; colormap hsv; title(重建输出相位); xlabel(像素); ylabel(像素); % 拓扑荷验证计算相位围绕中心的环路积分 % 在输出相位图上选择一个以中心为圆心的圆环路径 [output_X, output_Y] meshgrid(1:size(phase_output,2), 1:size(phase_output,1)); center_output [floor(size(phase_output,2)/2), floor(size(phase_output,1)/2)]; radius_verify 50; % 像素半径 circle_mask ((output_X - center_output(1)).^2 (output_Y - center_output(2)).^2) radius_verify^2; phase_on_circle phase_output(circle_mask); % 由于相位是包裹的需要先解包裹再计算差分 phase_unwrapped unwrap(phase_on_circle); phase_integral sum(diff(phase_unwrapped)); estimated_charge phase_integral / (2*pi); fprintf(理论拓扑荷: %d\n, topological_charge); fprintf(通过重建相位环路积分估算的拓扑荷: %.4f\n, estimated_charge);关键点解析可视化布局使用subplot将关键步骤的结果并列显示便于对比。注意坐标轴的单位转换*1e3显示毫米和合适的色彩映射hsv用于相位hot或gray用于光强。viscircles这个函数需要Image Processing Toolbox可以方便地在频谱图上画出我们应用的滤波窗口直观展示滤波位置是否正确对准了1级衍射斑。拓扑荷验证这是检验模拟成功与否的定量步骤。我们从重建的phase_output中提取一个环形区域的相位值。由于相位是包裹的被限制在(-π, π]直接差分求和会因2π跳变而错误。unwrap函数可以沿着数据序列消除这些跳变得到连续的相位变化。然后计算相邻像素的相位差之和再除以2π就得到了环路积分的缠绕数即估算的拓扑荷。如果模拟正确这个值应该非常接近你设置的topological_charge整数。输出结果中的微小偏差如2.9998来自于数值计算误差和滤波引入的微小畸变这是可以接受的。4. 从模拟到实践关键参数影响与调试心法程序跑起来只是第一步理解每个参数如何影响结果并能在结果不理想时进行调试才是真正掌握这项技能的关键。下面分享一些我调试这类程序时积累的心得。4.1 采样与窗口避免“失真”的源头模拟的准确性首先取决于采样是否满足奈奎斯特采样定理。对于我们的螺旋相位exp(i*l*θ)其空间频率与拓扑荷l和半径r有关。在靠近中心的位置方位角θ变化剧烈空间频率最高。经验法则为了充分采样一个相位周期2π内至少要有2个像素点。螺旋相位的空间频率f_θ ≈ l / (2πr)。在最小的有效半径r_min例如光束的束腰半径或你关心的最小特征尺度处要求dx ≤ (π * r_min) / |l|。例如l5,r_min0.5mm则dx ≤ 0.314mm。我们之前设置的L10mm,N1024dx≈0.0098mm远小于要求因此采样充足。但如果l很大或L很小导致dx大就可能出现采样不足导致螺旋条纹出现锯齿状混叠失真。在频谱图上你会看到高频分量折叠回低频区域。窗口尺寸L的权衡L太小高斯光束的边缘被生硬截断会在频谱中引入高频噪声。L太大在固定像素数N下dx变大可能导致采样不足。一个实用的方法是让高斯光束的振幅在模拟窗口边缘 (L/2) 处衰减到足够小例如exp(-(L/2)^2 / waist^2) 1e-3。根据这个条件可以反推所需的L。4.2 全息图编码消除“缝线”与优化衍射效率在3.3节生成全息图时我提到了相位跳变的问题。这里再深入一下。如果你直接用mod(phase_vortex phase_grating, 2*pi)而phase_vortex是atan2输出的(-π, π]值你会发现在θ ±π这条“缝线”上phase_vortex从-π跳变到π有一个2π的跃变。虽然mod运算以2π为模但-π φ_grating和π φ_grating取模2π后的结果在大多数情况下并不相等这就导致了一条灰度不连续的线称为“2π跳变线”。加载到SLM上这条线相当于一个额外的相位阶跃会散射光降低衍射效率并引入噪声。解决方案就是我代码中采用的phase_vortex_0_to_2pi mod(phase_vortex pi, 2*pi)。这相当于将atan2的(-π, π]范围整体平移π变成(0, 2π]再取模2π将其映射到[0, 2π)。这样原来的-π点变成了0π点变成了2π取模后也是0在θ ±π处相位值都是0实现了平滑衔接消除了跳变线。衍射效率纯相位全息图的衍射效率理论上限很高。但在模拟中衍射到1级的光强占比受到光栅调制度的影响。在我们的简单叠加模型中φ_hologram mod(φ_vortex 2π*fx*X, 2π)。这实际上是一个相位型闪耀光栅其衍射效率取决于相位调制的锯齿形状的完美程度。我们的mod运算产生了完美的锯齿波因此理论上可以将大部分光能衍射到1级。在实际SLM中由于像素离散化和相位调制深度限制效率会降低。4.3 滤波与重建如何获得“干净”的涡旋光斑4f系统模拟的核心是滤波。以下几个因素会直接影响重建质量滤波窗口的位置 (shift_pixels)必须精确对准1级衍射斑的中心。你可以先运行一次不带滤波的代码显示频谱图I_freq观察1级衍射斑在x方向偏离中心多少像素然后用这个观测值来设置shift_pixels而不是仅仅依靠理论计算。因为理论计算grating_frequency * L可能因为数值离散化而不完全精确。滤波窗口的大小 (filter_radius)窗口大小是一个权衡。窗口太小相当于只允许低频通过重建的光斑会丢失细节中心暗斑可能变大甚至消失环状结构变得模糊。窗口太大会包含来自其他衍射级如0级-1级的噪声或者包含频谱的高频噪声导致重建光斑背景不干净。通常窗口半径应略大于1级衍射斑的主瓣。你可以通过逐渐增大filter_radius并观察重建光斑的变化来确定一个最佳值当光斑环状结构清晰、中心暗斑明显且继续增大半径光斑不再明显改善时即可。滤波窗口的形状我们用了圆形窗口因为它符合透镜光瞳的形状。也可以使用方形窗口但方形窗口的锐利边缘会在空域引入振铃效应吉布斯现象。一种改进是使用切趾的圆形窗口如高斯型让窗口边缘平滑过渡到零可以减少振铃但会轻微降低分辨率。调试技巧如果重建出的光强图中心不是暗斑而是一个亮斑或者相位图不是清晰的螺旋条纹首先检查滤波窗口是否真的套住了1级衍射斑。其次检查全息图生成环节的相位跳变问题是否已解决。最后可以尝试输出中间变量比如滤波后的频谱abs(fftshift(U_freq_filtered)).^2看看是否只有一个亮斑1级被保留了下来。5. 超越基础程序扩展与高级应用探索掌握了基础模拟后这个程序可以作为一个平台进行许多有趣的扩展用于研究更复杂的现象或为实验做准备。5.1 模拟部分相干涡旋光束与大气湍流理想的涡旋光束是完全相干的。但在实际中如通过大气传输的光通信链路光束会受到湍流引起的随机相位扰动导致波前畸变拓扑荷可能发生模糊甚至跳变。我们可以在生成U_vortex后叠加一个随机相位屏来模拟湍流。一种简化的方法是使用符合Kolmogorov谱的随机相位分布% 生成Kolmogorov湍流相位屏简化模型 % 注意这是一个非常简化的模型用于演示概念 r0 0.1; % 弗里德参数表征湍流强度越小湍流越强 D L; % 光束口径 % 生成随机矩阵并赋予符合特定功率谱的权重 [fx, fy] meshgrid((-N/2:N/2-1)*(1/L), (-N/2:N/2-1)*(1/L)); freq_sq fx.^2 fy.^2; freq_sq(center_freq, center_freq) 1; % 避免除零 % Kolmogorov功率谱 PSD ~ f^(-11/3) phase_spectrum (freq_sq).^(-11/6); phase_spectrum(center_freq, center_freq) 0; % 去掉零频平均相位 % 生成随机相位并乘以功率谱 random_phase randn(N) 1i*randn(N); % 复高斯随机数 turbulence_phase_screen real(ifft2(fft2(random_phase) .* fftshift(sqrt(phase_spectrum)))); % 归一化并缩放相位扰动强度 turbulence_phase_screen turbulence_phase_screen / std(turbulence_phase_screen(:)) * (D/r0)^(5/6); % 应用到光束上 U_vortex_turbulent U_vortex .* exp(1i * turbulence_phase_screen);然后用这个U_vortex_turbulent去生成全息图并进行重建。你会观察到重建出的光强环变得破碎相位条纹扭曲、断裂。通过分析大量这样的随机实现可以统计拓扑荷探测的误码率评估湍流对OAM通信的影响。5.2 生成多拓扑荷复合全息图与模式复用轨道角动量复用的核心是同时生成或探测多个不同l的涡旋光束。我们可以在一个全息图上编码多个涡旋相位每个对应不同的光栅频率从而在空间上分离它们。% 生成承载两个不同拓扑荷的全息图 l1 1; l2 3; fx1 50; % 第一个OAM信道的光栅频率 fx2 -50; % 第二个OAM信道的光栅频率负号表示衍射方向相反 % 生成两个涡旋相位 phase_vortex1 angle(exp(1i * l1 * Theta)); phase_vortex2 angle(exp(1i * l2 * Theta)); phase_vortex1_0_to_2pi mod(phase_vortex1 pi, 2*pi); phase_vortex2_0_to_2pi mod(phase_vortex2 pi, 2*pi); % 方法一叠加法简单但效率非最优 % hologram_phase_multi mod(phase_vortex1_0_to_2pi 2*pi*fx1*X ... % phase_vortex2_0_to_2pi 2*pi*fx2*X, 2*pi); % 方法二复用多个计算全息图更常用 % 分别生成两个全息图然后选择相位值例如取和或取模 % 这里采用一种简单的复用将两个相位场相加后取模 hologram_phase_multi mod(phase_vortex1_0_to_2pi phase_vortex2_0_to_2pi ... 2*pi*fx1*X, 2*pi); % 注意这里只用一个光栅频率两个模式会衍射到同一方向但重叠。 % 更复杂的复用需要设计更精巧的编码算法如GS算法优化多通道全息图。在频谱面上你会看到两套衍射斑分别位于(fx1, 0)和(fx2, 0)附近。通过使用两个独立的滤波窗口可以分别重建出l1和l3的涡旋光束。这直接模拟了OAM复用的发射端。5.3 拓扑荷检测与模态分解算法模拟在接收端我们需要判断入射的涡旋光束携带的拓扑荷是多少。一种经典的方法是使用一个与之共轭的螺旋相位板进行匹配滤波。% 假设接收到一个未知拓扑荷 l_unknown 的涡旋光束 U_unknown % 我们使用一组已知拓扑荷 l_probe 的共轭相位板进行探测 l_probe_array -5:5; % 探测的拓扑荷范围 correlation_intensity zeros(size(l_probe_array)); for idx 1:length(l_probe_array) l_probe l_probe_array(idx); % 生成共轭相位板 exp(-i * l_probe * theta) probe_phase angle(exp(-1i * l_probe * Theta)); % 注意负号 probe_field exp(1i * probe_phase); % 假设振幅均匀 % 计算与待测场的重叠积分相关 overlap sum(sum(U_unknown .* conj(probe_field))); correlation_intensity(idx) abs(overlap).^2; end % 找出相关强度最大的探测拓扑荷 [~, max_idx] max(correlation_intensity); detected_charge l_probe_array(max_idx); fprintf(探测到的拓扑荷最可能为: %d\n, detected_charge);这个简单的匹配滤波模拟表明当l_probe等于l_unknown时exp(i*l_unknown*θ) * exp(-i*l_probe*θ) 1全场相位平坦重叠积分即傅里叶透镜后的焦点光强最大。通过扫描l_probe找到光强最大的位置即可判定l_unknown。在实际实验中这可以通过一个可编程的SLM动态加载不同l_probe的相位图来实现。通过将这些扩展模块集成到你的基础程序中你可以构建一个功能强大的涡旋光束仿真研究平台从基础生成到复杂环境模拟、从单通道到多路复用、从发射到探测全方位地深化对轨道角动量这一前沿领域的理解。本文还有配套的精品资源点击获取 MATLAB模拟涡旋光束：从轨道角动量原理到计算全息生成实践
简介本资源是一套面向光学工程、物理光学及研究生科研实践的MATLAB涡旋光束仿真工具聚焦涡旋相位生成、全息衍射建模与拓扑荷调控三大核心问题适用于光学信息处理、量子通信与微纳光操控等前沿方向的教学与实验设计。压缩包为1KB的RAR格式仅含1个主程序文件.m完整实现了基于贝塞尔函数的涡旋光束相位构建、全息光栅周期性调制建模、二维傅里叶变换fft2/ifft2驱动的一阶衍射场计算以及不同拓扑荷ℓ±1,±2,…下LG模式光束的对比可视化。已有1936人学习下载程序结构简洁、注释清晰可直接运行观察相位螺旋性、衍射光斑分布与拓扑荷对衍射级次能量分配的影响是理解涡旋光束动力学特性与全息调控机制的轻量级实操入口。1. 从一束“拧巴”的光说起涡旋光束的独特魅力如果你玩过手电筒或者观察过激光笔的光斑你会发现它们投射出的光斑通常是一个均匀的亮圆点。但自然界和实验室里存在一种“不走寻常路”的光——涡旋光束。想象一下你面前有一束光它的波前不是平坦的而是像螺旋楼梯或者开瓶器的螺纹一样绕着传播轴旋转前进。这种光在垂直于传播方向的横截面上中心点是一个暗点因为所有相位的波在这里相互抵消了而光强则呈一个明亮的环状分布。这个螺旋的“拧巴”程度在物理学上用一个整数l来量化它就是拓扑荷。拓扑荷为1意味着波前像右旋螺纹为-1则是左旋螺纹为2则拧得更紧相位变化更快。这束“拧巴”的光可不仅仅是物理学家实验室里的玩具。它的核心特性是携带轨道角动量。普通的光如平面波、高斯光束只有自旋角动量对应左旋或右旋圆偏振而涡旋光束多了一份“绕轴旋转”的角动量。这个特性让它成为了一个无比强大的信息载体。在光通信领域不同拓扑荷的涡旋光束可以作为相互正交的独立信道理论上可以无限叠加从而极大地提升通信容量这就是轨道角动量复用技术。在微观操控领域比如光学镊子利用涡旋光束的环状光强分布和角动量可以稳定地捕获并旋转微米甚至纳米尺度的粒子仿佛一双无形的手在操控微观世界。此外在显微成像、量子信息处理、激光加工等领域涡旋光束都展现出了独特的优势。那么如何“制造”并研究这种神奇的光束呢在实验上通常使用空间光调制器加载计算全息图来生成而在理论研究和教学验证中MATLAB成为了绝佳的模拟工具。它强大的矩阵运算和可视化能力让我们可以在电脑上“无中生有”清晰地看到涡旋光束的相位分布、光强图案并计算其拓扑荷为深入理解其物理本质和设计实际应用方案提供了低成本、高效率的途径。接下来我将带你一步步用 MATLAB 构建一个完整的涡旋光束模拟与全息生成程序并深入探讨其中的关键细节和容易踩坑的地方。2. 构建基石模拟程序的核心原理与数学模型在动手写代码之前我们必须搞清楚要在计算机里“构造”什么。模拟一束涡旋光束本质上是生成一个二维复数矩阵这个矩阵代表了光场的复振幅分布U(x, y) A(x, y) * exp(i * φ(x, y))。其中A(x, y)是振幅分布φ(x, y)是相位分布。对于最典型的拉盖尔-高斯涡旋光束其振幅就是高斯函数而相位则包含了涡旋的核心——螺旋相位项。2.1 螺旋相位与拓扑荷的数学表达涡旋光束的螺旋相位项是模拟的关键它在极坐标(r, θ)下表达最为简洁φ_spiral l * θ这里l就是拓扑荷可以是任意整数正、负或零。θ是方位角。当l ≠ 0时相位随角度θ从 0 到 2π 线性变化l * 2π。这意味着沿着一个围绕中心的闭合路径走一圈累积的相位变化是2πl。这正是拓扑荷的物理定义相位绕奇点中心暗点一周的缠绕数。在 MATLAB 的笛卡尔网格(X, Y)上我们需要计算每个像素点对应的方位角θ atan2(Y, X)。atan2函数返回的是(-π, π]区间内的角度值这直接生成了一个从-π跳变到π的相位分布即所谓的相位包裹现象。当我们用imagesc或imshow显示这个phase矩阵时就会看到熟悉的螺旋图案。l的绝对值决定了螺旋条纹的密度符号决定了旋转方向。2.2 全息图生成原理如何用计算全息“编码”涡旋光束模拟出涡旋光束的复振幅场后如何生成能用于空间光调制器加载的全息图呢最常用的是计算全息法特别是相位型全息图。因为大多数SLM只调制相位而不调制振幅。我们的目标是将目标光场涡旋光束的相位信息φ_target编码到一个可供SLM显示的灰度图中。一种经典方法是GS算法Gerchberg-Saxton算法它是一种迭代傅里叶变换算法可以在已知目标光强通常设为均匀和约束条件下如SLM仅能调制相位求解出最优的纯相位分布。但对于像涡旋光束这样相位形式明确的场我们常用更简单的闪耀光栅法或相位叠加法。相位叠加法非常直观全息图的相位φ_hologram由两部分组成一部分是承载涡旋信息的相位φ_vortex另一部分是起到偏折作用的闪耀光栅相位φ_gratingφ_hologram mod(φ_vortex φ_grating, 2π)φ_grating 2π * (fx * X fy * Y)其中(fx, fy)是光栅的频率决定了衍射光的方向。mod(..., 2π)运算确保了相位值在[0, 2π)之间这正好对应SLM上一个周期如0到2π对应灰度0到255。将这个相位分布φ_hologram线性映射到0-255的灰度区间就得到了8位灰度全息图。注意这里有一个关键细节。atan2给出的φ_vortex范围是(-π, π]而mod(..., 2π)运算期望输入是[0, 2π)或(-∞, ∞)。直接使用atan2的结果进行mod运算在θ接近-π的跳变边界处可能会产生不连续的人工条纹。更稳健的做法是使用angle(exp(1i * l * θ))函数它直接返回(-π, π]的相位且兼容复数运算或者先计算l*θ再使用mod(l*theta, 2*pi)将其转换到[0, 2π)区间。2.3 模拟观测傅里叶变换与衍射传播生成了全息图相位分布后我们如何在模拟中“看到”它产生的涡旋光束呢这需要模拟光通过一个光学系统后的传播。最常用的方法是角谱衍射理论。在MATLAB中我们可以利用快速傅里叶变换来高效实现。基本思路是将全息图表示的相位板视为一个复振幅透射率函数t(x,y) exp(i * φ_hologram)。当一束平面波模拟为均匀振幅垂直入射后出射光场U0 t(x,y)。为了观察在某个距离z处的衍射场我们对U0做二维傅里叶变换到频域乘以角谱传递函数H(fx, fy) exp(i * 2π * z * sqrt(1/λ² - fx² - fy²))然后再做逆傅里叶变换回空域得到Uz。Uz的模平方就是观察面上的光强分布I(x,y) |Uz|²。在实际编写用于教学的简化程序时我们常常使用更直接的傅里叶变换法来模拟一个4f滤波系统将全息图放在输入面经过第一个透镜进行傅里叶变换到频谱面在频谱面上我们通常只取其中一个衍射级通过滤波再经过第二个透镜进行逆傅里叶变换得到输出面的像。在这个像面上我们应该能看到清晰的涡旋光束光强图一个亮环和相位图螺旋条纹。这种方法的代码更简洁物理图像清晰非常适合原理演示。3. 手把手实现MATLAB核心代码逐行解析理论清晰后我们进入实战环节。下面我将分模块构建一个完整的模拟程序并解释每一行代码的意图和注意事项。3.1 环境初始化与参数设置任何模拟的第一步都是定义“世界”的尺度。我们需要设定模拟的分辨率、物理尺寸和波长。%% 1. 参数设置 clear; close all; clc; % 清空环境好习惯 % 物理参数 lambda 632.8e-9; % 波长单位米这里用He-Ne激光的典型波长 k 2 * pi / lambda; % 波数 % 模拟窗口参数 N 1024; % 像素数推荐使用2的幂次便于FFT且速度快 L 0.01; % 模拟窗口的物理尺寸单位米 (10mm x 10mm) dx L / N; % 空间采样间隔 dy dx; x linspace(-L/2, L/2, N); % 生成坐标轴 y x; [X, Y] meshgrid(x, y); % 生成二维网格 [Theta, R] cart2pol(X, Y); % 转换为极坐标Theta范围[-pi, pi] % 涡旋光束参数 topological_charge 3; % 拓扑荷 l可以改为任意整数测试 beam_waist 0.001; % 高斯光束的束腰半径单位米 (1mm)关键点解析N1024选择2的幂次如51210242048能最大化FFT算法的效率。太小如128会导致采样不足图像粗糙太大如4096会显著增加计算时间对于教学演示1024是个平衡点。L0.01这个窗口尺寸需要根据光束的束腰和预期观察的衍射范围来调整。如果L太小光束可能被截断太大则采样点过于稀疏。通常L应比光束主要部分如3倍束腰大一些。cart2pol这个函数非常关键它一次性计算了网格上每一点的半径R和角度Theta。注意Theta的范围是[-π, π]这直接对应atan2(Y,X)的结果。3.2 生成理想涡旋光束场接下来我们构造一个理想的拉盖尔-高斯模式的涡旋光束。%% 2. 生成理想涡旋光束拉盖尔-高斯近似 % 高斯振幅部分 A_gaussian exp(-R.^2 / beam_waist^2); % 螺旋相位部分 % 方法1直接使用 angle(exp(1i * l * theta))得到范围在[-pi, pi]的相位 phase_vortex angle(exp(1i * topological_charge * Theta)); % 方法2使用 mod 运算得到范围在[0, 2*pi)的相位 % phase_vortex mod(topological_charge * Theta, 2*pi); % 注意当 Theta 从 -pi 跳变到 pi 时方法2在边界可能产生微小跳变方法1更稳健。 % 构建复振幅场 U_vortex A_gaussian .* exp(1i * phase_vortex); % 计算光强 I_vortex abs(U_vortex).^2;关键点解析A_gaussian这里用了最简单的高斯分布忽略了拉盖尔多项式项这模拟的是拓扑荷为l、径向指数p0的拉盖尔-高斯光束。对于p0的模式光强会是多环结构代码需要引入拉盖尔多项式laguerreL(p, abs(l), 2*R.^2/beam_waist^2)。phase_vortex的计算我强烈推荐使用angle(exp(1i * l * Theta))。虽然exp(1i*l*Theta)再取angle看起来多了一步但它完美处理了相位缠绕得到数学上严格且连续的(-π, π]相位分布避免了mod函数在-π边界可能引入的数值不连续。这种不连续在后续计算中可能被放大产生伪影。点乘.*注意是点乘这是MATLAB中对矩阵元素进行逐元素运算的运算符。*是矩阵乘法在这里用会报错。3.3 生成用于SLM的相位全息图现在我们将这个涡旋相位编码到一张全息图中。%% 3. 生成计算全息图相位型 % 添加闪耀光栅将涡旋光束衍射到特定方向与零级光分离 grating_frequency 50; % 光栅频率单位周期/米 % 生成线性相位光栅 phase_grating 2 * pi * grating_frequency * X; % 只沿X方向偏折 % 将涡旋相位与光栅相位叠加并取模2π % 使用 atan2 生成的 phase_vortex 范围是[-pi, pi]先将其平移至[0, 2pi) phase_vortex_0_to_2pi mod(phase_vortex pi, 2*pi); % 先加pi使范围变为[0, 2pi) hologram_phase mod(phase_vortex_0_to_2pi phase_grating, 2*pi); % 将相位映射到8位灰度图 (0-255) hologram_gray uint8(255 * hologram_phase / (2*pi)); % 另一种常见编码方式余弦型全息图振幅型可用于理论理解 % hologram_cosine 0.5 0.5 * cos(phase_vortex phase_grating); % hologram_gray_cos uint8(255 * hologram_cosine);关键点解析grating_frequency这个参数控制衍射角。频率越大衍射角越大。需要根据SLM的像素尺寸和光学系统的设计来选择合适的值确保衍射后的光斑能完整进入后续透镜或探测器。在模拟中可以调整它来观察衍射光斑在输出面的位置变化。phase_vortex_0_to_2pi这是一个重要的预处理步骤。因为phase_vortex来自angle()函数范围是(-π, π]。而mod(..., 2π)函数对负数取模的结果在[0, 2π)例如mod(-0.1, 2*pi)会得到一个接近2π的正数。如果直接将phase_vortex与phase_grating相加再取模在phase_vortex为负值的区域会与正值区域产生一个2π的跳变编码到灰度图上就是一条明显的黑线。因此先通过mod(phase_vortex pi, 2*pi)将其整体平移至[0, 2π)区间可以避免这个跳变生成更平滑的全息图。uint8将[0, 2π)的相位线性映射到[0, 255]的整数灰度值。这是大多数SLM驱动软件所期望的输入格式8位BMP或PNG图像。3.4 模拟光学系统与衍射重建我们使用一个简化的4f系统模型来模拟通过全息图后的光场重建。%% 4. 模拟4f滤波系统重建涡旋光束 % 假设平面波垂直入射全息图的复振幅透射率为 exp(i * hologram_phase) U_input exp(1i * hologram_phase); % 入射光为均匀振幅1 % 第一步傅里叶变换到频谱面 U_freq fft2(fftshift(U_input)); % fftshift 将零频移到中心便于观察 I_freq log(1 abs(fftshift(U_freq)).^2); % 对数显示增强对比度 % 在频谱面进行滤波选择1级衍射 % 首先找到零频中心位置 center_freq floor(N/2) 1; % 创建一个圆形滤波窗口选择1级衍射斑 filter_radius 30; % 滤波孔径半径像素 [FX, FY] meshgrid(1:N, 1:N); filter_mask ((FX - center_freq).^2 (FY - center_freq).^2) filter_radius^2; % 将滤波窗口移动到1级衍射斑的大致位置由光栅频率决定 shift_pixels round(grating_frequency * L / (1/dx)); % 根据光栅频率估算偏移像素数 filter_mask_shifted circshift(filter_mask, [0, shift_pixels]); % 应用滤波 U_freq_filtered U_freq .* fftshift(filter_mask_shifted); % 注意滤波在未fftshift的频域进行 % 第二步逆傅里叶变换回空域输出面 U_output ifft2(U_freq_filtered); I_output abs(fftshift(U_output)).^2; % 输出光强 phase_output angle(fftshift(U_output)); % 输出相位关键点解析fftshift与ifftshift这是FFT模拟光学衍射中最容易混淆的点。fft2默认输出的零频在矩阵角落。fftshift将其移到中心方便我们观察和操作频谱。但在进行滤波时滤波掩模必须与fft2输出后的未移位频谱对齐。所以流程是U_freq fft2(U_input)- 对U_freq应用移位后的滤波掩模filter_mask_shifted-U_output ifft2(U_freq_filtered)- 对U_output使用fftshift将结果移到中心显示。简单记忆fftshift只用于显示和生成与显示坐标对齐的掩模实际的滤波乘法在“原始”频域进行。filter_radius和shift_pixels滤波是关键步骤。filter_radius决定了我们选择多少高频信息太小会损失细节使光斑失真太大会引入其他级次的噪声。shift_pixels需要根据光栅频率精确计算shift N * dx * fx其中fx grating_frequency是光栅空间频率。这里N*dx L所以shift_pixels round(grating_frequency * L)。由于矩阵索引是整数需要四舍五入。如果偏移不准滤波窗口可能套不住1级衍射斑导致重建失败。对数显示log(1...)频谱的强度动态范围很大直接显示abs(U_freq).^2可能只能看到最亮的零级看不到较弱的涡旋光束级次。取对数可以压缩动态范围让所有级次都可见。3.5 结果可视化与拓扑荷验证最后我们将所有结果绘制出来并验证生成的涡旋光束是否具有正确的拓扑荷。%% 5. 可视化 figure(Position, [100, 100, 1400, 800]); % 子图1理想涡旋光束的相位和光强 subplot(2,3,1); imagesc(x*1e3, y*1e3, phase_vortex); % 坐标转换为毫米显示 axis image; colorbar; colormap hsv; % HSV色彩映射最适合显示相位周期 title([理想涡旋相位 (l, num2str(topological_charge), )]); xlabel(x (mm)); ylabel(y (mm)); subplot(2,3,4); imagesc(x*1e3, y*1e3, I_vortex); axis image; colorbar; colormap gray; title(理想涡旋光强 (高斯包络)); xlabel(x (mm)); ylabel(y (mm)); % 子图2计算全息图 subplot(2,3,2); imshow(hologram_gray); title(相位型计算全息图 (8位灰度)); % 子图3频谱面光强对数 subplot(2,3,3); imagesc(log(1 abs(fftshift(U_freq)).^2)); axis image; colorbar; colormap hot; title(频谱面光强 (对数刻度)); hold on; % 绘制滤波窗口位置 viscircles([center_freqshift_pixels, center_freq], filter_radius, Color, g, LineWidth, 1.5); hold off; % 子图4重建输出光强 subplot(2,3,5); imagesc(I_output); axis image; colorbar; colormap hot; title(重建输出光强); xlabel(像素); ylabel(像素); % 子图5重建输出相位 subplot(2,3,6); imagesc(phase_output); axis image; colorbar; colormap hsv; title(重建输出相位); xlabel(像素); ylabel(像素); % 拓扑荷验证计算相位围绕中心的环路积分 % 在输出相位图上选择一个以中心为圆心的圆环路径 [output_X, output_Y] meshgrid(1:size(phase_output,2), 1:size(phase_output,1)); center_output [floor(size(phase_output,2)/2), floor(size(phase_output,1)/2)]; radius_verify 50; % 像素半径 circle_mask ((output_X - center_output(1)).^2 (output_Y - center_output(2)).^2) radius_verify^2; phase_on_circle phase_output(circle_mask); % 由于相位是包裹的需要先解包裹再计算差分 phase_unwrapped unwrap(phase_on_circle); phase_integral sum(diff(phase_unwrapped)); estimated_charge phase_integral / (2*pi); fprintf(理论拓扑荷: %d\n, topological_charge); fprintf(通过重建相位环路积分估算的拓扑荷: %.4f\n, estimated_charge);关键点解析可视化布局使用subplot将关键步骤的结果并列显示便于对比。注意坐标轴的单位转换*1e3显示毫米和合适的色彩映射hsv用于相位hot或gray用于光强。viscircles这个函数需要Image Processing Toolbox可以方便地在频谱图上画出我们应用的滤波窗口直观展示滤波位置是否正确对准了1级衍射斑。拓扑荷验证这是检验模拟成功与否的定量步骤。我们从重建的phase_output中提取一个环形区域的相位值。由于相位是包裹的被限制在(-π, π]直接差分求和会因2π跳变而错误。unwrap函数可以沿着数据序列消除这些跳变得到连续的相位变化。然后计算相邻像素的相位差之和再除以2π就得到了环路积分的缠绕数即估算的拓扑荷。如果模拟正确这个值应该非常接近你设置的topological_charge整数。输出结果中的微小偏差如2.9998来自于数值计算误差和滤波引入的微小畸变这是可以接受的。4. 从模拟到实践关键参数影响与调试心法程序跑起来只是第一步理解每个参数如何影响结果并能在结果不理想时进行调试才是真正掌握这项技能的关键。下面分享一些我调试这类程序时积累的心得。4.1 采样与窗口避免“失真”的源头模拟的准确性首先取决于采样是否满足奈奎斯特采样定理。对于我们的螺旋相位exp(i*l*θ)其空间频率与拓扑荷l和半径r有关。在靠近中心的位置方位角θ变化剧烈空间频率最高。经验法则为了充分采样一个相位周期2π内至少要有2个像素点。螺旋相位的空间频率f_θ ≈ l / (2πr)。在最小的有效半径r_min例如光束的束腰半径或你关心的最小特征尺度处要求dx ≤ (π * r_min) / |l|。例如l5,r_min0.5mm则dx ≤ 0.314mm。我们之前设置的L10mm,N1024dx≈0.0098mm远小于要求因此采样充足。但如果l很大或L很小导致dx大就可能出现采样不足导致螺旋条纹出现锯齿状混叠失真。在频谱图上你会看到高频分量折叠回低频区域。窗口尺寸L的权衡L太小高斯光束的边缘被生硬截断会在频谱中引入高频噪声。L太大在固定像素数N下dx变大可能导致采样不足。一个实用的方法是让高斯光束的振幅在模拟窗口边缘 (L/2) 处衰减到足够小例如exp(-(L/2)^2 / waist^2) 1e-3。根据这个条件可以反推所需的L。4.2 全息图编码消除“缝线”与优化衍射效率在3.3节生成全息图时我提到了相位跳变的问题。这里再深入一下。如果你直接用mod(phase_vortex phase_grating, 2*pi)而phase_vortex是atan2输出的(-π, π]值你会发现在θ ±π这条“缝线”上phase_vortex从-π跳变到π有一个2π的跃变。虽然mod运算以2π为模但-π φ_grating和π φ_grating取模2π后的结果在大多数情况下并不相等这就导致了一条灰度不连续的线称为“2π跳变线”。加载到SLM上这条线相当于一个额外的相位阶跃会散射光降低衍射效率并引入噪声。解决方案就是我代码中采用的phase_vortex_0_to_2pi mod(phase_vortex pi, 2*pi)。这相当于将atan2的(-π, π]范围整体平移π变成(0, 2π]再取模2π将其映射到[0, 2π)。这样原来的-π点变成了0π点变成了2π取模后也是0在θ ±π处相位值都是0实现了平滑衔接消除了跳变线。衍射效率纯相位全息图的衍射效率理论上限很高。但在模拟中衍射到1级的光强占比受到光栅调制度的影响。在我们的简单叠加模型中φ_hologram mod(φ_vortex 2π*fx*X, 2π)。这实际上是一个相位型闪耀光栅其衍射效率取决于相位调制的锯齿形状的完美程度。我们的mod运算产生了完美的锯齿波因此理论上可以将大部分光能衍射到1级。在实际SLM中由于像素离散化和相位调制深度限制效率会降低。4.3 滤波与重建如何获得“干净”的涡旋光斑4f系统模拟的核心是滤波。以下几个因素会直接影响重建质量滤波窗口的位置 (shift_pixels)必须精确对准1级衍射斑的中心。你可以先运行一次不带滤波的代码显示频谱图I_freq观察1级衍射斑在x方向偏离中心多少像素然后用这个观测值来设置shift_pixels而不是仅仅依靠理论计算。因为理论计算grating_frequency * L可能因为数值离散化而不完全精确。滤波窗口的大小 (filter_radius)窗口大小是一个权衡。窗口太小相当于只允许低频通过重建的光斑会丢失细节中心暗斑可能变大甚至消失环状结构变得模糊。窗口太大会包含来自其他衍射级如0级-1级的噪声或者包含频谱的高频噪声导致重建光斑背景不干净。通常窗口半径应略大于1级衍射斑的主瓣。你可以通过逐渐增大filter_radius并观察重建光斑的变化来确定一个最佳值当光斑环状结构清晰、中心暗斑明显且继续增大半径光斑不再明显改善时即可。滤波窗口的形状我们用了圆形窗口因为它符合透镜光瞳的形状。也可以使用方形窗口但方形窗口的锐利边缘会在空域引入振铃效应吉布斯现象。一种改进是使用切趾的圆形窗口如高斯型让窗口边缘平滑过渡到零可以减少振铃但会轻微降低分辨率。调试技巧如果重建出的光强图中心不是暗斑而是一个亮斑或者相位图不是清晰的螺旋条纹首先检查滤波窗口是否真的套住了1级衍射斑。其次检查全息图生成环节的相位跳变问题是否已解决。最后可以尝试输出中间变量比如滤波后的频谱abs(fftshift(U_freq_filtered)).^2看看是否只有一个亮斑1级被保留了下来。5. 超越基础程序扩展与高级应用探索掌握了基础模拟后这个程序可以作为一个平台进行许多有趣的扩展用于研究更复杂的现象或为实验做准备。5.1 模拟部分相干涡旋光束与大气湍流理想的涡旋光束是完全相干的。但在实际中如通过大气传输的光通信链路光束会受到湍流引起的随机相位扰动导致波前畸变拓扑荷可能发生模糊甚至跳变。我们可以在生成U_vortex后叠加一个随机相位屏来模拟湍流。一种简化的方法是使用符合Kolmogorov谱的随机相位分布% 生成Kolmogorov湍流相位屏简化模型 % 注意这是一个非常简化的模型用于演示概念 r0 0.1; % 弗里德参数表征湍流强度越小湍流越强 D L; % 光束口径 % 生成随机矩阵并赋予符合特定功率谱的权重 [fx, fy] meshgrid((-N/2:N/2-1)*(1/L), (-N/2:N/2-1)*(1/L)); freq_sq fx.^2 fy.^2; freq_sq(center_freq, center_freq) 1; % 避免除零 % Kolmogorov功率谱 PSD ~ f^(-11/3) phase_spectrum (freq_sq).^(-11/6); phase_spectrum(center_freq, center_freq) 0; % 去掉零频平均相位 % 生成随机相位并乘以功率谱 random_phase randn(N) 1i*randn(N); % 复高斯随机数 turbulence_phase_screen real(ifft2(fft2(random_phase) .* fftshift(sqrt(phase_spectrum)))); % 归一化并缩放相位扰动强度 turbulence_phase_screen turbulence_phase_screen / std(turbulence_phase_screen(:)) * (D/r0)^(5/6); % 应用到光束上 U_vortex_turbulent U_vortex .* exp(1i * turbulence_phase_screen);然后用这个U_vortex_turbulent去生成全息图并进行重建。你会观察到重建出的光强环变得破碎相位条纹扭曲、断裂。通过分析大量这样的随机实现可以统计拓扑荷探测的误码率评估湍流对OAM通信的影响。5.2 生成多拓扑荷复合全息图与模式复用轨道角动量复用的核心是同时生成或探测多个不同l的涡旋光束。我们可以在一个全息图上编码多个涡旋相位每个对应不同的光栅频率从而在空间上分离它们。% 生成承载两个不同拓扑荷的全息图 l1 1; l2 3; fx1 50; % 第一个OAM信道的光栅频率 fx2 -50; % 第二个OAM信道的光栅频率负号表示衍射方向相反 % 生成两个涡旋相位 phase_vortex1 angle(exp(1i * l1 * Theta)); phase_vortex2 angle(exp(1i * l2 * Theta)); phase_vortex1_0_to_2pi mod(phase_vortex1 pi, 2*pi); phase_vortex2_0_to_2pi mod(phase_vortex2 pi, 2*pi); % 方法一叠加法简单但效率非最优 % hologram_phase_multi mod(phase_vortex1_0_to_2pi 2*pi*fx1*X ... % phase_vortex2_0_to_2pi 2*pi*fx2*X, 2*pi); % 方法二复用多个计算全息图更常用 % 分别生成两个全息图然后选择相位值例如取和或取模 % 这里采用一种简单的复用将两个相位场相加后取模 hologram_phase_multi mod(phase_vortex1_0_to_2pi phase_vortex2_0_to_2pi ... 2*pi*fx1*X, 2*pi); % 注意这里只用一个光栅频率两个模式会衍射到同一方向但重叠。 % 更复杂的复用需要设计更精巧的编码算法如GS算法优化多通道全息图。在频谱面上你会看到两套衍射斑分别位于(fx1, 0)和(fx2, 0)附近。通过使用两个独立的滤波窗口可以分别重建出l1和l3的涡旋光束。这直接模拟了OAM复用的发射端。5.3 拓扑荷检测与模态分解算法模拟在接收端我们需要判断入射的涡旋光束携带的拓扑荷是多少。一种经典的方法是使用一个与之共轭的螺旋相位板进行匹配滤波。% 假设接收到一个未知拓扑荷 l_unknown 的涡旋光束 U_unknown % 我们使用一组已知拓扑荷 l_probe 的共轭相位板进行探测 l_probe_array -5:5; % 探测的拓扑荷范围 correlation_intensity zeros(size(l_probe_array)); for idx 1:length(l_probe_array) l_probe l_probe_array(idx); % 生成共轭相位板 exp(-i * l_probe * theta) probe_phase angle(exp(-1i * l_probe * Theta)); % 注意负号 probe_field exp(1i * probe_phase); % 假设振幅均匀 % 计算与待测场的重叠积分相关 overlap sum(sum(U_unknown .* conj(probe_field))); correlation_intensity(idx) abs(overlap).^2; end % 找出相关强度最大的探测拓扑荷 [~, max_idx] max(correlation_intensity); detected_charge l_probe_array(max_idx); fprintf(探测到的拓扑荷最可能为: %d\n, detected_charge);这个简单的匹配滤波模拟表明当l_probe等于l_unknown时exp(i*l_unknown*θ) * exp(-i*l_probe*θ) 1全场相位平坦重叠积分即傅里叶透镜后的焦点光强最大。通过扫描l_probe找到光强最大的位置即可判定l_unknown。在实际实验中这可以通过一个可编程的SLM动态加载不同l_probe的相位图来实现。通过将这些扩展模块集成到你的基础程序中你可以构建一个功能强大的涡旋光束仿真研究平台从基础生成到复杂环境模拟、从单通道到多路复用、从发射到探测全方位地深化对轨道角动量这一前沿领域的理解。本文还有配套的精品资源点击获取