延迟自外差法测量窄线宽激光器线宽的仿真与拟合

简介：本资源是一套面向本科及硕士阶段科研学习的激光线宽仿真拟合工具，基于Matlab实现DSH（Delay Self-Heterodyne）干涉信号建模与线宽参数反演，适用于光学测量、激光器性能评估及信号处理等研究场景。压缩包共20个文件，含11个核心Matlab函数（.m）、2个PDF技术文档（含仿真原理与使用说明）、2个PNG示意图（含工具界面与流程图）、3个文本说明文件（含Readme与License）以及2个CSV数据样例，整体体积仅1.83MB，结构紧凑、即装即用。已有61人下载学习，配套文档清晰标注运行环境（兼容Matlab 2014a/2019a）并附带完整仿真结果，降低初学者调试门槛；用户可直接调用工具箱模块开展DSH信号生成、噪声建模、频谱估计与非线性拟合全流程分析，亦可结合general_routines与conventional_fitting子模块拓展对比实验。 调了一台标称线宽几十kHz的窄线宽激光器，我第一反应是拿光谱仪直接看线宽。结果屏幕上那条谱线比我的拇指还粗——光谱仪分辨率到顶也就0.01nm，换算到1550nm波段大约是1.25GHz，和kHz级别的真实线宽差了四五个数量级，根本没法看。后来改用光电探测器加频谱仪，又发现直流附近全是弛豫振荡峰和低频噪声，想要的那一点线宽信息早就被淹没了。最后只能老老实实搭延迟自外差（DSH，Delayed Self-Heterodyne）测量链路。

但问题并没有在搭好链路之后结束。链路里每一根光纤、每一个转接头、AOM的驱动功率、探测器的带宽，甚至偏振态的漂移都在影响最终结果。数据出来了，我第一反应不是"测到了多少"，而是"这个数到底准不准"。所以我把整个DSH测量链路整体搬进了仿真：用已知参数的激光场生成拍频信号，再走一遍与实测完全相同的拟合算法，反过来检验程序的可信度。这篇文章就是这次仿真加拟合程序的完整复盘。

这里先说明一下，标题里的DSH指的是延迟自外差法，不是某些软件工具里的插件市场那种DSH。内容围绕三个问题展开：为什么要用DSH、怎么在仿真环境里生成拍频数据、以及拟合程序如何把线宽从频谱里反演出来。文中所有代码基于Python，参数配置和结果都有具体数值，正在做窄线宽激光器测量、或者想用仿真验证线宽反演算法的朋友可以直接参考。

1. 从"光谱仪读不出线宽"说起：DSH方案的提出与仿真必要性

1.1 频谱仪分辨率就是第一道坎

窄线宽激光器的线宽定义，通常指光场功率谱的半高全宽（FWHM）。kHz到MHz量级的线宽，如果直接拿光谱仪看，基本是死路一条。以1550nm波段为例，一台分辨率0.01nm的光谱仪，对应的频率分辨率大约是1.25GHz。激光器线宽是100kHz，谱线在屏幕上就是一整根像素条的宽度，连个轮廓都看不出来。就算分辨率做到0.001nm（很多实验室的顶级OSA也就这水平），也就对应125MHz，和100kHz之间依然隔着三个数量级。

那能不能用超高精度的法布里-珀罗干涉仪？可以，但扫描式FP腔对扫描电压、温漂和震动极敏感，而且自由光谱范围有限，标定也麻烦。更本质的问题在于，线宽这个参数是统计量，需要对光场的相位涨落做长时间观测，FP腔一次性只能看一小段光谱窗口，效率低，数据处理也繁琐。

DSH的思路是绕开"直接测量光频"这个难题，把光频域的测量转换到射频域。射频域的频谱分析仪分辨率做到1Hz都不难，关键是先让激光"自己和自己拍频"，把相位涨落信息降到一个可以直接分析的频率范围。

1.2 DSH的基本原理：让激光自己和自己"打架"

DSH的光路听起来并不复杂。

• 激光器输出被一个分束器分成两路。
• 一路经过一段延时光纤，延迟时间为τ。
• 另一路经过一个声光移频器（AOM），光频率整体平移一个量Ω，典型值是80MHz或200MHz。
• 两路光合束进入光电探测器，产生拍频信号。

假设激光电场写为：

E(t) = E0 · exp(j(2πν0t + φ(t)))

其中ν0是光频率，φ(t)是相位噪声。经过延迟和移频后，探测器上的拍频光电流（只关注射频项）正比于：

I_beat(t) ∝ cos(2πΩt + φ(t) - φ(t-τ))

