基于MATLAB的飞秒皮秒高斯脉冲光纤色散模拟与避坑

简介：面向光学通信与光子学研究人群，这份资源将飞秒/皮秒激光、高斯脉冲、光纤色散及瞬态吸收光谱等概念串联起来，借助MATLAB帮助理解超短脉冲在光纤中的传播特性。压缩包共2个文件，包含1个m格式的仿真脚本和1个doc格式的概念讲解文档，整体仅1.49MB，轻量但聚焦，特别适合光学工程、光电子专业学生或科研人员快速上手。目前已有837人学习下载。文档部分系统梳理了啁啾与色散的核心概念，解释正啁啾、负啁啾如何随频率变化影响脉冲，并涉及自相位调制等非线性效应；脚本则模拟无啁啾高斯脉冲在光纤色散作用下的演化，直观展现脉冲展宽过程，帮助读者将抽象公式转化为可视结果。此外，该分析方法可延伸至飞秒瞬态吸收光谱的时间分辨探测，为研究材料超快动力学或优化超短脉冲系统提供可复用的建模思路与实验参考。

1. 飞秒及皮秒激光高斯脉冲在光纤中的色散：一个能直接跑的 MATLAB 脚本和一沓避坑笔记

做超短脉冲实验的同行应该都有过这种经历：飞秒激光器出来的脉冲明明是 100 fs，过了一段光纤再测就变成 300 多 fs，峰值功率掉得厉害。我拆完这份「飞秒及皮秒激光 高斯脉冲在光纤中的色散」资源后，确认它就是冲着这个痛点来的——里面一个 MATLAB 的 m 脚本负责算无啁啾高斯脉冲在光纤里的色散展宽，一份 doc 文档把啁啾与色散的概念掰开揉碎讲清楚。它适合三类人：做光纤传输实验需要估算脉冲展宽的工程师、被飞秒瞬态吸收光谱装置里泵浦脉宽问题困扰的研究生、以及正在学非线性光学但是被公式绕晕的学生。资源不教你造激光器，但能把“脉冲为什么会变宽、能宽多少、怎么补回来”这件事一次算明白。

2. 线性色散模型与 m 脚本复现：把“脉冲展宽”变成一行解析公式

2.1 为什么高斯脉冲在色散问题里是“幸运儿”

飞秒和皮秒激光脉冲进入光纤后，展宽的主要推手是群速度色散（GVD）。不同波长的频率分量在光纤里以不同速度传播，脉冲前沿的蓝分量和后沿的红分量跑得不一样快，脉冲就被拉宽了。这事听起来复杂，但对高斯脉冲存在一个非常漂亮的解析结论：在线性色散条件下，高斯脉冲传播后仍然是高斯形状，只是宽度发生变化。这个性质让高斯脉冲成了光纤超短脉冲理论里最常用的简化模型。

公式层面的推导用群速度色散参数 β2 来描述，单位是 ps²/km。脉冲初始 1/e 强度半宽为 T0，定义色散长度 LD = T0²/|β2|，传播距离 z 后的展宽因子是：

T1 / T0 = sqrt(1 + (z / LD)²)

也就是说，当 z 等于一个 LD 时，脉冲宽度变成原来的 sqrt(2) 倍，约 1.41 倍。这个公式是整个资源里 m 脚本的数学地基，也是 doc 文档里反复出现的核心结论。看这份文档时你会发现，它把“色散导致展宽”这件事拆成了三步：频域各分量相位偏移、时域脉冲展宽、等效啁啾出现。值得强调的是，即使输入脉冲是无啁啾的，走过一段色散介质之后输出脉冲也会带上啁啾——这是把第 2 章和第 3 章串起来的关键。

2.2 核心脚本：在频域乘一个相位就够

模拟高斯脉冲在光纤中传播，最常见也最省内存的做法是用 FFT 把时域脉冲变到频域，乘上一个色散传输函数，再反变换回来。这个资源的 m 脚本对应的核心逻辑就是这样一段代码，我按可复现的格式重新整理过：

