4. 部分相干光场模拟：三种主流生成方法
好，咱们进入实战环节。前面几章讲了互强度函数的理论，那都是基础。但真正干活的时候，你面对的是“怎么把部分相干光造出来”这个问题。
我刚开始做相干成像仿真时，总习惯用完全相干光。结果呢？模拟出来的干涉条纹跟实验数据对不上。后来才意识到——真实光源哪有那么理想？激光器都有一定的谱宽和空间相干性。所以，学会生成部分相干光场，是连接理论与实验的必经之路。
目前业界主流方法有三种：随机相位屏法、模态分解法，以及直接数值积分法。这一章咱们重点讲前两个，第三个我会在实战篇里补充。
核心观点：部分相干光模拟的本质，就是用有限个“不完全相干”的子光场，去逼近真实的互强度函数。你不需要完美复现，只需要在统计意义上一致。
4.1 随机相位屏法：最直观的工程方法
这个方法的思想很简单：给一个完全相干的光场，叠加上随机相位扰动。你想想看，部分相干光不就是“相干性被破坏”的光吗？那我们就主动去破坏它。
具体怎么做呢？我一般分三步走：
- 生成初始相干光场——比如一个高斯光束的复振幅分布。
- 构造随机相位屏——用特定统计分布的随机数，生成一个二维相位矩阵。
- 叠加并统计平均——将相位屏乘到光场上，重复多次取平均。
这里有个关键点：随机相位屏的统计特性决定了相干度。如果你用均匀分布的随机相位，得到的是完全非相干光。如果你用高斯相关随机相位，就能控制相干长度。
我的经验：在实际项目中，我常用“高斯型相关函数”来生成相位屏。具体来说，相位屏的自相关函数设为 exp(-r²/ξ²)，其中 ξ 就是相干长度。这样生成的相位屏，物理意义明确，而且跟实验对得上。
来看一段核心代码。这是生成高斯相关随机相位屏的 Python 实现：
import numpy as np
def generate_phase_screen(N, dx, coherence_length):
"""
生成高斯相关随机相位屏
N: 网格点数
dx: 网格间距
coherence_length: 相干长度
"""
# 1. 生成白噪声
noise = np.random.randn(N, N)
# 2. 构造高斯滤波器
kx = np.fft.fftfreq(N, dx) * 2 * np.pi
ky = np.fft.fftfreq(N, dx) * 2 * np.pi
kx, ky = np.meshgrid(kx, ky)
# 高斯型功率谱
sigma_k = 1.0 / coherence_length
filter_h = np.exp(-(kx**2 + ky**2) / (2 * sigma_k**2))
# 3. 频域滤波
noise_fft = np.fft.fft2(noise)
phase_fft = noise_fft * filter_h
# 4. 反变换得到相位屏
phase_screen = np.real(np.fft.ifft2(phase_fft))
# 归一化到 [-π, π]
phase_screen = phase_screen / np.std(phase_screen) * np.pi
return phase_screen这段代码的核心逻辑是：白噪声 → 频域滤波 → 反变换。滤波器的宽度决定了相干长度。我习惯把相干长度设成光源尺寸的 1/10 到 1/5，这样模拟出来的部分相干光效果最自然。
注意：随机相位屏法有一个坑——单次实现不能代表部分相干光。你必须做 多次独立实现并取统计平均。我一般至少跑 100 次，如果要求高，500 次起步。次数太少，相干度曲线会抖动得很厉害。
4.2 模态分解法：更优雅的数学方法
如果说随机相位屏法是“暴力破解”，那模态分解法就是“精准打击”。它的理论基础是：任何部分相干光场都可以分解成一组正交相干模式的加权和。
这个思想来源于 Wolf 的相干模式分解理论。说白了，就是把互强度函数 J(r₁, r₂) 写成：
J(r₁, r₂) = Σ λₙ φₙ*(r₁) φₙ(r₂)
其中 λₙ 是特征值，φₙ 是正交模式。每个模式都是完全相干的，但模式之间互不相干。
这样做的好处是什么？你只需要生成有限个相干模式，叠加起来就是部分相干光。而且模式数量越少，相干度越高；模式数量越多，相干度越低。
| 对比项 | 随机相位屏法 | 模态分解法 |
|---|---|---|
| 计算效率 | 需要多次独立运行 | 一次计算，模式叠加 |
| 物理直观性 | 直观，容易理解 | 数学抽象，需理解特征分解 |
| 适用光源 | 任意光源 | 低相干度光源效果更好 |
| 精度控制 | 靠样本数量 | 靠模式截断阶数 |
| 我的推荐 | 工程仿真首选 | 理论分析首选 |
我举个例子。对于高斯-谢尔模型光束，它的相干模式有解析解——就是厄米-高斯函数。你只需要取前 M 个模式，加权求和就能得到任意相干度的光场。
def gsm_modes(x, y, w0, sigma_g, M):
"""
生成高斯-谢尔模型光束的前 M 个相干模式
x, y: 坐标网格
w0: 束腰半径
sigma_g: 相干长度
M: 模式数量
"""
# 计算特征值
beta = 1 / (4 * w0**2) + 1 / (2 * sigma_g**2)
gamma = np.sqrt(1 - (1/(2*sigma_g**2*beta))**2)
modes = []
eigenvalues = []
for m in range(M):
for n in range(M):
# 厄米-高斯模式
Hm = np.polynomial.hermite.hermval(
np.sqrt(2*beta) * x, [0]*m + [1])
Hn = np.polynomial.hermite.hermval(
np.sqrt(2*beta) * y, [0]*n + [1])
mode = Hm * Hn * np.exp(-beta * (x**2 + y**2))
mode = mode / np.sqrt(np.sum(np.abs(mode)**2))
# 特征值
lam = (1 - gamma) * gamma**(m+n)
modes.append(mode)
eigenvalues.append(lam)
return modes, eigenvalues这段代码里，gamma 控制了相干度。gamma 越接近 1，相干度越高，需要的模式越少。gamma 接近 0，相干度低，需要很多模式。
实战建议：我一般用模态分解法做理论验证，用随机相位屏法做工程仿真。如果时间紧，直接上随机相位屏法，100 次平均基本够用。如果要做论文级别的精度分析，模态分解法更可控。
4.3 两种方法的对比与选择
说了这么多，到底该用哪个？我总结几条经验：
- 场景一：模拟大气湍流——用随机相位屏法。因为湍流本身就是随机过程，物理上匹配。
- 场景二：分析部分相干成像系统——用模态分解法。因为你可以单独分析每个模式的传递特性。
- 场景三：快速原型验证——用随机相位屏法。代码简单，调参方便。
- 场景四：高精度理论计算——用模态分解法。误差可控，可解析推导。
我曾经在一个激光雷达项目中，需要模拟部分相干光在雾中的传播。一开始我用模态分解法，结果模式数量太多，算到崩溃。后来换成随机相位屏法，把雾的散射等效成随机相位扰动，效率提升了 10 倍。嗯，这就是工程思维——选对方法比算对结果更重要。
一个小技巧：如果你不确定选哪种，可以先用随机相位屏法跑一个快速版本，看看结果趋势。如果趋势合理，再决定要不要用模态分解法做精细化。我管这叫“先粗后精”策略，能省不少时间。
最后提醒一句：无论用哪种方法，都要验证生成的互强度函数是否与理论一致。我习惯的做法是：生成光场后，直接计算它的互强度函数，跟目标函数对比。如果偏差超过 5%，就调整参数或增加样本数。这个验证步骤，能帮你避免很多隐藏的 bug。