做结构健康监测和无损检测的同行应该都见过声发射波形试件受力、裂纹扩展、材料内部位错滑移时应变能瞬间释放产生一个持续时间很短的高频弹性波传感器拾取放大后屏幕上就是一段“突突”振铃衰减的波形。刚接触这行的人容易觉得声发射信号很难搞参数多、形态复杂实测信号里还混着噪声和反射。其实反过来想如果能在MATLAB里自己把单次声发射波形从头合成出来每个影响因素都能拆开看——中心频率怎么设、衰减快慢怎么写、噪声加到多大、包络长什么样一目了然。这篇内容就是我常用的做法用MATLAB模拟并绘制单个声发射波形示例的波形模板叫“Planck”包络采用(t/τ)e^(-t/τ)这种非对称脉冲形态完整流程从参数设定、信号合成到绘图全部走一遍。整套流程非常短适合信号处理入门者、声发射方向的研究生也适合想快速验证AE特征提取算法的工程师。1. 模拟前需要搞懂的几个关键概念1.1 声发射信号为什么是“振荡衰减”的形态声发射Acoustic Emission简称AE是材料内部因局部能量快速释放而产生瞬态弹性波的现象。金属材料里裂纹尖端滑移、复合材料里纤维断裂、压力容器里泄漏冲击都会产生AE事件。工程上拾取到的AE波形绝大多数是突发型burst信号幅值从噪声里快速升起到达峰值后迅速回落再经历若干次振荡慢慢衰减到零。这段波形本质上是被结构传递函数调制过的源信号但为了方便分析大家更习惯用一个标准形式去近似——高频正弦载波外面乘一条衰减包络。这种“指数衰减正弦振荡”的数学形态非常重要因为声发射特征参数都是从它身上提取的上升时间峰值前的时间、持续时间超过阈值的时长、振铃计数越过阈值振荡的次数、能量包络平方积分。我们在MATLAB里模拟波形其实就是在可控地制造这些参数。1.2 Planck型包络的数学定义与特性很多人看到“Planck”会想到量子论里的普朗克函数这里我借用这个名字指一种非对称脉冲包络。它的数学形式是h(t) (t / τ) · exp(-t / τ)其中τ是衰减时间常数。这个函数在t0时等于0没有突变到tτ时达到峰值峰值大小是e^(-1) ≈ 0.368之后按指数规律缓慢衰减。我常选择这种包络来模拟AE事件是因为真实突发型AE信号的物理过程通常是振动能量在极短时间里被激发到最大再逐渐耗散掉。这个“快上升、慢衰减”的形状比单纯的exp(-t/τ)更贴近实际因为纯指数包络在t0时幅值最大波形起始瞬间会有一个不自然的“台阶”而Planck包络从0平滑上升波形看起来干净真实。完整的信号模型是s(t) A · (t / τ) · exp(-t / τ) · sin(2π · f_c · t φ)其中A是峰值幅值fc是中心频率φ是初始相位。1.3 采样率与时间窗动手写代码前的必算项模拟前必须先确定采样率fs和时间窗T否则后面画出来的波形可能是假的。采样率要满足奈奎斯特采样定理即fs必须大于信号最高频率的两倍。对于窄带AE信号工程上更稳妥的建议是fs ≥ 10 × fc。本文示例给中心频率fc 30 kHz采样率我建议选500 kHz大约是中心频率的16.7倍。这个裕量一方面保证时域波形光滑便于包络和峰值读数另一方面也保证FFT频谱不会出现频谱混叠主频定位准确。时间窗T的选择取决于包络衰减速度。对τ0.3 ms来说包络在t5 ms时已经衰减到峰值幅值的百万分之几所以5 ms完全够用。对应采样点数为N fs × T 500,000 × 0.005 2500个点。点数不算多画图、计算FFT都很快这正是“单次波形模拟”的优势——改参数几乎零成本。提示如果让你模拟的信号中心频率提高到几百kHz记得按比例提高采样率。我习惯先算一个参数表把fs、T、fc、τ写好注释再动手写代码。2. MATLAB实操从零到一张完整的AE波形图2.1 先把参数表列清楚打开MATLAB编辑器新建一个脚本我习惯命名为ae_sim_planck.m。脚本最前面放一个参数区所有可调参数集中在这里加好注释。参数表如下参数取值含义fs500e3 Hz采样率满足fs ≥ 10×fcT5e-3 s信号总时长Nfs×T采样点数fc30e3 Hz中心频率tau0.3e-3 s包络衰减时间常数A1峰值幅值归一化phi0初始相位SNR_db25 dB加噪信噪比之所以把这些参数全部提到文件头部是因为后续做参数扫描时你只需要在一个地方改数字不需要翻代码逻辑。这也是我在工程脚本里的长期习惯参数区、信号合成区、绘图区三段式划分脚本几十次改动也不会乱。2.2 核心代码的逐段拆解第一段生成时间轴和包络代码是这样的% 参数区 fs 500e3; % 采样率 500 kHz T 5e-3; % 信号时长 5 ms N round(fs * T); % 采样点数 t (0:N-1) / fs; % 时间轴单位秒 % Planck 型 AE 信号参数 fc 30e3; % 中心频率 30 kHz tau 0.3e-3; % 衰减时间常数 0.3 ms A 1; % 峰值幅值 phi 0; % 初始相位 % 包络 env (t / tau) .* exp(-t / tau); env env / max(env); % 归一化让峰值幅度为 1这里有一个细节容易踩坑如果不做归一化原始包络在tτ处的峰值只有0.368而幅值参数A的全过程含义就会变得不直观。我把包络归一化到1后A就直接代表信号的峰值幅值调整A时心里非常有数。载波信号用sin还是cos区别只在于初始相位φ。对单次AE模拟来说φ0即可如果你要做多通道时延分析才需要考虑相位对齐问题。合成信号就是包络乘以载波再加高斯白噪声% 载波 carrier sin(2 * pi * fc * t phi); % 合成无噪信号 s A * env .* carrier; % 手动添加高斯白噪声按信噪比计算噪声功率 SNR_db 25; Ps mean(s.^2); Pn Ps / (10^(SNR_db / 10)); noise sqrt(Pn) * randn(size(s)); s_noisy s noise;为什么手动加噪而不直接用awgn函数因为awgn属于Communications Toolbox不是所有人都装了。手动加噪只用randn和基本运算任何MATLAB版本都能直接跑。计算过程就是先估计信号平均功率Ps再由信噪比定义求出噪声功率Pn Ps / 10^(SNR/10)最后生成对应标准差的高斯白噪声。这种写法通用、可控性强还能扩展成不同信噪比的多段信号。2.3 时域绘图把波形和包络画到同一张图上数据算完了绘图部分要做的就是把波形、上下包络、图例和坐标标注都安排好。我常用的绘图代码figure(Color, w); plot(t * 1e3, s_noisy, b, LineWidth, 0.8); hold on; plot(t * 1e3, env, r--, LineWidth, 1.5); plot(t * 1e3, -env, r--, LineWidth, 1.5); hold off; xlabel(时间 (ms)); ylabel(幅值); title(模拟声发射波形Planck 型包络); legend(AE信号噪声, 包络, 包络-负向); grid on; xlim([0, 5]);这里我把时间轴乘以1e3单位换成毫秒图上横轴显示的是0到5 ms直观好用。如果直接用秒坐标轴会显示0到0.005既不美观也不好读。把正负两个方向包络都画成红色虚线能一眼看出噪声对包络幅值的干扰程度。运行脚本后你会看到一条蓝色波形快速起振后逐渐衰减红色虚线正好包住振荡的上下边缘。这就是单人声发射波形的第一版可视化结果。我建议此时先不急着加噪声版本观察无噪信号s的波形确认包络峰值出现在t0.3 ms附近振荡频率肉眼数一下大约30 kHz。确认无误后再用加噪版本。3. 绘图进阶包络提取、频谱分析与参数对比3.1 用希尔伯特包络检查信号形态上一步画的包络是理论包络因为它是我们直接算出来的。但实际分析AE信号时我们往往只拿到波形数据需要从波形本身提取包络。最常用的方法就是希尔伯特变换把实信号变成解析信号取模就得到包络。MATLAB里一句话env_est abs(hilbert(s_noisy));注意hilbert函数属于Signal Processing Toolbox。如果没装这个工具箱可以退而求其次用“平方包络法”近似先对信号平方再用smoothdata做高斯平滑最后开方也能得到大致轮廓。不过实测下来希尔伯特包络更准对振荡型信号尤其好用。对加噪信号做希尔伯特提取时包络会有一点毛刺这是噪声导致的。解决思路分两步一是别用太小的时间窗比如这里5 ms的窗口希尔伯特包络依然稳定二是如果后续要做上升时间、持续时间、能量这些特征参数计算最好先对raw waveform做一个带通滤波滤掉带外噪声再提取包络。直接对高噪声信号提取包络测出来的上升时间会明显偏小因为噪声可能让包络提前触发阈值。3.2 用FFT验证主频是否正确波形画得再漂亮还得验证频域对不对。AE信号的特征频率是我们设定的中心频率30 kHzFFT应该在这个位置出现清晰的谱峰。代码如下L length(s_noisy); NFFT 2^nextpow2(L); % 补零到2的幂频谱更平滑 f fs / 2 * linspace(0, 1, NFFT/2 1); % 单边谱频率轴 S fft(s_noisy, NFFT); S_mag abs(S(1:NFFT/2 1)); % 取单边幅度 S_mag S_mag / max(S_mag); % 归一化 figure(Color, w); plot(f / 1e3, S_mag, b, LineWidth, 1.2); xlabel(频率 (kHz)); ylabel(归一化幅度); title(模拟AE信号的频谱); grid on; xlim([0, 100]);我一般用单边谱因为双边谱有一半是对称冗余分析主频没必要。NFFT补零不是增加真实频率分辨率只是让插值后的谱线更平滑真正的频率分辨率取决于原始时域长度这里T5 ms对应的分辨率是1/5 ms200 Hz30 kHz主频附近的峰完全够分辨。如果你的频谱峰值不在30 kHz附近最大概率是采样率设置出了问题或者载波频率参数写错。频谱验证是波形模拟里最值得做的自检项建议每次改参数后都跑一遍。3.3 参数扫描不同衰减时间常数与中心频率单次模拟只能看一条波形参数扫描能直观展示每个参数对波形的“雕塑力”。这里我给你一套快速对比的方法把tau和fc各准备几组值用for循环跑画在subplot里。tau_list [0.1e-3, 0.3e-3, 0.6e-3]; fc_list [15e3, 30e3, 60e3]; figure(Color, w); for i 1:3 tau tau_list(i); env_i (t / tau) .* exp(-t / tau); env_i env_i / max(env_i); s_i A * env_i .* sin(2 * pi * fc_list(i) * t); subplot(3, 1, i); plot(t * 1e3, s_i, b, LineWidth, 0.8); hold on; plot(t * 1e3, env_i, r--, LineWidth, 1.2); plot(t * 1e3, -env_i, r--, LineWidth, 1.2); hold off; xlim([0, 5]); ylim([-1.2, 1.2]); title(sprintf(tau%.1fms, fc%.0fkHz, tau*1e3, fc_list(i)/1e3)); grid on; end xlabel(时间 (ms));这段代码能跑出三行波形对比。第一行tau0.1 ms波形振铃几下就没了持续时间很短第三行tau0.6 ms振荡能拖到3 ms以后。这就是衰减时间常数对AE持续时间最直观的体现。如果你把fc提得更高比如60 kHz务必确认fs满足10倍采样关系否则时域波形会出现明显的“锯齿”失真。参数扫描是理解AE信号形态最好的实验方法比看十页公式都有用。4. 常见问题与避坑指南4.1 采样率不足波形变成“假的”我见过不少新手把fs设为40 kHz就去模拟30 kHz的AE信号结果画出来波形面目全非——频率不对、幅度忽大忽小。原因很简单40 kHz采样率对30 kHz信号只比奈奎斯特频率大10 kHz虽然没超过香农下限但可用的频带余量太薄波形重建效果极差。对于需要读峰值时序的应用强烈建议fs至少是中心频率的10倍。这里我的实践是模拟30 kHz信号就用500 kHz采样模拟150 kHz信号就用1.5 MHz采样。代价仅仅是数组变长对单次波形模拟来说完全不是负担。4.2 包络线毛刺与端点效应希尔伯特包络在信号两端经常出现“飞刺”显得包络突然抬高。这是解析信号处理在有限长度数据上的端点效应。解决办法并不复杂一是保留足够长的静默段让信号在起始和末尾都衰减到接近0二是在提取包络后对前后几十个点做简单截断或平滑处理。我个人在画图时会直接忽略端点两端约5%的区域因为AE特征分析本身也常常设置阈值不会把起始前和衰减后的微弱信号当有效数据。4.3 单位与标注细节这个坑虽小但容易让人困惑时间用秒还是毫秒频率用Hz还是kHz。我建议代码里全部用国际单位秒、Hz计算只有绘图时乘以1e3把横轴显示成毫秒频谱横轴除以1e3显示成kHz。这样既符合MATLAB习惯图上标注又友好。另外legend里的文字要写清楚区分“理论包络”和“估计包络”不然以后回看脚本很容易分不清哪条线是哪个来源。4.4 工具箱与版本兼容问题不同MATLAB版本对这个脚本没有兼容性问题核心运算全是fft、hilbert、plot、randn这些基础函数。唯一要注意的是plot语法在2014之后都是一致的旧版R2010之前略有差异但内容已经极少人用。新版MATLAB包括2026b打开这个脚本直接运行即可。如果你用envelope函数替代我这里的近似方法需要注意envelope属于Signal Processing Toolbox如果你用awgn加噪需要Communications Toolbox。我的原则是能用基础函数就不用工具箱函数保证脚本在纯净环境下也能跑通。这里整理一个常见问题速查表现象可能原因解决方法波形出现锯齿、频率失真采样率不够高fs设为10倍以上频率频谱主峰不在30 kHzfc参数写错或混叠检查参数重算NFFT包络两端飞刺希尔伯特端点效应加静默段、截断端点区域包络毛刺严重噪声太强先带通滤波再提取包络纵轴幅值超出想象未归一化包络加env env / max(env)报错找不到awgn缺少通信工具箱手动加噪sqrt(Pn)*randn(size(s))4.5 “Planck”模型还能怎么扩展单个波形模拟清楚之后扩展方向很多。可以把多个Planck波形按时间间隔拼接模拟连续多个AE事件可以让包络参数随机化生成一批波形用于训练神经网络分类器也可以加入声波传播时延的概念在传感器阵列模拟里用它做数据源。我个人实际使用中最喜欢拿它当“信号标准器”——在写AE特征提取程序振铃计数、上升时间、能量时先用模拟波形验证算法逻辑确保每个特征都能精确还原再去处理真实采集信号。这套工作流省掉的调试时间非常可观。最后分享一个小技巧把脚本开头参数区的SNR_db从25 dB改成60 dB跑一次就能看到无噪情况下的纯Planck波形再改成5 dB跑一次会看到信号几乎淹没在噪声里。这两条波形的对比比任何文字都能让你直观感受到信噪比对AE信号分析的影响。多存几个不同信噪比的输出图理解波形形态时非常有价值。