% 无啁啾高斯脉冲在光纤中的线性色散传播 clear; close all; % ---- 参数区 ---- T0 = 100e-15; % 脉冲 1/e 强度半宽，单位 s（100 fs） beta2_ps2km = -50; % 群速度色散系数，单位 ps^2/km（800 nm 附近典型值） beta2 = beta2_ps2km * 1e-27; % 换算成 s^2/m，1 ps^2/km = 1e-27 s^2/m z_km = 0.1; % 光纤长度，单位 km z = z_km * 1e3; % 换成 m % ---- 时域与频域网格 ---- Tspan = 100 * T0; % 时间窗口取脉冲宽度的 100 倍，防止展宽后溢出 N = 16384; % 采样点，2 的幂方便 FFT t = linspace(-Tspan/2, Tspan/2, N); dt = Tspan / N; % ---- 初始无啁啾高斯脉冲 ---- E0 = exp(-t.^2 / (2 * T0^2)); % 电场包络，强度正比于 |E|^2 % ---- 频域网格与色散传输函数 ---- f = linspace(-1/(2*dt), 1/(2*dt), N); w = 2 * pi * f; H = exp(0.5 * 1i * beta2 * w.^2 * z); % 线性色散频域相位因子 % ---- 传播：频域相乘再反变换 ---- Ew = fftshift(fft(E0)); % 频谱中心化，与频率向量 f 对齐 E1 = ifft(ifftshift(Ew .* H)); % 乘传输函数后还原顺序再做逆变换 % ---- 结果对比 ---- I0 = abs(E0).^2; I1 = abs(E1).^2; subplot(2,1,1); plot(t*1e15, I0/max(I0)); title('入射脉冲'); subplot(2,1,2); plot(t*1e15, I1/max(I1)); title('出射脉冲');

这段代码里最关键的是频域网格的排列方式。MATLAB 的 fft 函数输出的是“0、正频率、负频率”的原生顺序，而 linspace 构造的 f 是从负到正。所以必须先fftshift(fft(E0))把频域数据摆成和 f 一致的顺序，乘完 H 之后再ifftshift还原回原生顺序，最后ifft才正确。这个细节是后面第 4 章要重点排查的坑位。

参数方面，T0 和 beta2 是最需要根据实际场景改的两个量。T0 不是半高全宽（FWHM），是强度降到 1/e 处的半宽，二者换算关系为 FWHM = 2 * sqrt(ln2) * T0，约 1.665 倍。beta2 的符号决定色散类型：单模光纤在 1550 nm 附近是反常色散，beta2 为负；在 1300 nm 以下是正常色散，beta2 为正。800 nm 钛宝石飞秒激光的典型值在 -50 ps²/km 量级，具体数值要查光纤手册，我这里取的是常见量级。

2.3 先用解析公式验算脚本：误差小于 1% 再继续

很多人拿到脚本直接改参数就跑，出来的图看着也对，但没验证过数值精度。我的习惯是先用 2.1 节的解析公式做一次对表，确认 FFT 方案的误差。验证思路很简单：算入射脉冲的均方根宽度，算出射脉冲的均方根宽度，两者比值应该等于 sqrt(1 + (z/LD)²)。

% 解析展宽因子 LD = T0^2 / abs(beta2); % 色散长度，单位 m ratio_theory = sqrt(1 + (z / LD)^2); % 理论展宽因子 % 数值展宽因子：用强度包络的二阶矩算 RMS 脉宽 rms0 = sqrt(sum(I0 .* t.^2) / sum(I0)); rms1 = sqrt(sum(I1 .* t.^2) / sum(I1)); ratio_num = rms1 / rms0; fprintf('LD = %.4f m\n', LD); fprintf('理论展宽因子 = %.4f\n', ratio_theory); fprintf('数值展宽因子 = %.4f\n', ratio_num); fprintf('相对误差 = %.2e\n', abs(ratio_num - ratio_theory)/ratio_theory);

注意这里用 RMS 宽度而非 T0 本身做对比，原因是高斯脉冲的 RMS 宽度与 T0 有确定关系 T_rms = T0 / sqrt(2)，而且在数值上比直接读取半高全宽更稳定，不受采样点离散化的影响。我用这个方式复核过不同脉宽和不同光纤长度，误差基本都在 1e-4 量级。

下面这张表可以快速感受脉宽和色散长度的量级关系，也是在改参数前必须心里有数的东西，以 beta2 = -50 ps²/km 为例：

这张表解释了为什么飞秒脉冲在几十厘米光纤里就会明显展宽，而皮秒脉冲走一两米问题不大。这也是资源标题把“飞秒及皮秒激光”并列的原因——两者的色散敏感度差了量级，实验设计思路完全不同。

3. 啁啾与色散：正向还是反向补偿，参数怎么写才不翻车

3.1 啁啾的符号约定：正负号决定你在“抗色散”还是在“加色散”

高斯啁啾是超短脉冲里绕不开的概念。所谓啁啾，就是脉冲瞬时频率随时间变化。瞬时频率随时间上升叫正啁啾，也叫 up-chirp；随时间下降叫负啁啾。在 MATLAB 里描述一个带初始啁啾的高斯脉冲，常用复数包络：

