企业官网建设流程全解析

简介：这份仿真脚本围绕高斯光束在大气湍流中的传输展开，面向光学工程、大气光学与自由空间光通信方向的学生和研究人员，用于分析光强闪烁、相位畸变等核心问题，也可作为课程设计或科研预研的起点。包内共1个m脚本，大小仅1KB，属于轻量级可运行源码，包含高斯光束初始参数设置、大气湍流参数定义、湍流屏生成及光强/相位统计计算等模块，适合快速复现Rytov近似与Kolmogorov湍流模型下的典型传播过程。资源已有1394人浏览学习。通过运行该脚本，可直观比较有无湍流时光强分布的变化，观察相位结构函数，理解Cn²谱和湍流强度对光束扩展、光强闪烁的影响；其简洁结构也便于二次开发，为后续扩展涡旋光束、多屏传输或自适应光学补偿等研究提供基础参考。

1. 高斯光束在大气湍流里的传输仿真：先搞清楚要解决什么问题

做激光通信、激光雷达或者自由空间光路的人，多半会遇到同一个场景：发射端出来的高斯光束，在实验室里光斑圆润、能量分布干净，但一到室外走个一两公里，光斑就开始抖、开始闪，接收面上的中心光强忽高忽低，甚至出现破碎的环状结构。这不是光学平台没调平，也不是激光器不稳定，而是大气湍流在起作用——温度起伏导致空气折射率随机波动，光束相位被一遍遍“揉搓”。与其花大价钱做外场实验反复试错，不如先把这套随机介质建模成可复现的仿真流程：把连续的大气信道切成多层相位屏，让高斯光束一层层穿过去，在接收面统计光强分布、光斑漂移和闪烁指数。

这份高斯光束大气湍流传输仿真资源，核心就是帮你搭起这样一套数值实验环境，覆盖光强分布、湍流相位屏生成、分步传播和参数调优。适合三类人：刚接触大气光学仿真、想把理论公式落成可运行代码的研究生；做激光链路预算需要预估闪烁衰减的工程师；以及想验证自适应光学方案但不想上来就砸硬件成本的项目组。接下来从数学模型讲起，一路走到代码实现和参数调优。

2. 湍流信道的数学模型：从折射率起伏到相位屏逼近

2.1 高斯光束的复振幅：束腰、瑞利长度与波前曲率

任何传输仿真第一步都是把光源表达清楚。基模高斯光束在自由空间中的复振幅分布写成：

E(r, z) = E0 * (w0/w(z)) * exp(-r²/w(z)²) * exp(-i*k*z - i*k*r²/(2*R(z)) + i*zeta(z))

其中 w0 是束腰半径，w(z) = w0sqrt(1+(z/zR)²) 是传播到 z 处的光斑半径，zR = piw0²/lambda 是瑞利长度，R(z) = z*(1+(zR/z)²) 是等相面曲率半径。这几个参数决定了光束在自由空间里的扩散速度：束腰越小，瑞利长度越短，光束发散越快；波长越长，衍射越明显。

实际仿真中不直接对 E(r,z) 表达式逐点采样，因为湍流扰动会让相位项变得极其复杂。更常见的做法是只初始化发射面的复振幅场 E(r,0)，后续传播全部交给衍射数值算法。初始化时注意：如果你要模拟的是准直光束，那么发射面波前曲率 R(0) 为无穷大，相位项为零；如果要模拟聚焦光束，则要在相位里预置一个球面波因子。

2.2 大气湍流的统计描述：Kolmogorov 谱与 von Kármán 谱

大气湍流对光束的影响体现在折射率的随机起伏上。Kolmogorov 湍流理论给出折射率结构函数服从 2/3 次方律，对应的功率谱密度为：

Phi_n(kappa) = 0.033 * Cn² * kappa^(-11/3)

其中 Cn² 是折射率结构常数，单位 m^(-2/3)，代表湍流强度；kappa 是空间波数。这个形式在数学上简洁，但存在两个奇点：kappa 趋于零时谱密度发散，kappa 趋于无穷时总能量发散。所以工程仿真里我更习惯用 von Kármán 谱，它引入了内尺度 l0 和外尺度 L0：