这个式子里的关键，是相位差φ(t) - φ(t-τ)。激光器的相位噪声是随机的，延迟τ越长，这两个时刻的相位差涨落就越大。拍频信号的功率谱，由这个相位差的时间统计决定。当延迟时间远大于激光器的相干时间τc时，拍频频谱会呈现洛伦兹线型，其FWHM正好是激光器线宽的2倍。

"自己和自己拍"的好处在于：两路光的频率在拍频瞬间相减，光频的绝对抖动被消掉了，剩下的纯粹是相位噪声统计。这就是为什么它能把一个GHz量级的问题降到一个几百kHz、在射频域可以精细测量的问题。

1.3 为什么先做仿真再上实测

搭好DSH链路之后，我很快发现一个尴尬处：实测链路上的变量太多了。光纤长度差一点，偏振态稍微漂移，AOM驱动功率波动，探测器的带宽滚降，每一个因素都会改变拍频谱的形状。数据出来了，拟合出了线宽，但你没法判断这个结果是激光器真实线宽，还是测量链路引入的系统偏差。

仿真环境里没有这个问题。激光器的线宽、延迟时间、移频量、采样率、噪声基底全部由我设置，拟合结果的"真值"是已知的。这样做的价值有三个：

• 先验证拟合算法本身有没有算错，避免拿着错误程序去处理实测数据。
• 评估测量参数（尤其是延迟时间）是否满足要求。仿真里可以随意改τ，看看不同配置下的结果偏差。
• 提前暴露测量中容易被忽视的坑，比如频谱泄漏、噪声耦合、采样长度不够。这些在实测里会让人排查到崩溃，但在仿真里一眼就能看清。

所以这个zip里的程序，本质上是一个"测量链路仿真器"加"线宽反演算法验证器"。先让仿真数据通过，再去处理真实拍频数据，心里有底得多。

2. 仿真信号怎么生成：相位噪声、延迟拍频与频谱估计的参数逻辑

2.1 相位噪声建模：维纳过程为什么够用

生成仿真激光场，核心是生成相位噪声φ(t)。对于单模激光器，如果频率噪声是白噪声，那么相位是一个维纳过程，也就是随机游走。表达式为：

φ(t + Δt) = φ(t) + sqrt(2πΔνΔt) · n

其中Δν是激光线宽的FWHM，n是标准正态分布随机数。用这个模型生成的光场，其功率谱恰好是洛伦兹线型，线宽正好等于Δν。这个性质让"生成已知线宽的仿真数据"变得非常直接：我设置Δν=50kHz，生成出来的数据理论上就是50kHz，然后就可以去检验拟合程序能不能把它反演出来。

注意这里使用复数基带表示。实际光频率ν0高达193THz，但仿真完全不需要处理这么高的频率，把光载波去掉，只保留相位随时间的变化就够了。基带信号中，光的幅度不变，相位做随机游走，这就是相干性退化的全部体现。

2.2 时域信号生成与采样参数选择

采样率和信号长度的选择是整个仿真能不能跑通的关键。拍频信号的中心频率是Ω，也就是AOM移频量，所以采样率至少要满足奈奎斯特条件。如果用80MHz的AOM，采样率取fs=400MHz比较合理，给信号留出至少5倍的带宽余量。

延迟时间的实现要转成整数采样点数。假设采样率400MHz，延迟τ=20μs，对应的采样点数是：

n_tau = int(20e-6 × 400e6) = 8000

也就是说，当前时刻的相位减去8000个采样点之前的相位，就得到了延迟自外差所需的相位差。这个操作在代码里就是数组移位。

生成拍频光电流的核心代码不长。我用一个数组存相位噪声，做一次移位，再乘上AOM的余弦项：

import numpy as np fs = 400e6 # 采样率 400 MHz T = 2e-3 # 信号总时长 2 ms N = int(fs * T) # 总采样点数 t = np.arange(N) / fs delta_nu = 50e3 # 目标激光线宽 50 kHz tau = 20e-6 # 延迟时间 20 us Omega = 80e6 # AOM移频 80 MHz n_tau = int(tau * fs) rng = np.random.default_rng(42) # 相位噪声：维纳过程 dphi = np.sqrt(2 * np.pi * delta_nu / fs) * rng.standard_normal(N) phi = np.cumsum(dphi) # 延迟后的相位 phi_delayed = np.concatenate([np.zeros(n_tau), phi[:-n_tau]]) # 拍频光电流 beat = np.cos(2 * np.pi * Omega * t + phi - phi_delayed)