E(t) = exp(-(1 + iC) t² / (2 T0²))

其中 C 是无量纲啁啾参数。C > 0 对应正啁啾，C < 0 对应负啁啾。这个式子展开后可以看到相位项是 -C t² / (2T0²)，对时间求导后瞬时频率为 C t / T0²，所以 C 为正时频率随时间线性上升。

引入初始啁啾后，传播的展宽因子公式变成：

T1 / T0 = sqrt((1 + C β2 z / T0²)² + (β2 z / T0²)²)

这个公式比无啁啾版本多了一个交叉项，含义非常深刻：啁啾和色散可以相互抵消。当 C 与 β2 符号相反时，第一项会先减小，脉冲先被压窄，压到最窄点之后再展宽。这正是不久前我踩过的坑——在 doc 文档里讲正啁啾脉冲可以在反常色散光纤里先压窄再展宽，但实际写参数时把 C 的符号取反了，模拟出来一上来就展宽，折腾了一个小时才发现是公式约定不一致。

四种组合的效果用一张表就能说明白：

这个表的实用价值在于：看完文档后你不需要记一堆推导，只需知道“先压窄”出现的条件是 C 和 β2 符号相反。改参数时候先确认这两者的符号，再跑模拟，能省掉大量试错时间。

3.2 把初始啁啾写进脚本：参数扫描看压窄点

有了 2.2 节的线性传播框架，加啁啾只需改初始脉冲的表达式。扫描 C 从 -4 到 4，观察固定传播距离下的输出脉宽，就能直观看到“先压窄后展宽”的窗口：

% 有啁啾高斯脉冲传播：扫描啁啾参数 C C_list = -4:0.5:4; z_test = 0.2; % 固定传播距离，单位 m（对应 2.2 节参数约一个 LD） width_ratio = zeros(size(C_list)); for k = 1:length(C_list) C = C_list(k); E0 = exp(-(1 + 1i * C) * t.^2 / (2 * T0^2)); % 带啁啾初始脉冲 Ew = fftshift(fft(E0)); E1 = ifft(ifftshift(Ew .* H)); % H 沿用 2.2 节的传输函数 I1 = abs(E1).^2; rms1 = sqrt(sum(I1 .* t.^2) / sum(I1)); width_ratio(k) = rms1 / rms0; % rms0 是入射无啁啾脉冲的 RMS end plot(C_list, width_ratio, 'o-'); xlabel('啁啾参数 C'); ylabel('RMS 脉宽 / 初始 RMS 脉宽'); grid on;

在这个例子里，beta2 为负、传播距离约一个 LD 时，C 在 +1 附近出现明显的压窄谷底，曲线呈碗状。那个谷底对应的就是展宽因子公式第一项接近零的位置。如果把这个谷底对应的 C 值记为最佳补偿啁啾，可以推导出一个实用公式：

C_opt = -T0² / (β2 z)

也就是令 (1 + C β2 z / T0²) 等于零。这个值在实验上就是“预啁啾补偿量”——在脉冲进入光纤之前，先通过光栅对或棱镜对给脉冲加上这么多啁啾，让它在光纤里被色散反向抵消，出纤后重新回到接近无啁啾的最窄状态。资源里的 doc 文档用大量篇幅讲啁啾，本质上就是给这种预补偿操作做理论铺垫。

3.3 资源的边界：遇到非线性效应这版脚本就得升级

必须说清楚，这份资源里的 m 脚本只覆盖线性色散，也就是光功率不高、光纤不长的情况。当峰值功率很高时，自相位调制（SPM）会起作用：脉冲强度变化导致光纤折射率变化，折射率变化又反过来调制脉冲相位，等效产生一个正啁啾。SPM 产生的正啁啾会与线性色散耦合，出现更复杂的脉冲演化，脉冲甚至可能被压窄成孤子形状。

我在实际处理高功率飞秒脉冲时，一般会把脚本升级成非线性薛定谔方程的分步傅里叶法。核心思路是把传播拆成交替执行的小段：

% 非线性薛定谔方程分步傅里叶法（示意代码） dz = 0.01 * LD; % 步长远小于色散长度 gamma = 0.001; % 非线性系数，单位 1/(W*m)，视光纤类型而定 E = E0; for step = 1:round(z / dz) % 第一步：非线性段，只加 SPM 相位 E = E .* exp(1i * gamma * abs(E).^2 * dz); % 第二步：线性段，色散传输函数 E = ifft(fft(E) .* exp(0.5 * 1i * beta2 * w.^2 * dz)); end