Phi_n(kappa) = 0.033 * Cn² * exp(-kappa²/kappa_l²) / (kappa² + kappa_0²)^(11/6)

kappa_l = 5.92/l0 对应内尺度截止，kappa_0 = 2*pi/L0 对应外尺度低频饱和。外尺度决定了湍流涡旋的最大尺寸，直接影响光束漂移的幅度；内尺度影响高频相位起伏，对闪烁指数的高频成分有贡献。做仿真时建议直接用 von Kármán 谱，逼真程度和数值稳定性都比纯 Kolmogorov 谱好。

除了谱模型，还有两个工程上常用的导出量。第一个是 Fried 相干长度 r0，它表示大气湍流相干孔径的大小：

r0 = 0.185 * (lambda² / (Cn² * L))^(3/5)

r0 越小湍流越强。第二个是 Rytov 方差 sigma_R² = 1.23 * Cn² * k^(7/6) * L^(11/6)，它用来判断闪烁强弱：小于 0.3 是弱起伏，0.3 到 1 是中等起伏，大于 1 进入强起伏饱和区。这两个量在后续参数设计里会反复用到，是判断仿真结果合理性的标尺。

2.3 为什么用相位屏：把连续随机介质切成离散薄层

真实大气沿传输路径是连续随机介质，光束每走一毫米都会受到相位扰动，但数值仿真无法对连续介质逐点采样。相位屏法的核心思想是：把整段传输路径 L 分成 Nz 段，每段长度 dz = L/Nz。在每一段的末端放一个薄相位屏，这个屏只改变光束的相位、不改变振幅，屏与屏之间按自由空间衍射传播。数学上是对随机介质做了一次“分段冻结”近似。

这个近似的物理依据是：湍流涡旋的特征时间远大于光波穿过单个分段时间，所以可以把每一段内的湍流看成静止的。只要 dz 选取得当，相位屏之间的互相关性可以忽略，多层屏就能很好地逼近真实连续扰动。相位屏的相位分布由湍流功率谱决定，生成方法下一章展开。需要注意的是，相位屏的统计特性必须与它代表的这一段路径长度匹配，每层屏的相位方差和 Cn²*dz 成正比。

3. 仿真主循环实现：相位屏生成与分步传播

3.1 整体仿真框架：初始化、循环、统计三段式

把仿真拆成三个模块，代码结构会清晰很多：初始化模块负责建立网格、生成光源复振幅、预计算传播核；主循环负责 Nz 次交替执行自由传播和相位屏扰动；统计模块负责记录每个接收面的光强分布，供后续计算闪烁指数等指标。下面这份 Python 代码用 numpy 实现，逻辑同样可以平移到 MATLAB。

3.2 相位屏生成：FFT 方法加次谐波补偿

相位屏生成是整套仿真里最容易出问题的地方。直接做法是对功率谱做逆傅里叶变换：生成频域复随机场，乘以谱密度的平方根，再 IFFT 回空间域。但 FFT 网格的采样间隔有限，低频分量（大尺度涡旋）采样不足，导致生成的光束漂移偏小，光斑破碎程度不够。

import numpy as np def generate_phase_screen(N, dx, L0, l0, Cn2, dz, seed=None): """ 生成 von Kármán 谱相位屏 N: 网格边长（像素数） dx: 空间采样间隔（m） L0: 外尺度（m） l0: 内尺度（m） Cn2: 折射率结构常数（m^-2/3） dz: 该相位屏代表的路径长度（m） """ if seed is not None: np.random.seed(seed) # 频域网格 kx = np.fft.fftfreq(N, d=dx) * 2 * np.pi ky = np.fft.fftfreq(N, d=dx) * 2 * np.pi KX, KY = np.meshgrid(kx, ky) K2 = KX**2 + KY**2 # von Kármán 谱 kappa_l = 5.92 / l0 kappa_0 = 2 * np.pi / L0 Phi = 0.033 * Cn2 * np.exp(-K2 / kappa_l**2) / (K2 + kappa_0**2)**(11/6) # 频域随机场，保留厄米特对称性保证实值输出 delta = np.sqrt(2 * np.pi / (N * dx)**2) * np.sqrt(Phi) cn = np.random.randn(N, N) + 1j * np.random.randn(N, N) phase_screen_f = cn * delta # 保证 IFFT 结果是实数 phase_screen_f = np.fft.fftshift(phase_screen_f) phase_screen = np.fft.ifft2(phase_screen_f).real * N**2 # 次谐波低频补偿 sub_harmonics = add_subharmonics(N, dx, L0, l0, Cn2, dz) return phase_screen + sub_harmonics