关键参数之间的匹配关系必须提前算清楚。

• 信号总时长T决定FFT的频率分辨率Δf=1/T。要分辨50kHz线宽，最好将分辨率控制在5kHz左右，也就是T至少200μs。取T=2ms，分辨率0.5kHz，这样拟合时线型轮廓不会被分辨率过度抹平。
• 延迟时间τ必须远大于激光器相干时间τc=1/(πΔν)。对于50kHz线宽，τc≈6.4μs（以FWHM=50kHz计算，τc = 1/(π×50e3) ≈ 6.4μs），20μs延迟大约是3倍相干时间，已经接近可接受，但还不是理想状态。后面章节会专门说这个比值对结果的影响。
• 采样率若取得太低，拍频信号会出现混叠，峰形会变成"折叠"的样子，这个在频谱图上非常难排查。

2.3 FFT频谱计算的细节：加窗、分段平均与频段截取

生成拍频信号之后，做FFT得到功率谱，这一步看似简单，但有几个细节直接决定拟合质量。

第一个细节是加窗。矩形窗的旁瓣在频域里只衰减约13dB，对洛伦兹这种重尾谱线来说，泄漏出来的尾巴会把拟合结果污染得很厉害。用汉宁窗，旁瓣衰减能到31dB左右，代价是主瓣展宽约一倍。对于50kHz线宽，使用65536点FFT时主瓣展宽约12kHz，相比线宽仍然足够小，可以接受。

第二个细节是谱平均。直接用整段数据做一次FFT，噪声方差极大，谱线边缘毛刺很重，拟合时这些毛刺会引入额外误差。我用Welch方法做分段平均，每段65536点，重叠50%，将2ms信号切成约60段。频率分辨率约6.1kHz，但谱线平滑度大幅提升，拟合稳定性明显好于单次FFT。

from scipy.signal import welch freqs, psd = welch(beat, fs=fs, nperseg=65536, noverlap=32768, window='hann', scaling='density')

第三个细节是频段截取。拟合时不需要把整个0到200MHz的频谱都交给拟合器，只取峰值附近±2MHz的区间就够了。截取范围太大会把远处的噪声基底也纳入拟合，太小又会截断洛伦兹的尾巴，导致拟合参数偏移。符号上，先定位峰值索引，然后从频谱中切出以峰值中心、宽度约2~3MHz的数据段。

3. 拟合代码的核心：洛伦兹线型、2倍关系与参数初值

3.1 为什么拍频谱是洛伦兹，以及那个绕不开的2倍关系

仿真的相位噪声模型是维纳过程，这决定了拍频信号的功率谱是洛伦兹线型。理论上的推导并不复杂：相位差φ(t)-φ(t-τ)的统计性质由激光器线宽Δν决定，当延迟τ足够远大于相干时间时，拍频谱的形状为洛伦兹函数，且其FWHM等于2Δν。

这个2倍关系是DSH方法里最容易翻车的地方。很多人拟合出了拍频FWHM，直接把这个数当激光器线宽报出去，结果比真实线宽大了整整一倍。因为拍频过程相当于对光谱做了一次"自卷积"，谱线宽度自然扩展了2倍。

还有一点要注意半高全宽和半高半宽的定义。我用的拟合函数是：

S(f) = A / (1 + ((f - f0) / gamma)^2) + B

其中gamma是洛伦兹的半高半宽（HWHM），FWHM = 2×gamma。拍频FWHM = 2×激光线宽FWHM，这两个2倍关系一抵消，拟合出的gamma值在数值上恰好等于激光器的FWHM线宽Δν。也就是说，代码里不需要额外除以2，只要确认好定义，直接输出gamma就是待测激光线宽。这个容易记混，我每次跑完都会故意看一眼输出和设定值是否一致，用仿真来校验。

3.2 拟合流程与代码骨架

整个拟合程序按下面几步走：

1. 读入频谱数据（频率数组、功率谱数组）。
2. 在截取的窗口内找到峰值位置，作为f0初值。
3. 用峰值幅度减本底噪声，作为A初值。
4. 在谱线两侧找半高位置，估算gamma初值。
5. 调用scipy.optimize.curve_fit做非线性最小二乘拟合。
6. 输出拟合参数、协方差、相对误差。

核心拟合函数和主流程如下：