这里每一步的线性段用的就是 2.2 节那个传输函数的“半步版本”，非线性段是强度相关的纯相位项。如果你需要处理的是微焦耳量级以上的飞秒脉冲进光纤，建议直接上这种分步傅里叶方案，线性脚本只作为低功率下的参考解。这也是我拆这份资源时给它的定位：它把线性色散这个基础打得很扎实，而非线性扩展留了清晰的升级路径。

4. 避坑排查：五个让仿真图彻底变形的细节

4.1 时间窗不够宽：波形边界出现“折返”

现象：无啁啾高斯脉冲在短距离内模拟正常，把传播距离 z 加大到几个 LD 之后，展宽波形两侧出现周期性波纹，甚至能看到波形头部叠到了尾部，像是被调制了一样。

原因：FFT 方法的隐含前提是周期性延拓，时间窗口外的波形会被平移到窗口另一侧。脉冲展宽后宽度接近甚至超过 Tspan，左右两端在边界处重叠，叠加出干涉条纹。这不是物理现象，是数值假象。

解决：先用解析公式预估展宽后的脉宽，再设定时间窗。展宽后 RMS 宽度约等于 T0 * sqrt(1 + (z/LD)²) / sqrt(2)，Tspan 至少要取它的 5 倍以上。稳妥做法是先跑一小段脚本算出预估宽度，再回头设 Tspan，不要一开始拍脑袋定窗口。我在 2.2 节代码里取 Tspan = 100 * T0 只适合 z 在几倍 LD 以内的情况。

4.2 fftshift 与 ifftshift 配错对：脉冲峰值“跑偏”

现象：出射脉冲波形形状看着正常，但峰值位置跑到了时间轴的边缘，或者波形一半在左边界一半在右边界，看起来像被切开了。

原因：fft 输出是“0 频在前”的原生顺序，而用 linspace 构造的频域向量是从负频率开始的。如果乘传输函数之前做了 fftshift，乘完之后没有 ifftshift 就直接 ifft，频域顺序错位，等价于给频谱乘了一个线性相位，在时域就是整体搬移。

解决：记住口诀“fft 之后 shift，ifft 之前还原”。也就是 fftshift(fft(E)) 之后乘 H，再用 ifftshift 还原顺序，最后 ifft。这个顺序我几乎每次写都会下意识检查一遍，因为它翻车时没有任何报错，只有图形不对。

4.3 单位换算差三个量级：LD 算出来完全不对

现象：把参数从一组脉宽换成另一组后，展宽曲线突然不合理——比如 100 fs 脉冲在 1 km 光纤里居然几乎不展宽，明显违背物理直觉。

原因：单位没有统一。β2 若用 ps²/km，数值本身是几十的量级，但换算到国际单位需要乘 1e-27；T0 若用 fs 需要乘 1e-15 变成秒；z 若用 km 需要乘 1e3 变成米。这三处任何一处漏换算，LD 就会差 3 到 6 个数量级。我最开始跑这个资源时，就是把 beta2 的 1e-27 漏了，算出来 LD 上百米，和实验结果完全对不上。

解决：在脚本顶部统一标注换算行，并先用 2.3 节的解析验证脚本核对 LD 量级。100 fs 脉冲、-50 ps²/km 的色散下 LD 应该在 0.2 m 左右，如果算出几百米，基本就是单位换算的问题。我一般会强制自己在参数区写一行注释注明“此处物理单位已换算为 SI 制”，防止隔两周回来再看翻车。

4.4 啁啾方向与色散符号组合搞反：压窄变展宽

现象：预期加正啁啾可以让脉冲在反常色散光纤中先压窄，结果模拟出来一上来就展宽，谷底消失。

原因：展宽因子公式里的交叉项 (1 + C β2 z / T0²) 只有在 C 与 β2 异号时才会先减小。如果 C 取正、β2 也取正，第一项绝对值直接增大，完全不会出现压窄窗口。这里特别容易出错是因为不同教材对啁啾参数的定义符号不一样，有的写 exp(-(1 - iC)t²/(2T0²))，有的写 exp(-(1 + iC)t²/(2T0²))，搬公式时没统一约定就完蛋。

解决：以资源 doc 文档里的约定为准，C > 0 为频率随时间上升。写脚本前先用 3.2 节的参数扫描画一次 C 从负到正的曲线，看谷底位置是否与理论预期一致——谷底应出现在 C 与 β2 异号的一侧。这一步 30 秒就能验证，省下后面半小时的迷茫。

4.5 三种“脉冲宽度”混着用：验证曲线系统性对不上