这段代码的关键逻辑在两步：第一步用 FFT 方法生成主相位屏，随机场乘以谱密度的平方根再变换回空间域；第二步调用add_subharmonics补低频。delta里的sqrt(2*pi/(N*dx)^2)是功率谱到频域振幅的归一化系数，不同教材写法略有差异，但本质是把连续谱密度离散化到频域网格上。fftshift是为了把零频移到中心，保证 IFFT 结果的相位分布正确。

add_subharmonics的实现不展开写完整代码，核心思路是：把频域低频区按 3x3、9x9 等更细的网格重新采样，用同样谱公式生成小尺寸随机场，再放大插值叠加回主相位屏。这一步能显著增强光束漂移的真实感，强烈建议保留。

3.3 分步传播：角谱法主循环

光束在两层相位屏之间走自由空间衍射。工程上有两种选择：菲涅尔衍射积分和角谱法。角谱法的好处是单步传播没有近场限制，只要网格采样满足条件，任意距离都能算。

def angular_spectrum_propagate(E, lamb, dx, dz): """ 角谱法自由空间传播 E: 输入复振幅场（N x N） lamb: 波长（m） dx: 空间采样间隔（m） dz: 传播距离（m） """ N = E.shape[0] k = 2 * np.pi / lamb # 频域坐标 fx = np.fft.fftfreq(N, d=dx) fy = np.fft.fftfreq(N, d=dx) FX, FY = np.meshgrid(fx, fy) F2 = FX**2 + FY**2 # 传递函数，注意 evanescent 波截止 H = np.exp(-1j * np.pi * lamb * dz * F2) H[F2 > (1/lamb)**2] = 0 # 正变换-相乘-逆变换 Ef = np.fft.fft2(E) Ef = Ef * H return np.fft.ifft2(Ef)

主循环就是一个 Nz 次的循环：每走一段 dz，先做自由传播，再乘一个相位屏。注意相位屏是以复指数形式施加的：E = E * np.exp(1j * phase_screen)。

def turb_sim_main(N=512, dx=0.002, lamb=1.55e-6, L=1000, Nz=10, Cn2=1e-14, L0=10.0, l0=0.01, seed=42): """ 高斯光束经大气湍流传输的主循环 """ # 初始化发射面高斯光束，束腰 w0=0.02m x = np.arange(-N/2, N/2) * dx X, Y = np.meshgrid(x, x) r2 = X**2 + Y**2 w0 = 0.02 E0 = np.exp(-r2 / w0**2) E = E0.copy() z = 0.0 dz = L / Nz rng_seed = seed for i in range(Nz): E = angular_spectrum_propagate(E, lamb, dx, dz) # 生成该段对应的相位屏 ps = generate_phase_screen(N, dx, L0, l0, Cn2, dz, seed=rng_seed+i) E = E * np.exp(1j * ps) z += dz # 可在此处记录光强 I = |E|^2 I = np.abs(E)**2 return I

主循环里我习惯对每一段相位屏用不同的随机种子，避免屏与屏之间出现统计相关。generate_phase_screen里传入的dz是这段路径代表的湍流强度累计值，层数越多每层越薄，相位屏方差越小，但总扰动守恒。完成循环后，I就是接收面的光强分布，可以直接画等高线图，也可以计算质心、光斑半径和闪烁指数。

有一点务必注意：相位屏施加的位置在每段传播之后，而不是之前。先传播完一段自由空间，再打上这一段积累的相位扰动，这个顺序对应物理过程的因果性。反过来会让能量分布失真。

4. 网格与湍流参数怎么设：采样条件、Rytov 方差与分段数

4.1 采样间隔和网格大小：混叠是第一道门槛

