简介本资源是一套面向本科及硕士阶段科研学习的激光线宽仿真拟合工具基于Matlab实现DSHDelay 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级别的真实线宽差了四五个数量级根本没法看。后来改用光电探测器加频谱仪又发现直流附近全是弛豫振荡峰和低频噪声想要的那一点线宽信息早就被淹没了。最后只能老老实实搭延迟自外差DSHDelayed 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其中Δν是激光线宽的FWHMn是标准正态分布随机数。用这个模型生成的光场其功率谱恰好是洛伦兹线型线宽正好等于Δν。这个性质让生成已知线宽的仿真数据变得非常直接我设置Δν50kHz生成出来的数据理论上就是50kHz然后就可以去检验拟合程序能不能把它反演出来。注意这里使用复数基带表示。实际光频率ν0高达193THz但仿真完全不需要处理这么高的频率把光载波去掉只保留相位随时间的变化就够了。基带信号中光的幅度不变相位做随机游走这就是相干性退化的全部体现。2.2 时域信号生成与采样参数选择采样率和信号长度的选择是整个仿真能不能跑通的关键。拍频信号的中心频率是Ω也就是AOM移频量所以采样率至少要满足奈奎斯特条件。如果用80MHz的AOM采样率取fs400MHz比较合理给信号留出至少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的频率分辨率Δf1/T。要分辨50kHz线宽最好将分辨率控制在5kHz左右也就是T至少200μs。取T2ms分辨率0.5kHz这样拟合时线型轮廓不会被分辨率过度抹平。延迟时间τ必须远大于激光器相干时间τc1/(πΔν)。对于50kHz线宽τc≈6.4μs以FWHM50kHz计算τc 1/(π×50e3) ≈ 6.4μs20μ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, fsfs, nperseg65536, noverlap32768, windowhann, scalingdensity)第三个细节是频段截取。拟合时不需要把整个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是洛伦兹的半高半宽HWHMFWHM 2×gamma。拍频FWHM 2×激光线宽FWHM这两个2倍关系一抵消拟合出的gamma值在数值上恰好等于激光器的FWHM线宽Δν。也就是说代码里不需要额外除以2只要确认好定义直接输出gamma就是待测激光线宽。这个容易记混我每次跑完都会故意看一眼输出和设定值是否一致用仿真来校验。3.2 拟合流程与代码骨架整个拟合程序按下面几步走读入频谱数据频率数组、功率谱数组。在截取的窗口内找到峰值位置作为f0初值。用峰值幅度减本底噪声作为A初值。在谱线两侧找半高位置估算gamma初值。调用scipy.optimize.curve_fit做非线性最小二乘拟合。输出拟合参数、协方差、相对误差。核心拟合函数和主流程如下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]), maxfev10000 ) 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必须大于0A必须大于0B不设负值。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τc约6.4μs拟合线宽相对误差2μs0.31约35kHz-30%5μs0.78约43kHz-14%10μs1.56约47kHz-6%20μs3.13约49.2kHz-1.6%50μs7.81约50.1kHz0.2%注意这是我当前仿真参数下的典型结果具体数值会因采样长度和拟合窗口改变而略有浮动但趋势非常稳定τ/τ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这一特定情况下好用我拿不同线宽、不同延迟时间跑了一组完整矩阵。这里统一采用T2ms、fs400MHz、Welch分段平均65536点、汉宁窗、拟合窗口±2MHz。设定线宽延迟ττ/τc拟合线宽相对误差1kHz100μs0.31约0.73kHz-27%1kHz500μs1.57约0.94kHz-6%10kHz50μs1.57约9.4kHz-6%10kHz200μs6.28约10.0kHz0.1%50kHz20μs3.13约49.2kHz-1.6%50kHz50μs7.81约50.1kHz0.2%100kHz20μs6.28约100.2kHz0.2%1MHz10μs31.4约998kHz-0.2%三组结论很清楚。线宽越宽越容易测准。因为同样延迟时间下宽线宽的激光器相干时间更短τ/τc比值更大。1kHz以下线宽需要非常长的延迟。按50μs延迟、τ/τc6以上的标准反推测1kHz线宽需要延迟约190μs对应光纤长度约39km。很多实验室没有这个条件所以超窄线宽测量要么依赖更大的延迟池要么需要用零差法配合锁相环。仿真和实测相比最大的优势就是可以在几分钟内把这种边界试出来不用真的去拉几十公里光纤。5.2 为什么仿真结果不能直接等于实测结果仿真数据是干净的没有偏振衰落、没有探测器饱和、没有AOM驱动不稳、没有光纤温度漂移。实测数据里这些因素一个都不会少。所以仿真的结论不能说成这个程序就是准的而应该说在理想条件下程序的反演算法本身是准的实测偏差主要来自测量链路的物理限制程序层面的误差已经被标定过了。这句话听起来像套话但实际上非常重要。我见过不少人把仿真验证后的拟合程序直接拿去处理实测数据结果发现线宽总比预期小一截然后开始怀疑程序。其实不是程序错了而是实测链路的延迟时间不够、噪声基底耦合严重或者谱线被振动的旁瓣干扰。程序只是忠实反映了输入数据里的信息数据本身被污染了再准的拟合也救不回来。5.3 扩展方向从白噪声到1/f噪声、从洛伦兹到Voigt很多窄线宽激光器的频率噪声并不是纯白噪声而是混有1/f分量。这种情况下光谱形状不再是纯洛伦兹低频段会被1/f噪声拉起成更尖的形状整体线型更像Voigt轮廓。如果你拿纯洛伦兹去拟合这类数据得到的gamma会落在线型的谷底和峰顶之间既不是洛伦兹分量也不是高斯分量解释起来很尴尬。仿真的扩展也简单。在相位噪声生成时除了白噪声驱动项再加一个低频分量比如用指数自相关的颜色噪声叠加到相位中。此时就不该再强行用洛伦兹拟合而是选用伪Voigt函数或Voigt函数拟合参数本文还有配套的精品资源点击获取 延迟自外差法测量窄线宽激光器线宽的仿真与拟合
简介本资源是一套面向本科及硕士阶段科研学习的激光线宽仿真拟合工具基于Matlab实现DSHDelay 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级别的真实线宽差了四五个数量级根本没法看。后来改用光电探测器加频谱仪又发现直流附近全是弛豫振荡峰和低频噪声想要的那一点线宽信息早就被淹没了。最后只能老老实实搭延迟自外差DSHDelayed 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其中Δν是激光线宽的FWHMn是标准正态分布随机数。用这个模型生成的光场其功率谱恰好是洛伦兹线型线宽正好等于Δν。这个性质让生成已知线宽的仿真数据变得非常直接我设置Δν50kHz生成出来的数据理论上就是50kHz然后就可以去检验拟合程序能不能把它反演出来。注意这里使用复数基带表示。实际光频率ν0高达193THz但仿真完全不需要处理这么高的频率把光载波去掉只保留相位随时间的变化就够了。基带信号中光的幅度不变相位做随机游走这就是相干性退化的全部体现。2.2 时域信号生成与采样参数选择采样率和信号长度的选择是整个仿真能不能跑通的关键。拍频信号的中心频率是Ω也就是AOM移频量所以采样率至少要满足奈奎斯特条件。如果用80MHz的AOM采样率取fs400MHz比较合理给信号留出至少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的频率分辨率Δf1/T。要分辨50kHz线宽最好将分辨率控制在5kHz左右也就是T至少200μs。取T2ms分辨率0.5kHz这样拟合时线型轮廓不会被分辨率过度抹平。延迟时间τ必须远大于激光器相干时间τc1/(πΔν)。对于50kHz线宽τc≈6.4μs以FWHM50kHz计算τc 1/(π×50e3) ≈ 6.4μs20μ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, fsfs, nperseg65536, noverlap32768, windowhann, scalingdensity)第三个细节是频段截取。拟合时不需要把整个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是洛伦兹的半高半宽HWHMFWHM 2×gamma。拍频FWHM 2×激光线宽FWHM这两个2倍关系一抵消拟合出的gamma值在数值上恰好等于激光器的FWHM线宽Δν。也就是说代码里不需要额外除以2只要确认好定义直接输出gamma就是待测激光线宽。这个容易记混我每次跑完都会故意看一眼输出和设定值是否一致用仿真来校验。3.2 拟合流程与代码骨架整个拟合程序按下面几步走读入频谱数据频率数组、功率谱数组。在截取的窗口内找到峰值位置作为f0初值。用峰值幅度减本底噪声作为A初值。在谱线两侧找半高位置估算gamma初值。调用scipy.optimize.curve_fit做非线性最小二乘拟合。输出拟合参数、协方差、相对误差。核心拟合函数和主流程如下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]), maxfev10000 ) 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必须大于0A必须大于0B不设负值。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τc约6.4μs拟合线宽相对误差2μs0.31约35kHz-30%5μs0.78约43kHz-14%10μs1.56约47kHz-6%20μs3.13约49.2kHz-1.6%50μs7.81约50.1kHz0.2%注意这是我当前仿真参数下的典型结果具体数值会因采样长度和拟合窗口改变而略有浮动但趋势非常稳定τ/τ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这一特定情况下好用我拿不同线宽、不同延迟时间跑了一组完整矩阵。这里统一采用T2ms、fs400MHz、Welch分段平均65536点、汉宁窗、拟合窗口±2MHz。设定线宽延迟ττ/τc拟合线宽相对误差1kHz100μs0.31约0.73kHz-27%1kHz500μs1.57约0.94kHz-6%10kHz50μs1.57约9.4kHz-6%10kHz200μs6.28约10.0kHz0.1%50kHz20μs3.13约49.2kHz-1.6%50kHz50μs7.81约50.1kHz0.2%100kHz20μs6.28约100.2kHz0.2%1MHz10μs31.4约998kHz-0.2%三组结论很清楚。线宽越宽越容易测准。因为同样延迟时间下宽线宽的激光器相干时间更短τ/τc比值更大。1kHz以下线宽需要非常长的延迟。按50μs延迟、τ/τc6以上的标准反推测1kHz线宽需要延迟约190μs对应光纤长度约39km。很多实验室没有这个条件所以超窄线宽测量要么依赖更大的延迟池要么需要用零差法配合锁相环。仿真和实测相比最大的优势就是可以在几分钟内把这种边界试出来不用真的去拉几十公里光纤。5.2 为什么仿真结果不能直接等于实测结果仿真数据是干净的没有偏振衰落、没有探测器饱和、没有AOM驱动不稳、没有光纤温度漂移。实测数据里这些因素一个都不会少。所以仿真的结论不能说成这个程序就是准的而应该说在理想条件下程序的反演算法本身是准的实测偏差主要来自测量链路的物理限制程序层面的误差已经被标定过了。这句话听起来像套话但实际上非常重要。我见过不少人把仿真验证后的拟合程序直接拿去处理实测数据结果发现线宽总比预期小一截然后开始怀疑程序。其实不是程序错了而是实测链路的延迟时间不够、噪声基底耦合严重或者谱线被振动的旁瓣干扰。程序只是忠实反映了输入数据里的信息数据本身被污染了再准的拟合也救不回来。5.3 扩展方向从白噪声到1/f噪声、从洛伦兹到Voigt很多窄线宽激光器的频率噪声并不是纯白噪声而是混有1/f分量。这种情况下光谱形状不再是纯洛伦兹低频段会被1/f噪声拉起成更尖的形状整体线型更像Voigt轮廓。如果你拿纯洛伦兹去拟合这类数据得到的gamma会落在线型的谷底和峰顶之间既不是洛伦兹分量也不是高斯分量解释起来很尴尬。仿真的扩展也简单。在相位噪声生成时除了白噪声驱动项再加一个低频分量比如用指数自相关的颜色噪声叠加到相位中。此时就不该再强行用洛伦兹拟合而是选用伪Voigt函数或Voigt函数拟合参数本文还有配套的精品资源点击获取