from scipy.optimize import curve_fit def lorentzian(f, A, f0, gamma, B): return A / (1 + ((f - f0) / gamma) ** 2) + B # 粗定位峰值 peak_idx = np.argmax(psd_masked) f0_guess = freqs_masked[peak_idx] B_guess = np.median(psd_masked[:50]) # 用窗口远端当作本底 A_guess = psd_masked[peak_idx] - B_guess # 粗估gamma：找半高位置，取半宽 half_level = (psd_masked[peak_idx] + B_guess) / 2 left_side = psd_masked[:peak_idx] right_side = psd_masked[peak_idx:] left_half_idx = peak_idx - np.argmin(np.abs(left_side[::-1] - half_level)) right_half_idx = peak_idx + np.argmin(np.abs(right_side - half_level)) gamma_guess = (freqs_masked[right_half_idx] - freqs_masked[left_half_idx]) / 2 # 拟合 popt, pcov = curve_fit( lorentzian, freqs_masked, psd_masked, p0=[A_guess, f0_guess, gamma_guess, B_guess], bounds=([0, f0_guess - 1e6, 1e3, 0], [np.inf, f0_guess + 1e6, 1e7, np.inf]), maxfev=10000 ) A_fit, f0_fit, gamma_fit, B_fit = popt print(f"设定线宽: {delta_nu/1e3:.2f} kHz") print(f"拟合线宽: {gamma_fit/1e3:.2f} kHz")

3.3 初值估计的意义和拟合边界

curve_fit本质上是迭代优化，初值给得不好，可能收敛到局部极值，甚至直接不收敛。经验上，A、B、f0这三个参数都比较容易估计，gamma是最敏感的。如果gamma初值比真实值大很多，拟合器可能把洛伦兹尾巴当成本底，把A和B调大来"解释"数据，得到明显偏大的gamma。所以上面的半高估算是值得做的，不能图省事直接写死一个固定值。

拟合边界也要设置好。gamma必须大于0，A必须大于0，B不设负值。f0的边界设定在初始猜测±1MHz内，防止拟合器在数据里乱跳。这里所有边界都是物理意义上的合理约束，而不是随意限制。如果初始窗口选得太宽，比如±10MHz，拟合时本底B的耦合会非常强，A和B相互竞争。用±2MHz或±3MHz通常比较合适。

加权拟合也是可选的优化项。如果频谱数据经过了Welch平均，各频点的方差近似服从卡方分布，理想加权方式是乘以对应频点的谱值倒数。实际操作中我用curve_fit的sigma参数传入1/sqrt(psd)做简单加权，可以让拟合更注重峰顶和近峰区域，而不是被远端噪声拉扯。不过这个不是必须的，没有加权也能拿到不错的结果，加权主要用于噪声偏大的数据。

4. 仿真暴露的三个真坑：延迟不足、频谱泄漏、噪声耦合

4.1 延迟时间不足导致的系统性偏窄

这是DSH仿真中最容易忽视的问题。前面说过，延迟时间τ要远大于相干时间τc，拍频谱才是干净的洛伦兹，FWHM才会严格等于2Δν。如果τ不够长，相位差φ(t)-φ(t-τ)还没充分演化到统计平衡，拍频谱的形状就会偏离洛伦兹，峰值区域变得更锐利，测出来的线宽会系统性偏窄。

我特意做了一组对比仿真，固定激光线宽50kHz不变，改变延迟时间，观察拟合结果的变化趋势。实际算下来：

注意这是我当前仿真参数下的典型结果，具体数值会因采样长度和拟合窗口改变而略有浮动，但趋势非常稳定：τ/τc小于1时，误差在10%以上；到了3倍左右，误差降到2%以内；超过5倍就基本可以忽略了。

这个结论对实测的指导意义很大。假如要测一台线宽10kHz的激光器，τc约32μs，理论上至少需要300μs的延迟，对应光纤长度约60km。实际工程中很多实验室只有10km或20km的延迟光纤，此时τ/τc只有3到6，测出来的线宽会偏窄几个百分点。仿真在这里的价值就是预先告诉你：你手里的光纤长度够不够，测出来的数该不该信。

4.2 洛伦兹重尾和频谱泄漏的组合拳

洛伦兹函数的尾巴按1/f²衰减，衰减得比高斯慢得多。这带来一个麻烦：如果FFT没有加窗，主瓣能量会因为矩形窗的卷积而向两侧泄漏，在谱线的远端形成明显的旁瓣结构。单次FFT时这些旁瓣像锯齿一样，会干扰本底估计，进而让B参数在拟合时偏离真实值。

我在不加窗和加汉宁窗两种条件下分别跑了同一组数据。不加窗时，2MHz窗口两侧的频谱底大约比真实本底高出一个数量级，拟合出的B明显偏大，Gamma从50kHz漂到53kHz左右，有6%的正偏差。加窗后，旁瓣被压下去，两侧频谱回到了接近噪声底的水平，Gamma回归正常。