网格参数选错，仿真结果就毫无意义。空间采样间隔 dx 必须足够小，能分辨光斑的精细结构；但如果 dx 过小，网格数 N 固定时计算区域太小，光束发散后撞到边界产生伪影。一个常用的约束是相邻像素间的相位差要小于 pi，否则角谱法会混叠。经验公式是 dx 要满足：

dx > lamb * L / N

这个条件保证角谱传递函数的最大相位变化不超过 pi。比如波长 1.55 um，距离 1 km，N=512，那么 dx 要大于 3 mm。但 3 mm 采样间隔对束腰 2 cm 的光束来说又略大——光斑中心只有大概 6~7 个采样点，光强分布细节不够。这时就提高 N 到 1024，把 dx 降到 1.5 mm，同时满足防混叠和分辨率。

网格大小的第二层考量是边界。高斯光束自由传播 1 km 后，光斑半径大概从 2 cm 扩到几十 cm，加上湍流引起的扩散和漂移，计算区域至少要是最终光斑半径的 3~4 倍。否则光瓣撞到边界再绕回来，接收面光强会出现干涉条纹，这就是常说的“边界翻车”。

4.2 湍流强度参数：Cn²、r0 与 Rytov 方差的联动

Cn² 的取值随高度和气象条件变化极大。近地面白天强烈日照下可达 1e-13 甚至更高，夜晚或阴天降到 1e-15 量级，高空平流层可以低到 1e-17。做仿真前先想清楚你模拟的场景，再定 Cn²。

表：不同条件下推荐参数组合

Rytov 方差大于 1 后进入强起伏区，相位屏法依然可用，但闪烁指数会出现饱和效应——不再随距离线性增长。这时如果仿真结果里闪烁指数远高于 1，先别怀疑代码，再看看理论上的饱和规律。

4.3 分段数 Nz：多少层相位屏才算够

分段数直接影响仿真精度，也影响计算量。太少的相位屏无法体现湍流的纵向随机性，光束只被几层屏“猛拍”，光强分布出现明显的人工条纹；太多层相位屏生成和传播计算量翻倍。

我常用的判定方法是看单层相位屏代表的路径上，Rytov 方差增量 d_sigma_R² = 1.23 * Cn² * k^(7/6) * dz^(11/6) 是否远小于 1。一般要求 d_sigma_R² < 0.05，这样每层屏只产生小扰动，多层叠加逼近连续介质。对应到 1 km 路径、1e-14 的 Cn²，Nz = 10 就够；但如果 Cn² = 1e-13，同样距离建议 Nz 至少取 20。另一个经验是 Nz 不要小于 5，少于 5 层的仿真结果基本不具备统计意义。

分段方式可以等距划分，也常见按湍流强度非均匀划分：近地面段湍流强、网格更密，高层段稀疏。但初始调试时先用等距，跑通了再优化。

5. 常见问题与避坑：5条血泪经验

5.1 光斑漂移量偏小，光强分布过于“干净”

现象：加了相位屏后，接收面光斑仍然圆润，只是轻微变形，没有观察到明显的位置抖动。

原因：相位屏的低频成分不足。FFT 方法生成相位屏时，频域网格最小波数受限于 2pi/(Ndx)，对应最大涡旋尺寸只有 N*dx/pixel 量级。外尺度 L0 = 10 m 的湍流涡旋在 1 m 量级的计算区域内根本表达不出来，低频能量缺失直接导致光束大角度漂移没有被模拟到。

解决：加入次谐波补偿，把低频区细分采样。我在generate_phase_screen里保留了add_subharmonics这个钩子，就是为这个准备的。次谐波阶数取 3~4 层即可，再多计算量上去了收益很小。

5.2 长距离传播后能量不守恒，光强整体衰减

现象：同样的初始能量，自由传播 1 km 后总能量明显下降，甚至出现环状伪影。

原因：角谱法的传递函数 H = exp(-ipilambdzF2) 对高频是纯相位调制，理论上不会损失能量。但如果 F2 超过 (1/lamb)^2，传递函数变成指数衰减的倏逝波项。代码里虽然把 H 置零了，但置零本身也截断了能量。另一个更常见的坑是采样不满足 dx > lamb*L/N，导致高频混叠折返。