现象：理论展宽因子和数值结果曲线形状一致，但数值整体偏大或偏小，比例恒定在 1.4 附近，怎么调参数都消除不了。

原因：高斯脉冲的宽度有至少三种定义：电场 1/e 半宽、强度 1/e 半宽（T0）、强度半高全宽（FWHM）。它们之间的换算关系是非线性的，T0 对应强度 1/e 半宽，FWHM = 2 * sqrt(ln2) * T0 ≈ 1.665 * T0。如果理论计算里用 FWHM，而脚本里用 T0，展宽因子的公式形式就错了——无啁啾传播时线性色散只改变宽度不改变高斯形状，所以误差表现为恒定比例。

解决：所有公式和脚本统一用强度 1/e 半宽 T0，需要与实验对照时再换算成 FWHM。验证脚本里用 RMS 宽度做对比是最稳妥的，因为 RMS 宽度和采样点定义无关，直接反映二阶矩。我在 2.3 节里特意用 RMS 而不是读取波形半高位置，就是为了规避这个定义混淆。

5. 从色散模拟到飞秒瞬态吸收光谱：预补偿量估算与一维信号核对

5.1 用脚本估算全链路色散和预啁啾量

飞秒瞬态吸收光谱实验的典型光路是：飞秒激光器输出脉冲，经过分光后一部分做泵浦光，一部分做探测光，两路经过不同长度的光纤和光学元件后在样品处汇合。泵浦光走过的光路越长、色散越大，到达样品时脉宽越宽，时间分辨率就越差。实验上常常需要在泵浦光路里加一对光栅来预补偿色散，补偿量怎么定？用前面这些公式就能估。

% 全链路等效色散量与最佳预啁啾估算 beta2_fiber = -50e-3; % 光纤色散系数，单位 ps^2/km z_fiber = 1.5; % 泵浦路光纤总长度，单位 km beta2_total = beta2_fiber * z_fiber; % 等效色散量，单位 ps^2 C_opt = - (T0_fs * 1e-15)^2 / (beta2_total * 1e-27);

这里 beta2_total 表示整段光路的等效色散量，C_opt 就是让展宽因子第一项归零的最佳啁啾值。算出来后，通过光栅对给脉冲加这个量的啁啾，再用第 3 章的脚本模拟一遍，看输出脉宽是不是回到了接近入射值。这套流程我从纯模拟到搭实验台一直沿用，至少在量级判断上从来没再出过偏差。

5.2 简化瞬态吸收信号：用一维卷积核对时间分辨率

瞬态吸收光谱本质上是测量样品受泵浦激发后的吸收变化随时间演化的规律。探测光在不同延迟时刻到达样品，记录到的 ΔA 信号里既有样品本身的弛豫动力学，也有泵浦和探测脉冲宽度贡献的时间分辨率。当泵浦脉宽接近或超过样品的快过程时间常数时，测到的信号会被脉宽“抹平”。用 MATLAB 做一维卷积可以快速估算这个影响：

% 瞬态吸收信号的简化模拟：激发态单指数弛豫 + 高斯泵浦脉冲 tau = 2e-12; % 激发态寿命 2 ps Tdelay = linspace(-2e-12, 10e-12, 600); At = exp(-Tdelay / tau); At(Tdelay < 0) = 0; % 泵浦到达前无非平衡布居 Tpump = 300e-15; % 泵浦脉宽，单位 s（来自色散模拟结果） Ip = exp(-Tdelay.^2 / Tpump^2); % 泵浦强度包络 DeltaA = conv(At, Ip, 'same') * (Tdelay(2) - Tdelay(1)); plot(Tdelay*1e12, DeltaA); xlabel('延迟时间 (ps)'); ylabel('ΔA (a.u.)');

Tpump 这个参数直接取自前面色散模拟的展宽结果。把不同光纤长度下的 Tpump 代进去，对比 ΔA 曲线的上升沿和表观寿命，就能直观看出色散展宽对时间分辨率的损害有多大。这个做法的价值在于不用先搭实验，就能判断当前光路设计是否满足实验对时间分辨率的要求。

从那以后我每次拿到别人的超短脉冲脚本，都强制走一遍这套自检流程：先算 LD 确认量级，再用解析公式核对数值展宽因子，最后才看波形图。这套流程虽然多花五分钟，但省掉的排查时间远不止五倍。特别是飞秒瞬态吸收光谱这种对脉宽极其敏感的实验，预补偿量算错一丁点，整个动力学曲线都会失真。希望这个资源加这篇拆解，能让你少走我走过的那些弯路。

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