这类问题在实测中几乎一定会遇到，因为实测数据的动态范围更复杂，探测器噪声、放大器的残余响应都会和泄漏交织在一起。仿真给我的经验是：永远不要不加窗就直接拟合；如果数据质量差，窗口外的裙边仍然很高，那就先做一个高通或带通滤波，把谱线附近的频段分离出来，再去拟合。

4.3 本底噪声耦合导致A、B、gamma参数互咬

洛伦兹拟合成败的关键往往不在峰顶端，而在本底B的处理。四个参数A、B、f0、gamma之间是强耦合的：A抬升会把峰拉高，B抬升会把整体垫高，gamma拉伸会把线型变宽，三种效果在频谱上可能产生相似的视觉结果。如果噪声很大，拟合器很容易通过调大B和gamma来同时"解释"噪声和真实线型，导致线宽偏大。

这个问题在仿真里可以精确控制地验证。我给拍频信号加了不同幅度的白噪声，结果噪声幅度从0升到信号幅度的1%时，拟合gamma从50kHz漂到55kHz；升到5%时直接飘到70kHz以上。这说明噪声本底对拟合结果的影响是显著的。

对策有三种：

• 拟合前用窗口两侧远离峰的区域做B的初值估计，并且在拟合时把B的边界限制在初值附近一个很窄的范围内，不让它自由漂移。
• 使用加权拟合，噪声大的频点权重小，避免远端噪声主导参数更新。
• 如果本底确实平坦且已知，就把B固定为常数，只拟合A、f0、gamma三个参数。这个最直接，也最稳。

我在最终的仿真程序里做了灵活性设计：config参数里有fit_fix_b这个开关，默认关掉，让B自由拟合；如果实测数据本底明显平坦，就打开固定B模式。这看似很小的设计，实际在处理低SNR数据时非常有用。

5. 验证边界：不同线宽量级下拟合程序的可靠区间

5.1 一组完整的验证结果

为了确认程序不是只在50kHz这一特定情况下好用，我拿不同线宽、不同延迟时间跑了一组完整矩阵。这里统一采用T=2ms、fs=400MHz、Welch分段平均65536点、汉宁窗、拟合窗口±2MHz。

三组结论很清楚。

• 线宽越宽，越容易测准。因为同样延迟时间下，宽线宽的激光器相干时间更短，τ/τc比值更大。
• 1kHz以下线宽需要非常长的延迟。按50μs延迟、τ/τc=6以上的标准反推，测1kHz线宽需要延迟约190μs，对应光纤长度约39km。很多实验室没有这个条件，所以超窄线宽测量要么依赖更大的延迟池，要么需要用零差法配合锁相环。
• 仿真和实测相比，最大的优势就是可以在几分钟内把这种边界试出来，不用真的去拉几十公里光纤。

5.2 为什么仿真结果不能直接等于实测结果

仿真数据是"干净"的，没有偏振衰落、没有探测器饱和、没有AOM驱动不稳、没有光纤温度漂移。实测数据里这些因素一个都不会少。所以仿真的结论不能说成"这个程序就是准的"，而应该说"在理想条件下，程序的反演算法本身是准的；实测偏差主要来自测量链路的物理限制，程序层面的误差已经被标定过了"。

这句话听起来像套话，但实际上非常重要。我见过不少人把仿真验证后的拟合程序直接拿去处理实测数据，结果发现线宽总比预期小一截，然后开始怀疑程序。其实不是程序错了，而是实测链路的延迟时间不够、噪声基底耦合严重，或者谱线被振动的旁瓣干扰。程序只是忠实反映了输入数据里的信息，数据本身被污染了，再准的拟合也救不回来。

5.3 扩展方向：从白噪声到1/f噪声、从洛伦兹到Voigt

很多窄线宽激光器的频率噪声并不是纯白噪声，而是混有1/f分量。这种情况下，光谱形状不再是纯洛伦兹，低频段会被1/f噪声拉起成更"尖"的形状，整体线型更像Voigt轮廓。如果你拿纯洛伦兹去拟合这类数据，得到的gamma会落在线型的谷底和峰顶之间，既不是洛伦兹分量也不是高斯分量，解释起来很尴尬。

仿真的扩展也简单。在相位噪声生成时，除了白噪声驱动项，再加一个低频分量，比如用指数自相关的颜色噪声叠加到相位中。此时就不该再强行用洛伦兹拟合，而是选用伪Voigt函数或Voigt函数，拟合参数

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