解决：先跑一遍零湍流基准测试——把 Cn² 设成 0，看自由传播的光斑半径是否吻合 w(z) 理论值，总能量是否守恒到浮点误差范围。如果能量掉了，优先检查 dx 和 N 是否满足防混叠条件，其次把传播分成更多步子。

5.3 闪烁指数不随距离增长，甚至反向下降

现象：增大传播距离 L，接收面光强起伏反而变小，与 Rytov 理论的增长趋势矛盾。

原因：相位屏数量不足。当 Nz 固定时，增大 L 意味着每层屏代表更长的 dz，单层扰动过大，光线在每一层都被剧烈偏折，等效于人为引入了强散射，出现了与真实物理不符的“伪饱和”。

解决：按 4.3 节的判据倒推 Nz。L 从 1 km 增大到 2 km，其他参数不变，Nz 也要跟着翻倍。我在调试时会在主循环里打印单层 Rytov 方差增量，一旦超过 0.05 就自动报警，提醒自己加密分段。

5.4 两次运行结果差异巨大，无法复现

现象：固定所有参数，只改随机种子，接收面光斑形状和峰值位置完全不同，有时连光强量级都有差异。

原因：这不是 bug，而是湍流的随机本质。单次仿真只是蒙特卡洛的一次抽样，只有多次独立运行做系综平均才具备统计意义。很多人第一次跑仿真后看到结果不稳定，就以为代码有问题，其实结论恰恰相反——结果对种子敏感才说明相位屏的随机性生效了。

解决：在统计阶段跑 50~100 次独立种子，对每次的闪烁指数、质心位置做平均。频率更高的是“滚动相位屏”技巧：用一个大尺寸相位屏，每次按窗口滑动取子区域，等效一次运行多条独立路径，省掉重复初始化光源的时间。

5.5 相位屏出现明显的网格纹理或周期伪影

现象：接收面光强分布里有规律性的横竖条纹或重复图案，看起来像衍射光栅叠加在光斑上。

原因：FFT 相位屏隐含周期性延拓假设，大尺度涡旋跨越计算边界时会从另一侧折返回来，形成周期伪影。尤其当外尺度 L0 与计算区域尺寸可相比时，这种伪影非常明显。

解决：三个办法叠加使用。第一，计算区域大于外尺度 L0 的 2~3 倍；第二，在相位屏四周乘以一个余弦窗函数，把边界区域的相位扰动压到零；第三，接受边界区域的失真，统计指标只在中心区域计算。第三个办法在工程上最常用，因为边界伪影对中心光斑的影响有限，不值得为它翻倍计算量。

6. 从单帧到统计量：闪烁指数与光束质量评价

单帧光强分布看得再多，也只是“一张照片”。真正对工程有意义的，是透过几十上百次独立传输提取出的统计量。我最先算的是闪烁指数——接收光强的归一化方差，定义是 sigma_I² = <I²>/² - 1。弱湍流下它应接近 Rytov 方差，强湍流下出现饱和，不会无限增长。这个指标直接对应通信链路的信噪比预算，是判断仿真结果物理合理性的第一道关卡。

其次是光束质心漂移。每次传输后算一次光强质心位置 (xc, yc)，统计它在水平方向的方差。这个方差与外尺度 L0 密切相关，L0 设得太小漂移会偏小。我一般会把仿真算出的质心漂移均方根值，和经典公式 0.97 * (lambda*L/r0)^(1/3) * L 做个对比，偏差在 20% 以内就说明低频补偿做对了。

光束质量评价建议用 M² 因子。接收面光强场做二阶矩计算，得到光斑半径 w(z)，再对多重 z 位置拟合双曲线，提取 M² = pi * w0_m * theta / lambda，其中 w0_m 是拟合出的束腰，theta 是远场发散角。湍流会让 M² 从 1 增大到 2~5 不等，这个数字比单看光斑半径直观得多。拟合时注意不要用强湍流状态下的几个点硬套，多取 5~8 个距离位置，拟合稳定度更好。

最后建议你做一个“零湍流自检”脚本：把 Cn² 设为 0，跑完整流程，确认光斑随距离的扩展符合自由空间高斯光束公式，闪烁指数为 0。从那以后我每次跑新参数，都强制先过一遍这个自检，再开湍流——这一招帮我拦下了至少三回网格参数配错的翻车。希望帮到你。

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