★ FEATURED ARTICLE 简介这份资源是面向通信工程、电子信息工程等专业本科生与研究生的信号与系统综合实验报告围绕连续时间LTI系统展开时域与频域特性分析及滤波器设计。内容涵盖稳定系统的零极点设计、直接/级联/并联三种实现结构的信号流图与仿真对比、参数微扰对零极点分布与稳定性的敏感性分析、周期信号频谱特性、LTI滤波器设计与音频信号频谱分析等核心环节并配有Matlab代码实现与结果讨论。资源包共1个docx文件约2.62MB内含完整实验目的、原理推导、步骤说明与结果分析结构清晰可直接用于课堂教学或自主学习。目前已有1303人学习下载。读者可借此掌握系统函数零极点与频率特性、稳定性之间的关系理解不同实现结构对参数变化的敏感度差异并积累Matlab仿真与音频信号处理的实操经验。1. 连续时间 LTI 系统综合实验从时域卷积到频域滤波的一条完整链路很多同学做连续时间 LTI 系统实验时习惯把时域分析和频域分析拆成两个互不相干的脚本一个画冲激响应一个画幅频特性交完报告就完事。真正到了滤波器设计环节却发现前面算出来的东西根本接不上——截止频率对不上、卷积结果和理论值差一个尺度、仿真出来的波形怎么看都不像低通。问题不在 Matlab 用得熟不熟而在于没有把「时域卷积—频域相乘—滤波器实现」当成一条链路来跑。这个综合实验要解决的就是这件事用一套可复现的 Matlab 流程把连续时间 LTI 系统的冲激响应、卷积、频率响应、滤波特性串起来最后落到一个能实际跑出结果的滤波器设计上。适合正在做信号与系统课程实验的本科生也适合想快速把理论公式变成可运行代码的工程师。下面按「先立住理论、再动手复现、最后避坑」的顺序展开。2. 时域特性分析冲激响应与卷积积分的 Matlab 落地2.1 为什么先做时域LTI 系统的「身份证」是冲激响应连续时间 LTI 系统的核心结论只有一条系统对任意输入的响应等于输入与冲激响应的卷积。也就是说只要拿到 h(t)这个系统就被完全确定了。时域分析的价值在于它不依赖变换域的任何假设直接给出物理上可解释的结果——输入信号经过系统后每个时刻的输出是过去所有输入的加权叠加。在 Matlab 里做时域分析绕不开两个问题一是连续信号怎么离散化二是卷积积分怎么用数值方法近似。常见做法是用足够小的采样间隔 dt 把连续信号采样成序列然后用数值积分或离散卷积近似。这里的关键参数是 dt 和观测时间窗口 T它们直接决定结果的精度和计算量。我一般会先确定信号的有效带宽或变化最快的时间尺度然后取 dt 小于该尺度的十分之一。比如信号最高频率分量约 10 Hz那 dt 取 0.001 s 就比较稳妥。观测窗口 T 要覆盖冲激响应衰减到可忽略的时间否则卷积结果尾部会被截断。2.2 用数值积分实现卷积代码与参数说明下面这段代码实现了一个通用的一阶系统冲激响应和卷积计算。系统用微分方程描述先求冲激响应再对给定输入做卷积。% 连续时间LTI系统时域分析冲激响应与卷积 clear; clc; close all; % 参数设置 dt 0.001; % 采样间隔决定数值精度 T 5; % 总观测时间 t 0:dt:T; % 时间向量 % 定义系统一阶系统 tau*dy/dt y x tau 0.5; % 时间常数 h (1/tau) * exp(-t/tau); % 冲激响应 h(t) % 定义输入信号矩形脉冲 x double(t 1 t 2); % 数值卷积近似连续卷积 y dt * conv(x, h); % 乘以dt是连续卷积的数值近似 y y(1:length(t)); % 截取与t等长的部分 % 绘图 figure; subplot(3,1,1); plot(t, x); title(输入信号 x(t)); xlabel(t/s); subplot(3,1,2); plot(t, h); title(冲激响应 h(t)); xlabel(t/s); subplot(3,1,3); plot(t, y); title(输出 y(t)x(t)*h(t)); xlabel(t/s);这段代码的逻辑是先构造时间向量和冲激响应再构造输入最后用conv做离散卷积并乘以 dt 来近似连续卷积积分。dt越小近似越准但计算量越大T要足够长让 h(t) 衰减到接近零否则卷积尾部会失真。tau是系统时间常数它决定了冲激响应的衰减速度也决定了系统带宽——tau 越小系统响应越快带宽越宽。运行后你会看到输出 y(t) 是一个被一阶系统「平滑」过的矩形脉冲上升沿和下降沿都变成了指数形状。这就是时域分析最直观的价值你能直接看到系统对输入做了什么。2.3 时域分析里最容易翻车的三个参数第一个是 dt。很多人直接取 0.01 甚至 0.1结果高频成分被严重混叠卷积出来的波形和理论值差很远。判断方法很简单把 dt 减半再跑一次如果结果明显变化说明 dt 还不够小。第二个是卷积后的截取位置。conv输出的长度是length(x)length(h)-1直接取前length(t)个点相当于只保留了卷积结果的前半段。正确做法是取中心附近与 t 对应的部分或者把时间轴对齐后截取。上面代码里用y(1:length(t))是一种简化处理实际做的时候要检查卷积结果的峰值位置是否合理。第三个是冲激响应的幅度尺度。连续系统的冲激响应和离散序列之间差一个 dt 因子如果忘了乘 dt输出幅度会差几个数量级。这个坑几乎每个人都踩过血泪经验就是先用手算的简单例子验证幅度再跑复杂信号。3. 频域特性分析频率响应、Bode 图与频谱的对应关系3.1 从时域到频域频率响应到底在描述什么频率响应 H(jω) 是冲激响应 h(t) 的傅里叶变换它描述的是系统对不同频率正弦输入的幅度缩放和相位偏移。换句话说如果你给系统输入一个频率为 ω 的正弦波输出还是同频率的正弦波只是幅度乘以 |H(jω)|相位加上 ∠H(jω)。这就是 LTI 系统在频域最本质的性质。在 Matlab 里频率响应可以用解析法求也可以用数值法求。解析法适合简单系统比如一阶系统 H(jω)1/(1jωτ)数值法适合复杂系统或实验数据。常见做法是用freqs函数直接求模拟滤波器的频率响应或者用 FFT 对冲激响应做数值变换。我一般会先用解析法画出理论曲线再用数值法验证两者对上之后再往下做滤波器设计。这样做的原因是滤波器设计阶段会大量依赖频率响应的形状如果频域分析本身不可靠后面全是空中楼阁。3.2 用 freqs 和 FFT 交叉验证频率响应下面代码同时用解析法和数值法计算一阶系统的频率响应并画在一起对比。% 频域分析解析法 vs 数值法 clear; clc; close all; tau 0.5; w logspace(-1, 2, 500); % 频率范围0.1到100 rad/s % 方法一解析法 H_analytic 1 ./ (1 1j*w*tau); % 方法二数值法对冲激响应做FFT dt 0.001; t 0:dt:10; h (1/tau) * exp(-t/tau); N 2^nextpow2(length(h)); H_fft fft(h, N) * dt; % 乘以dt近似连续傅里叶变换 w_fft (0:N-1) * (2*pi/(N*dt)); % 只取正频率部分对比 w_fft_pos w_fft(1:N/2); H_fft_pos H_fft(1:N/2); % 绘图对比 figure; subplot(2,1,1); semilogx(w, 20*log10(abs(H_analytic)), b, LineWidth, 1.5); hold on; semilogx(w_fft_pos, 20*log10(abs(H_fft_pos)), r--, LineWidth, 1.5); legend(解析法, FFT数值法); title(幅频特性); ylabel(幅度/dB); grid on; subplot(2,1,2); semilogx(w, angle(H_analytic)*180/pi, b, LineWidth, 1.5); hold on; semilogx(w_fft_pos, angle(H_fft_pos)*180/pi, r--, LineWidth, 1.5); legend(解析法, FFT数值法); title(相频特性); ylabel(相位/度); xlabel(频率/(rad/s)); grid on;这段代码的关键在于H_fft fft(h, N) * dt这一步。FFT 默认是对离散序列做变换要近似连续傅里叶变换必须乘以 dt。另外N取 2 的幂次是为了 FFT 计算效率nextpow2自动找到不小于序列长度的最小 2 的幂。频率轴w_fft的构造要注意2*pi/(N*dt)是频率分辨率乘以索引得到实际频率。运行后你会看到两条曲线在低频段几乎重合高频段 FFT 结果会有轻微偏差这是数值截断和混叠造成的。如果偏差很大优先检查 dt 是否足够小、观测时间是否足够长。3.3 频域分析中三个容易忽略的细节第一个是频率轴的单位。freqs用的是 rad/s而很多实验报告要求 Hz换算关系是 ω2πf。如果画图时忘了换算截止频率会差 6.28 倍滤波器设计直接跑偏。第二个是 FFT 的对称性。实信号的频谱是共轭对称的正频率和负频率幅度相同、相位相反。画图时只取正频率部分是对的但要注意索引范围是1:N/2还是1:N/21差一个点在某些频率上会造成明显误差。第三个是数值傅里叶变换的截断效应。如果 h(t) 在观测窗口结束时还没衰减到零相当于给信号加了一个矩形窗频谱会出现波纹。解决办法是增大 T或者给 h(t) 加一个平滑窗函数。我一般会先画 h(t) 的尾部确认它已经接近零再往下做。4. 滤波器设计从频率响应指标到可实现的传递函数4.1 滤波器设计的本质用有理函数逼近理想频率响应理想低通滤波器在截止频率处是垂直下降的但这样的频率响应对应的冲激响应是 sinc 函数非因果且无限长物理上不可实现。所以实际滤波器设计都是在做一件事用一个可实现的有理传递函数 H(s) 去逼近理想特性同时满足通带波纹、阻带衰减、过渡带宽度等指标。常见的逼近方法有 Butterworth、Chebyshev、Elliptic 等。Butterworth 通带最平坦但过渡带最宽Chebyshev 过渡带更窄但有通带波纹Elliptic 过渡带最窄但通带和阻带都有波纹。选哪种取决于你的指标要求如果对通带平坦度要求高选 Butterworth如果要求过渡带尽量窄选 Elliptic。在 Matlab 里模拟滤波器设计主要用butter、cheby1、cheby2、ellip这几个函数它们返回分子分母系数再用freqs验证频率响应。4.2 用 butter 设计低通滤波器并验证指标下面代码设计一个截止频率 100 Hz、阶数可调的低通滤波器并验证其幅频特性是否满足要求。% 模拟低通滤波器设计Butterworth clear; clc; close all; fc 100; % 截止频率单位Hz wc 2*pi*fc; % 转换为rad/s order 4; % 滤波器阶数 % 设计模拟低通滤波器 [b, a] butter(order, wc, s); % s表示模拟滤波器 % 计算频率响应 f logspace(0, 3, 1000); % 1Hz到1000Hz w 2*pi*f; H freqs(b, a, w); % 绘图验证 figure; subplot(2,1,1); semilogx(f, 20*log10(abs(H)), LineWidth, 1.5); hold on; yline(-3, r--, LineWidth, 1); % -3dB参考线 xline(fc, g--, LineWidth, 1); % 截止频率参考线 title(低通滤波器幅频特性); ylabel(幅度/dB); grid on; subplot(2,1,2); semilogx(f, angle(H)*180/pi, LineWidth, 1.5); title(相频特性); ylabel(相位/度); xlabel(频率/Hz); grid on; % 验证-3dB点 [H_fc, ~] freqs(b, a, wc); fprintf(截止频率处幅度%.4f dB\n, 20*log10(abs(H_fc)));这段代码的核心是butter(order, wc, s)其中order是滤波器阶数wc是截止角频率s表示设计模拟滤波器而不是数字滤波器。阶数越高过渡带越窄但相位非线性也越严重实际实现时对元件精度要求也越高。运行后你会看到幅频曲线在 100 Hz 处正好下降 3 dB这就是 Butterworth 滤波器的定义。如果 -3dB 点不在 100 Hz检查wc是否用了 rad/s 而不是 Hz。这个单位坑几乎每届学生都会踩一次。4.3 阶数选择与过渡带指标的权衡滤波器阶数不是越高越好。阶数每增加 1过渡带衰减斜率大约增加 6 dB/倍频程但同时也带来两个问题一是相位非线性加剧群延迟在截止频率附近变化剧烈二是对元件容差更敏感实际电路实现时更容易振荡。我一般会先根据过渡带指标估算所需阶数。比如要求 100 Hz 截止、200 Hz 处衰减至少 40 dBButterworth 的衰减公式是 20n·log10(f/fc) dB代入得 20n·log10(2)≈6n dB要 40 dB 就需要 n≈7。然后取 n7 或 8 跑一下看实际曲线是否满足再决定是否换 Chebyshev 来降低阶数。如果对相位线性度有要求比如要传输脉冲信号那就要考虑 Bessel 滤波器它的通带群延迟最平坦但过渡带最宽。选型没有绝对优劣只有是否匹配你的指标。5. 避坑与排查连续时间 LTI 系统实验里最常见的五个翻车点5.1 卷积结果幅度不对先查 dt 有没有乘现象时域卷积出来的输出幅度比理论值大几百倍或小几百倍。原因连续卷积的数值近似是y dt * conv(x, h)很多人直接写conv(x, h)忘了乘 dt。解决先用手算的简单例子验证幅度比如 x 和 h 都是矩形时卷积结果是三角形峰值等于重叠面积用这个来校准。5.2 频率响应曲线在高频段翘起来混叠或截断现象FFT 算出来的幅频特性在高频段不降反升。原因要么 dt 太大导致高频混叠要么 h(t) 截断太早导致频谱泄漏。解决把 dt 减半再跑一次如果曲线改善说明是混叠把 T 加倍再跑一次如果改善说明是截断。两个都试直到曲线稳定。5.3 滤波器截止频率对不上Hz 和 rad/s 混用现象设计 100 Hz 低通结果 -3dB 点在 628 Hz。原因butter的截止频率参数默认是 rad/s直接填了 100 而不是 2π×100。解决统一在 rad/s 下计算画图时再换算成 Hz。养成习惯所有频率变量加后缀_hz或_w避免混用。5.4 阶数太高导致数值不稳定检查极点位置现象高阶滤波器比如 10 阶以上的频率响应出现异常尖峰或凹陷。原因高阶传递函数的极点可能非常靠近虚轴数值计算时精度不够。解决用二阶节级联SOS形式代替直接形式Matlab 里用zp2sos转换。或者降低阶数换用 Chebyshev 或 Elliptic 逼近。5.5 时域和频域结果对不上检查变换对是否一致现象时域卷积得到的输出做 FFT 后和「输入频谱×频率响应」对不上。原因时域和频域用了不同的 dt 或不同的观测窗口。解决时域和频域分析必须用同一套时间参数先确定 dt 和 T再分别做时域和频域计算最后交叉验证。如果对不上优先检查两边的 dt 是否一致。6. 进阶技巧用零极点分布快速判断系统特性6.1 零极点图比频率响应更直观的系统「指纹」频率响应画出来是一条曲线但零极点图是几个点。别小看这几个点它们决定了系统的一切极点位置决定稳定性零点位置决定阻带特性零极点之间的距离决定频率响应的形状。对于连续时间 LTI 系统极点必须在左半平面系统才稳定极点越靠近虚轴对应频率处的增益越高零点越靠近虚轴对应频率处的增益越低。我一般会在设计滤波器之前先画零极点图确认极点都在左半平面再往下做。如果发现极点在右半平面不用往下跑了系统本身就不稳定。6.2 用 zplane 和 pzmap 快速检查稳定性下面代码设计一个滤波器并画出零极点图同时判断稳定性。% 零极点分析稳定性快速判断 clear; clc; close all; fc 100; wc 2*pi*fc; order 4; [b, a] butter(order, wc, s); % 求零极点 [z, p, k] tf2zp(b, a); % 画零极点图 figure; zplane(z, p); title(零极点分布); % 判断稳定性 if all(real(p) 0) fprintf(系统稳定所有极点都在左半平面\n); else fprintf(系统不稳定存在右半平面极点\n); end % 打印极点位置 fprintf(极点位置\n); disp(p);这段代码用tf2zp从传递函数系数求零极点zplane画图然后检查所有极点的实部是否小于零。Butterworth 滤波器的极点均匀分布在左半平面的半圆上阶数越高极点越靠近虚轴这也是高阶滤波器数值敏感的原因。6.3 从零极点图反推频率响应的三个经验规则第一个规则极点离虚轴越近该频率处增益越高。如果极点正好在虚轴上该频率处增益无穷大这就是振荡器。第二个规则零点离虚轴越近该频率处增益越低。如果零点在虚轴上该频率处增益为零这就是陷波器。第三个规则零极点成对出现时如果零点在极点左侧高频段衰减更快如果零点在极点右侧高频段可能出现增益回升。这三个规则不用算频率响应就能快速判断系统的大致形状。我一般在设计滤波器时先用零极点图确认大致形状再用freqs精确验证。这样比直接调参数快得多也不容易陷入「改一个参数跑一次仿真」的循环。最后说一个我自己的习惯每次做完滤波器设计都会把零极点图、幅频特性、相频特性、冲激响应四张图放在一起看。如果四张图讲的故事一致这个设计基本就靠谱了如果哪张图看起来不对劲优先查那张图对应的参数。这个习惯帮我省了很多后悔药。希望帮到你。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助？
