简介针对信号与线性系统综合实验课程这份报告类资源面向通信工程、电子信息工程等专业已完成基础知识学习的本科生与研究生。内容围绕连续时间LTI系统展开以设置零极点的方式设计稳定的4至6阶系统分别用直接、级联、并联三种结构建流图并借助Matlab仿真时域冲激响应、阶跃响应与幅频、相频特性通过微小改变系数对比不同结构在零极点分布、稳定性与响应上的敏感性差异同时构建周期信号分析频谱设计LTI滤波器完成滤波验证并延伸到音频信号的获取、输出与频谱分析。资源包为单个docx文档约2.62MB含实验目的、原理推导、步骤、完整Matlab代码与结果分析可直接作为课程实验报告模板或自主学习参考。目前已有1303人学习适合需要对照代码复现仿真、理清零极点与频率特性关系、补齐频域卷积定理等知识点的读者。1. 从 s 平面上的一对极点说起这套连续时间 LTI 系统综合实验要交付什么通信、电子信息这类专业的综合实验里连续时间 LTI 系统的时域频域特性分析与滤波器设计几乎绕不开。题目看着温和给一个 4~6 阶系统在 s 平面上摆几对零极点用 Matlab 把冲激响应、阶跃响应、幅频相频画出来再设计一个滤波器把 cos(t)cos(2t) 里的某个分量压下去。真正动手才发现极点实部差 0.01谐振峰就能挪出一大截分母系数手抄错一位稳定性结论直接反过来零极点位置记反整段滤波结论全部作废。这套方案面向的是要交一份完整综合实验报告的人——从系统函数的构造、三种实现结构的等价性与敏感性对比到截断频谱的误差解释再到陷波型滤波器的零极点设计和音频信号的采集分析一条线走完每一段都留下能复现的代码和能写进结论的判据而不是画几张图就交差。2. 需求拆解与工具选型符号计算与数值仿真两条腿怎么分工2.1 三层任务结构从系统分析到频谱验证再到音频应用整套实验其实是三层递进的交付物每一层的产物形式不一样验收口径也不一样。第一层是基础部分设计一个 4~6 阶的稳定连续系统其中至少含一对靠近虚轴的共轭复数极点用直接型、级联型、并联型三种结构分别实现画出信号流图对每种结构做时域和频域仿真再做参数微扰比较三种结构对系数变化的敏感程度。这一层交付的是「一组对照图 一条敏感性结论」结论必须落到哪种结构在有限字长下更稳。第二层是拓展部分构造一个含多个频率分量的周期信号截取不同长度做频谱分析把截断后的谱和理论谱摆在一起解释差异然后用 s 平面零极点法设计滤波器滤掉其中一个分量验证滤波效果。这一层交付的是「误差解释 滤波前后对比」核心是把频域卷积定理讲通。第三层是应用设计用声卡录一段音频存成 wav读回来画波形和频谱。这一层交付的是「一段可回放的音频 一张能读出基频和共振峰的频谱图」顺带把采样定理落到实处。三层之间不是并列关系。第一层把系统函数、零极点、稳定性这条线打通第二层把这条线搬到信号侧用窗谱去解释泄漏第三层把连续时间的直觉换成离散的 FFT稍不留神就会在频率轴上翻车。写报告时如果按这三层组织章节读者也就是批改的人能顺着你的推导走而不是在一堆孤立截图里找逻辑。需要提醒的是三层任务的代码风格要统一都用脚本都清空工作区都显式给参数。中途换变量名、换符号变量声明方式最后拼报告时会非常痛苦。2.2 符号计算与数值函数的边界ilaplace/fourier 与 freqs/step 各管哪一段Matlab 里做信号与系统仿真有两条技术路线选错了会白折腾很久。符号路线用syms声明符号变量用ilaplace求拉普拉斯反变换拿冲激响应用fourier求傅里叶变换拿频谱用ezplot或fplot画图。优点是表达式干净能直接写进报告的公式栏缺点是含heaviside的符号表达式在求fourier时经常返回一堆带pi、dirac的东西画出来还不是你想要的样子。数值路线用tf、freqs、impulse、step、lsim输入是系数向量输出是数值数组。优点是稳定、快、不会出符号怪缺点是要自己管频率轴和时间轴采样点选错了图形会失真。我的做法是两条腿都留着推导和公式展示走符号路线图形和定量指标走数值路线最后交叉验证一次。这里有个理论前提值得单独说一句——对因果稳定系统冲激响应的傅里叶变换就等于系统函数在 sjω 上的取值也就是 H(jω)H(s)|_{sjω}。这个等式成立的条件是收敛域包含虚轴等价于所有极点都在左半平面。所以你完全可以跳过fourier(h(t))直接做变量替换syms s w real Hs (2.8*s 1.99) / ((s^2 0.2*s 0.03) * (s 1) * (s 2)); Hw subs(Hs, s, 1j*w); % 稳定系统的频响 系统函数在 sjw 处取值 mag abs(Hw); ph angle(Hw) * 180/pi; fplot(mag, [-6 6]); % 幅频响应这段代码的关键在第二行subs把符号变量s整体替换成1j*w不做积分也不做变换速度快且不会引入dirac项。[-6 6]是频率扫描范围一般取谐振峰位置的三到五倍就够看清形状。如果系统不稳定有右半平面极点这个替换在数学上不成立必须先说明前提再画图否则结论是错的。反过来涉及阶跃响应、任意输入响应这类需要数值卷积的场合符号法就不好使了得切回lsim或step。判断标准很简单要公式用符号要曲线上的具体数值用数值函数。2.3 系统函数的参数取值分母系数、极点位置与稳定性的自洽性检查系统函数一旦写定后面所有结论都挂在它身上所以参数必须在动手画图之前自查一遍。原始设计给的是 H(s)(2.8s1.99)/((s²0.2s0.01)(s1)(s2))。这里有个必须提前发现的问题二次因子 s²0.2s0.01 的判别式是 0.2²−4×0.010恰好为零也就是说它给的是 s−0.1 的二重实极点而不是题目要求的「一对靠近虚轴的共轭复数极点」。想拿到共轭复极点常数项必须大于 0.01取 0.03 时判别式为 0.04−0.12−0.08极点变成 −0.1±j0.1414实部仍是 −0.1贴着虚轴正好满足靠近虚轴的要求。第二个问题是分母系数别手抄。把因式乘开和硬编码一串小数这两件事很容易对不上。让 Matlab 自己展开c 0.03; % 二次因子常数项0.01 才能得到共轭复极点 a2 [1 0.2 c]; % s^2 0.2s c a conv(a2, conv([1 1], [1 2])); % 再乘 (s1)(s2)得到 5 个分母系数 b [2.8 1.99]; % 分子 2.8s 1.99 p roots(a); % 极点 z roots(b); % 零点 max(real(p)) % 小于 0 才算稳定这是硬判据conv做多项式乘法roots求根这两行是整个实验的地基。max(real(p))0就是稳定性判据的直接实现所有极点实部为负系统稳定只要有一个实部大于零或者虚轴上出现二重及以上极点就不稳定。这行判断建议直接打印出来贴在报告里比任何文字描述都有说服力。顺便把几个后面要反复用的量算出来直流增益 H(0)1.99/0.06≈33.2说明这是个低频增益很高的系统零点在 s−0.7107是左半平面实零点极点分别是 −0.1±j0.1414、−1、−2。极点对离虚轴只有 0.1阻尼比约 0.58在 ω≈0.1414 附近会顶出一个谐振峰。这个峰值位置是后面判断参数微扰有没有改变系统本质的参照点必须先记下来。3. 总体设计三种实现结构与参数敏感性的对照方案3.1 直接型、级联型、并联型的系数转换链路同一个系统函数有三种实现结构数学上完全等价工程上差别很大。实验要求把三种都搭出来做对照转换靠三个函数完成。tf2sos把传递函数拆成二阶基本节的级联返回一个 n×6 的矩阵和一个总增益 Gsos2tf是它的逆运算tf2zp求零极点增益形式zp2tf反变换回去。这条链路走通了三种结构的系数就都到手了[sos, G] tf2sos(b, a); % 直接型 - 级联型sos 每行是 [b0 b1 b2 a0 a1 a2] [b2, a2] sos2tf(sos, G); % 级联型 - 直接型用于验证转换无损 [z, p, k] tf2zp(b, a); % 直接型 - 零极点增益型 [b3, a3] zp2tf(z, p, k); % 零极点增益型 - 直接型并联结构的基础sos矩阵的每一行对应一个二阶节前三个数是分子系数后三个是分母系数。G是各节之间的总增益画信号流图时必须把它标在输入端否则级联结构的幅频响应会整体差一个常数倍——这是报告里常见的扣分点图看着对数值对不上。b2/a2和原始的b/a理论上数值相等可以直接用max(abs(b-b2))验证量级应该在 1e-12 以下如果差得明显说明转换过程里丢了增益或者阶数没配平。zp2tf得到的b3/a3用于构造并联结构把每个极点对应的部分分式展开出来各支路独立并联输出相加。三种结构在信号流图上的区别画图时体现为直接型是一条长链所有系数耦合在一起级联型是若干二阶节串起来节与节之间只有信号传递并联型是若干支路并排各自独立再求和。图形本身不难难的是在图上把反馈系数标对——直接型的反馈支路跨越整个链条级联型每个二阶节内部自己闭环这是后面解释敏感性的全部依据。3.2 零极点分布到幅频响应的映射规则实验报告里最容易被追问的一句是为什么这个滤波器长这样。答案全在零极点相对虚轴的位置上有几条定量的经验规则可以直接用。极点越靠近虚轴幅频响应在该极点虚部对应的频率处峰值越高、越尖相频在该处下降越快零点越靠近虚轴幅频在该零点虚部对应的频率处凹陷越深相频在该处上升越快。零极点彼此靠近或者关于虚轴对称成对出现对幅频的影响会互相抵消曲线趋于平坦零极点都远离虚轴整条幅频曲线也平坦。用这套规则去读 H(s)(2.8s1.99)/((s²0.2s0.03)(s1)(s2))极点对 −0.1±j0.1414 离虚轴 0.1会在 ω≈0.14 处顶一个峰极点 −1 和 −2 离得远只贡献缓慢的滚降零点在 −0.71比极点对离虚轴远得多压不住这个峰。所以整体是低频高增益 低频段一个谐振峰 高频衰减的形状本质上是个带谐振的低通。这套规则也是设计滤波器时的正推工具想让某个频率分量过不去就在该频率对应的 s±jω 处放一对零点想让某个分量突出就在 s−σ±jωσ 小处放一对极点。第 4.3 节会具体用一次。3.3 微扰实验的设计改哪个系数、改多少、看什么指标微扰不是随便改一个数看看图变没变得设计成可控的对照实验否则结论写不出因果。改哪个系数直接型改分母多项式的高次项或常数项系数因为直接型的极点位置是高次方程的根任何一个系数的扰动都会传播到所有极点级联型改某一个二阶节的 a1 或 a2理论上只影响那一对极点并联型改某一条支路的极点系数只影响该支路。改多少题目建议在万分之一到 5% 之间。我的经验是分两档做一档 0.1%模拟有限字长下的量化误差一档 2%模拟元件容差两档放在一起才能看出小扰动是否已经导致结构失效。只做一档结论没有梯度。看什么指标不要只看图变没变要量化。至少记录四个量——极点最大位移 max|Δp|、谐振峰频率偏移 Δω、阶跃响应超调量变化、稳态值变化。做一个记录表把三种结构在同样扰动下的四项指标填进去结构扰动系数max|Δp|谐振峰偏移超调量变化稳态值变化直接型a4 ±0.1%待测待测待测待测级联型第二节 a2 ±0.1%待测待测待测待测并联型第一支路 a1 ±0.1%待测待测待测待测表里的待测由你自己的脚本填不要抄别人的数。填完之后那条结论自然就出来了扰动量相同时直接型的极点位移通常最大、谐振峰偏移最明显级联和并联的响应基本锁在被扰动的那一对零极点附近。原因是直接型的高阶多项式根对系数极其敏感系数只在第六位小数上动一下根可能已经跑出可观的距离级联和并联把高阶问题拆成了若干二阶问题每个二阶节的根只依赖自己的三个系数扰动被局部化了。这句话写进报告时最好配一张四宫格的pzmap对照图未扰动、直接型扰动后、级联型扰动后、并联型扰动后。四张图的坐标范围必须一致否则视觉上会骗人。4. 关键模块实现时域响应、截断频谱与滤波器落地4.1 冲激响应与阶跃响应的双轨求解冲激响应是系统函数反变换的直接结果用符号法拿表达式用数值法拿曲线两边对一下。syms s t Hs (2.8*s 1.99) / ((s^2 0.2*s 0.03) * (s 1) * (s 2)); ht ilaplace(Hs); % 冲激响应表达式默认 t 0 fplot(ht, [0 60]); % 只画 t 0 的部分别乘 heavisideilaplace返回的表达式本身已经隐含了因果性不需要再乘heaviside(t)。乘了之后符号引擎在 t0 处会遇到跳变画出来的图可能只剩半截或者出现一条竖线属于自己给自己找麻烦。[0 60]这个区间是按时间常数选的最靠近虚轴的极点实部是 −0.1时间常数 10 秒画 6 个时间常数足够看到衰减到稳态。阶跃响应同理符号法用ilaplace(Hs/s)数值法用stepb [2.8 1.99]; a [1 3.2 2.63 0.49 0.06]; % 由 conv 展开得到与上面 Hs 的分母一致 sys tf(b, a); step(sys, 60); % 数值阶跃响应与符号结果对照a这五个系数必须和符号表达式里乘开的结果一致[1 3.2 2.63 0.49 0.06]就是 (s²0.2s0.03)(s1)(s2) 展开的系数。数值法和符号法画在同一条曲线上如果形状差得远八成是系数对不上先回去查第 2.3 节。阶跃响应的稳态值可以用终值定理快速核对lim(t→∞)y(t)lim(s→0)s·(Hs/s)H(0)33.2。图上最后稳定在 33 附近说明整条链路没问题。这个核对花十秒钟能省掉后面半小时的排查。4.2 截断长度对频谱的影响与频域卷积定理的验证拓展部分的核心是让人亲眼看到截断就是乘矩形窗。理论频谱是四根冲激ω±1、±2一旦截断频谱就变成理论谱和窗谱的卷积除以 2π冲激被抹开成 Sa 形状这就是泄漏。先构造信号并做数值积分近似r 0.02; % 时间步长对应采样率 50 Hz够覆盖到 2 rad/s T 10; % 截取时长秒 t -T:r:0; % 截取区间长度 T f cos(t) cos(2*t); % 待分析信号 N 200; W 2*pi; % 频率扫描半宽 k -N:N; w k * W / N; % 频率轴共 2N1 个点 F r * f * exp(-1j*w.*t).; % 数值积分1×L 乘 L×(2N1) mag abs(F); ph angle(F) * 180/pi;第四行的矩阵维度是这个脚本最容易翻车的地方。f是 1×L 的行向量exp(-1j*w.*t)是 (2N1)×L 的矩阵因为w.是列向量、t是行向量转置之后变成 L×(2N1)再和f相乘才得到 1×(2N1) 的频谱。如果写成f*exp(-1j*t.*w)维度对不上直接报错如果是老版本容忍了某种广播得到的数也是错的。r是积分权重等价于 dt少了它幅度会差 50 倍。截取长度做三组T5、T10、T20。注意最低频率分量是 ω1周期 2π≈6.28 秒T5 已经小于一个周期属于不满足约束的截取——正好留着当反面教材用。三组结果的差别会非常直观T20 时窗谱主瓣宽度约 4π/T≈0.63ω1 和 ω2 两根谱线分得很开峰值位置准T5 时主瓣宽度约 2.5两根谱线糊成一团峰值位置还会偏移相位谱上出现明显的线性斜坡。这个斜坡来自截取区间不以原点对称相当于给信号乘了一个 e^{-jωT/2} 的相位因子幅度谱不受影响相位谱整体倾斜。做报告时把这一点写清楚比只贴三张图强得多。对称截取比如 −10 到 10相位谱就没有这个斜坡可以拿来做对照证明斜坡确实来自时间平移而不是计算错误。4.3 陷波型滤波器的零极点设计与卷积求输出滤波器设计部分给的系统函数是 H(s)(s²−0.2s1.97)/(s²0.2s1.01)。先把零极点算准这里非常容易记反。b [1 -0.2 1.97]; a [1 0.2 1.01]; z roots(b); % 零点0.1 ± j1.4在右半平面 p roots(a); % 极点-0.1 ± j1.0在左半平面零点实部是 0.1落在右半平面说明这是个非最小相位系统极点实部 −0.1贴着虚轴。零点虚部 1.4、极点虚部 1.0两点在频率轴上错开形成的形状是ω≈1.0 处因极点靠近虚轴而隆起ω≈1.4 处因零点靠近虚轴而凹陷ω 很大时分子分母同阶幅频趋于 1。把这条幅频算出来看数w 0:0.005:4; H freqs(b, a, w); mag abs(H); ph angle(H) * 180/pi; H1 polyval(b, 1j*1) / polyval(a, 1j*1); % ω1 处的复增益 H2 polyval(b, 1j*2) / polyval(a, 1j*2); % ω2 处的复增益 abs(H1), abs(H2)freqs直接给连续系统的频率响应输入是分子分母系数和频率向量。polyval那两行是手算校验把 sjω 代进去求复数值abs拿幅度、angle拿相位。跑出来 |H(j1)|≈4.95|H(j2)|≈0.69比值约 7.2 倍。也就是说输入 f(t)cos(t)cos(2t) 经过这个系统之后ω1 的分量被放大近 5 倍ω2 的分量被压到 0.69两个分量的幅度比从 1:1 变成约 7:1视觉上基本只剩低频那一个。用lsim拿输出波形sys tf(b, a); t 0:0.01:60; f cos(t) cos(2*t); y lsim(sys, f, t); idx t 40; % 跳过瞬态只取稳态段 y steady y(idx);lsim对任意输入求响应本质是数值卷积比符号法求ilaplace(Hs*Fs)稳得多。t 40这个截断点是按瞬态时长估的极点实部 −0.1时间常数 10 秒四倍时间常数后瞬态衰减到 2% 以下可以认为进入稳态。稳态段的 FFT 应该只在 ω1 和 ω2 附近有峰且 ω1 的峰高明显大于 ω2这就把滤波效果从定性描述变成了可测量的数。如果想让 ω2 被压得更狠把零点挪到 ±j2 就行此时 H(s)(s²4)/((s0.5)²1) 这一类形式会在 ω2 处形成真正的零点理论上完全滤除该分量。这是报告里可以加的一段改进讨论能体现设计意识而不只是照抄参数。4.4 音频信号的采集、回放与频谱分析应用部分把前面所有连续时间的直觉换成离散处理坑最多。rec audiorecorder(8000, 16, 1); % 8 kHz 采样、16 bit 量化、单声道 recordblocking(rec, 3); % 录 3 秒阻塞式录完才继续 y getaudiodata(rec); % 取回 [-1,1] 区间的双精度序列 audiowrite(speech.wav, y, 8000); % 存盘采样率必须和录制时一致 [y, fs] audioread(speech.wav); % 读回fs 应为 8000 sound(y, fs); % 回放 N 2^nextpow2(length(y)); % 取 2 的幂次FFT 最快 Y fft(y, N) / length(y); % 归一化幅度才有物理意义 fax (0:N-1) * fs / N; % 频率轴单位 Hz plot(fax(1:N/2), 20*log10(abs(Y(1:N/2)) eps));audiorecorder三个参数分别是采样率、位深、通道数8 kHz 对语音够用奈奎斯特上限 4 kHz覆盖语音主要能量但对音乐明显不够会听出闷。recordblocking是阻塞调用程序会停在那里等 3 秒别以为是卡死了。fft之后除以length(y)是归一化不做这一步幅度值随信号长度变化换一段录音图就完全不一样没法比较。频率轴用(0:N-1)*fs/N生成只画前一半是因为实信号的频谱共轭对称后一半是冗余的。20*log10换成 dB 是为了让弱分量也看得见加eps防止对零取对数出-Inf。看频谱时关注三件事基频在哪语音通常在 100~300 Hz、共振峰在哪元音的前三个峰大致在 500 Hz、1500 Hz、2500 Hz 附近、高频能量衰减得多快。这三点写进报告比画出了频谱图有意义得多。如果频谱上出现明显的高频折叠或者一片糊先查采样率设置和 FFT 点数别怀疑信号本身。5. 避坑与排查这套实验里最容易翻车的五个地方5.1 分母系数和因式形式对不上稳定性结论直接失效现象roots(a)算出来的极点和手算的因式根不一致pzmap上有极点跑到右半平面但系统明明是稳定的。原因分母系数是手工从因式乘开抄下来的抄错一位或者用了一个和因式不匹配的旧数组。多项式系数对根的影响是全局的常数项错 0.01根可能移动 0.1。解决一律用conv展开不用手抄a_check conv([1 0.2 0.03], conv([1 1], [1 2])); % 应得 [1 3.2 2.63 0.49 0.06] max(abs(a - a_check)) % 差值应在 1e-12 量级conv嵌套两次就能处理四阶再多就换成poly配合根向量。这个自检花两行代码是整个实验最划算的保险。5.2 二次因子判别式为零拿到的不是共轭复极点现象题目要求一对靠近虚轴的共轭复数极点但pzmap上显示的是实轴上的两个重合点幅频响应也没有明显的谐振峰。原因s²0.2s0.01 的判别式恰好等于零两个根都是 −0.1退化成二重实极点。这种情况下系统仍然稳定但和题目要求的共轭复极点不是一回事后面解释谐振现象时会对不上。解决把常数项往上调让判别式小于零。取 0.03 得到 −0.1±j0.1414实部仍是 −0.1满足靠近虚轴取 0.0101 得到的虚部太小峰不明显不好看。改完记得同步更新分母系数并重跑稳定判据。报告里把这次修正写进去反而是个加分项说明你验算过。5.3 相频响应画出一条直线或者报错问题出在符号变量没定义现象angle(Fw)报 Undefined function or variable Fw或者对含heaviside的表达式求fourier后画出来是一条奇怪的水平线。原因前一段代码里变量名是Hw这一段写成了Fw或者fourier处理含阶跃函数的符号表达式时返回值不含你预期的频率信息。解决变量名统一或者在每段开头clear后重新定义更重要的是对稳定因果系统跳过fourier直接用subs(Hs, s, 1j*w)做变量替换。这样得到的表达式干净abs和angle都能正常作用。如果一定要用fourier先确认表达式里没有heaviside或者把冲激响应用数值法算出来再做fft。5.4 截断频谱的幅度差了 50 倍或者矩阵维度报错现象频谱图的纵轴数值明显不合理或者执行到F r*f*exp(-j*t*w)时报 Matrix dimensions must agree。原因前者是漏了积分步长r把数值积分当成了求和后者是行向量和列向量的方向没对齐t和w的位置决定了结果矩阵的形状。解决固定写法F r * f * exp(-1j*w.*t).;其中f是行向量w.是列向量t是行向量中间那个矩阵是 (2N1)×L转置后和f相乘得到 1×(2N1)。如果实在记不住直接换成fft加fftshift再手动构造频率轴虽然是对离散序列做变换但至少维度不会错只要注意补零和分辨率的关系就行。5.5 滤波器零极点记反整段滤波结论站不住现象报告里写零点 0.1±j、极点 −0.1±1.4j但画出来的幅频曲线在 ω1.4 处是峰而不是谷和文字描述矛盾。原因直接套用了别人的结论没有自己roots一遍。实际算出来零点是 0.1±j1.4极点是 −0.1±j1.0虚部数值刚好对调了结论就全反了。解决任何零极点结论都从roots(b)和roots(a)出发先算再说。另外注意零点实部是正的位于右半平面这意味着系统是非最小相位的——幅频响应和把零点镜像到左半平面的最小相位系统完全一样但相位响应不同。这一点如果能在报告里指出来说明你真的看懂了零极点图而不只是画了个图。6. 验收与答辩把图、数、结论对上号的检查方法实验报告交上去之前我习惯做一遍三对齐检查图上的特征点、代码里的数、结论里的话三者必须能一一对应。这套检查不花什么时间但能挡掉大部分追问。第一步查系统函数。把roots(a)的输出、pzmap上点的位置、报告里写的极点坐标放在一起三处必须完全一致。尤其注意共轭复极点的虚部0.1414 和 0.14 在图上分不出来在数值上差 1%如果后面要讨论谐振峰频率这个差别会被放大。第二步查稳定性。max(real(p))的符号、文字结论里的稳定/不稳定、阶跃响应是否收敛到有限值这三者必须相互印证。三者只要有一个对不上先停下来查别往下写。第三步查滤波效果。|H(j1)| 和 |H(j2)| 的数值、输出波形里两个分量的幅度比、频谱图上两个峰的高度比这三者应该大致吻合。输出波形的幅度比可以这样估稳态段里找相邻的极大值低频分量周期 2π、高频分量周期 π包络的起伏深度直接反映比值。如果算出来是 7 倍但图上看着只有 2 倍八成是稳态段截早了瞬态还没衰减完。答辩时高频出现的几个问题答案其实都在前面几节里。问为什么级联结构对参数不敏感答二阶基本节的根只依赖本节三个系数扰动被局部化问截断为什么会让谱线变宽答时域乘矩形窗等于频域卷积窗谱主瓣宽度和截取长度成反比问这个滤波器为什么在 ω1.4 处凹陷答零点 0.1±j1.4 靠近虚轴距离只有 0.1形成谷值。这三个问题答得出来基本就稳了。还有一个容易被忽略的点所有图的坐标范围和坐标轴标签。四宫格对照图的横纵坐标必须统一否则视觉对比毫无意义频率轴要标单位rad/s 还是 Hz时间轴要标秒dB 和线性幅度不能混着画。这些细节不加分但错了会扣分。我把这套检查做成一张表每次交报告前逐项打勾检查项判据通过标准分母系数conv展开 vs 硬编码差值 1e-12极点位置roots(a)vs 报告文字完全一致稳定性max(real(p))vs 结论符号一致直流增益H(0) vs 阶跃稳态值误差 1%频响峰值freqs扫描 vs 理论估算频率偏差 5%滤波比|H(j1)|/|H(j2)| vs 波形幅度比同量级坐标轴所有对照图范围一致、有单位这张表最初是因为一次翻车才建起来的当时四张pzmap的坐标范围各不相同直接型扰动后的极点看着几乎没动实际上已经跑出画面外了结论写反了还没察觉。从那以后我每次画对照图都强制走一遍统一坐标范围的流程宁可多写两行xlim/ylim也不相信自己的眼睛。希望这套流程帮到你。本文还有配套的精品资源点击获取 
长河君 · 2026/10/28 5:32:13 · 阅读 1.2k 连续时间LTI系统Matlab实验：零极点、稳定性与滤波器设计
简介针对信号与线性系统综合实验课程这份报告类资源面向通信工程、电子信息工程等专业已完成基础知识学习的本科生与研究生。内容围绕连续时间LTI系统展开以设置零极点的方式设计稳定的4至6阶系统分别用直接、级联、并联三种结构建流图并借助Matlab仿真时域冲激响应、阶跃响应与幅频、相频特性通过微小改变系数对比不同结构在零极点分布、稳定性与响应上的敏感性差异同时构建周期信号分析频谱设计LTI滤波器完成滤波验证并延伸到音频信号的获取、输出与频谱分析。资源包为单个docx文档约2.62MB含实验目的、原理推导、步骤、完整Matlab代码与结果分析可直接作为课程实验报告模板或自主学习参考。目前已有1303人学习适合需要对照代码复现仿真、理清零极点与频率特性关系、补齐频域卷积定理等知识点的读者。1. 从 s 平面上的一对极点说起这套连续时间 LTI 系统综合实验要交付什么通信、电子信息这类专业的综合实验里连续时间 LTI 系统的时域频域特性分析与滤波器设计几乎绕不开。题目看着温和给一个 4~6 阶系统在 s 平面上摆几对零极点用 Matlab 把冲激响应、阶跃响应、幅频相频画出来再设计一个滤波器把 cos(t)cos(2t) 里的某个分量压下去。真正动手才发现极点实部差 0.01谐振峰就能挪出一大截分母系数手抄错一位稳定性结论直接反过来零极点位置记反整段滤波结论全部作废。这套方案面向的是要交一份完整综合实验报告的人——从系统函数的构造、三种实现结构的等价性与敏感性对比到截断频谱的误差解释再到陷波型滤波器的零极点设计和音频信号的采集分析一条线走完每一段都留下能复现的代码和能写进结论的判据而不是画几张图就交差。2. 需求拆解与工具选型符号计算与数值仿真两条腿怎么分工2.1 三层任务结构从系统分析到频谱验证再到音频应用整套实验其实是三层递进的交付物每一层的产物形式不一样验收口径也不一样。第一层是基础部分设计一个 4~6 阶的稳定连续系统其中至少含一对靠近虚轴的共轭复数极点用直接型、级联型、并联型三种结构分别实现画出信号流图对每种结构做时域和频域仿真再做参数微扰比较三种结构对系数变化的敏感程度。这一层交付的是「一组对照图 一条敏感性结论」结论必须落到哪种结构在有限字长下更稳。第二层是拓展部分构造一个含多个频率分量的周期信号截取不同长度做频谱分析把截断后的谱和理论谱摆在一起解释差异然后用 s 平面零极点法设计滤波器滤掉其中一个分量验证滤波效果。这一层交付的是「误差解释 滤波前后对比」核心是把频域卷积定理讲通。第三层是应用设计用声卡录一段音频存成 wav读回来画波形和频谱。这一层交付的是「一段可回放的音频 一张能读出基频和共振峰的频谱图」顺带把采样定理落到实处。三层之间不是并列关系。第一层把系统函数、零极点、稳定性这条线打通第二层把这条线搬到信号侧用窗谱去解释泄漏第三层把连续时间的直觉换成离散的 FFT稍不留神就会在频率轴上翻车。写报告时如果按这三层组织章节读者也就是批改的人能顺着你的推导走而不是在一堆孤立截图里找逻辑。需要提醒的是三层任务的代码风格要统一都用脚本都清空工作区都显式给参数。中途换变量名、换符号变量声明方式最后拼报告时会非常痛苦。2.2 符号计算与数值函数的边界ilaplace/fourier 与 freqs/step 各管哪一段Matlab 里做信号与系统仿真有两条技术路线选错了会白折腾很久。符号路线用syms声明符号变量用ilaplace求拉普拉斯反变换拿冲激响应用fourier求傅里叶变换拿频谱用ezplot或fplot画图。优点是表达式干净能直接写进报告的公式栏缺点是含heaviside的符号表达式在求fourier时经常返回一堆带pi、dirac的东西画出来还不是你想要的样子。数值路线用tf、freqs、impulse、step、lsim输入是系数向量输出是数值数组。优点是稳定、快、不会出符号怪缺点是要自己管频率轴和时间轴采样点选错了图形会失真。我的做法是两条腿都留着推导和公式展示走符号路线图形和定量指标走数值路线最后交叉验证一次。这里有个理论前提值得单独说一句——对因果稳定系统冲激响应的傅里叶变换就等于系统函数在 sjω 上的取值也就是 H(jω)H(s)|_{sjω}。这个等式成立的条件是收敛域包含虚轴等价于所有极点都在左半平面。所以你完全可以跳过fourier(h(t))直接做变量替换syms s w real Hs (2.8*s 1.99) / ((s^2 0.2*s 0.03) * (s 1) * (s 2)); Hw subs(Hs, s, 1j*w); % 稳定系统的频响 系统函数在 sjw 处取值 mag abs(Hw); ph angle(Hw) * 180/pi; fplot(mag, [-6 6]); % 幅频响应这段代码的关键在第二行subs把符号变量s整体替换成1j*w不做积分也不做变换速度快且不会引入dirac项。[-6 6]是频率扫描范围一般取谐振峰位置的三到五倍就够看清形状。如果系统不稳定有右半平面极点这个替换在数学上不成立必须先说明前提再画图否则结论是错的。反过来涉及阶跃响应、任意输入响应这类需要数值卷积的场合符号法就不好使了得切回lsim或step。判断标准很简单要公式用符号要曲线上的具体数值用数值函数。2.3 系统函数的参数取值分母系数、极点位置与稳定性的自洽性检查系统函数一旦写定后面所有结论都挂在它身上所以参数必须在动手画图之前自查一遍。原始设计给的是 H(s)(2.8s1.99)/((s²0.2s0.01)(s1)(s2))。这里有个必须提前发现的问题二次因子 s²0.2s0.01 的判别式是 0.2²−4×0.010恰好为零也就是说它给的是 s−0.1 的二重实极点而不是题目要求的「一对靠近虚轴的共轭复数极点」。想拿到共轭复极点常数项必须大于 0.01取 0.03 时判别式为 0.04−0.12−0.08极点变成 −0.1±j0.1414实部仍是 −0.1贴着虚轴正好满足靠近虚轴的要求。第二个问题是分母系数别手抄。把因式乘开和硬编码一串小数这两件事很容易对不上。让 Matlab 自己展开c 0.03; % 二次因子常数项0.01 才能得到共轭复极点 a2 [1 0.2 c]; % s^2 0.2s c a conv(a2, conv([1 1], [1 2])); % 再乘 (s1)(s2)得到 5 个分母系数 b [2.8 1.99]; % 分子 2.8s 1.99 p roots(a); % 极点 z roots(b); % 零点 max(real(p)) % 小于 0 才算稳定这是硬判据conv做多项式乘法roots求根这两行是整个实验的地基。max(real(p))0就是稳定性判据的直接实现所有极点实部为负系统稳定只要有一个实部大于零或者虚轴上出现二重及以上极点就不稳定。这行判断建议直接打印出来贴在报告里比任何文字描述都有说服力。顺便把几个后面要反复用的量算出来直流增益 H(0)1.99/0.06≈33.2说明这是个低频增益很高的系统零点在 s−0.7107是左半平面实零点极点分别是 −0.1±j0.1414、−1、−2。极点对离虚轴只有 0.1阻尼比约 0.58在 ω≈0.1414 附近会顶出一个谐振峰。这个峰值位置是后面判断参数微扰有没有改变系统本质的参照点必须先记下来。3. 总体设计三种实现结构与参数敏感性的对照方案3.1 直接型、级联型、并联型的系数转换链路同一个系统函数有三种实现结构数学上完全等价工程上差别很大。实验要求把三种都搭出来做对照转换靠三个函数完成。tf2sos把传递函数拆成二阶基本节的级联返回一个 n×6 的矩阵和一个总增益 Gsos2tf是它的逆运算tf2zp求零极点增益形式zp2tf反变换回去。这条链路走通了三种结构的系数就都到手了[sos, G] tf2sos(b, a); % 直接型 - 级联型sos 每行是 [b0 b1 b2 a0 a1 a2] [b2, a2] sos2tf(sos, G); % 级联型 - 直接型用于验证转换无损 [z, p, k] tf2zp(b, a); % 直接型 - 零极点增益型 [b3, a3] zp2tf(z, p, k); % 零极点增益型 - 直接型并联结构的基础sos矩阵的每一行对应一个二阶节前三个数是分子系数后三个是分母系数。G是各节之间的总增益画信号流图时必须把它标在输入端否则级联结构的幅频响应会整体差一个常数倍——这是报告里常见的扣分点图看着对数值对不上。b2/a2和原始的b/a理论上数值相等可以直接用max(abs(b-b2))验证量级应该在 1e-12 以下如果差得明显说明转换过程里丢了增益或者阶数没配平。zp2tf得到的b3/a3用于构造并联结构把每个极点对应的部分分式展开出来各支路独立并联输出相加。三种结构在信号流图上的区别画图时体现为直接型是一条长链所有系数耦合在一起级联型是若干二阶节串起来节与节之间只有信号传递并联型是若干支路并排各自独立再求和。图形本身不难难的是在图上把反馈系数标对——直接型的反馈支路跨越整个链条级联型每个二阶节内部自己闭环这是后面解释敏感性的全部依据。3.2 零极点分布到幅频响应的映射规则实验报告里最容易被追问的一句是为什么这个滤波器长这样。答案全在零极点相对虚轴的位置上有几条定量的经验规则可以直接用。极点越靠近虚轴幅频响应在该极点虚部对应的频率处峰值越高、越尖相频在该处下降越快零点越靠近虚轴幅频在该零点虚部对应的频率处凹陷越深相频在该处上升越快。零极点彼此靠近或者关于虚轴对称成对出现对幅频的影响会互相抵消曲线趋于平坦零极点都远离虚轴整条幅频曲线也平坦。用这套规则去读 H(s)(2.8s1.99)/((s²0.2s0.03)(s1)(s2))极点对 −0.1±j0.1414 离虚轴 0.1会在 ω≈0.14 处顶一个峰极点 −1 和 −2 离得远只贡献缓慢的滚降零点在 −0.71比极点对离虚轴远得多压不住这个峰。所以整体是低频高增益 低频段一个谐振峰 高频衰减的形状本质上是个带谐振的低通。这套规则也是设计滤波器时的正推工具想让某个频率分量过不去就在该频率对应的 s±jω 处放一对零点想让某个分量突出就在 s−σ±jωσ 小处放一对极点。第 4.3 节会具体用一次。3.3 微扰实验的设计改哪个系数、改多少、看什么指标微扰不是随便改一个数看看图变没变得设计成可控的对照实验否则结论写不出因果。改哪个系数直接型改分母多项式的高次项或常数项系数因为直接型的极点位置是高次方程的根任何一个系数的扰动都会传播到所有极点级联型改某一个二阶节的 a1 或 a2理论上只影响那一对极点并联型改某一条支路的极点系数只影响该支路。改多少题目建议在万分之一到 5% 之间。我的经验是分两档做一档 0.1%模拟有限字长下的量化误差一档 2%模拟元件容差两档放在一起才能看出小扰动是否已经导致结构失效。只做一档结论没有梯度。看什么指标不要只看图变没变要量化。至少记录四个量——极点最大位移 max|Δp|、谐振峰频率偏移 Δω、阶跃响应超调量变化、稳态值变化。做一个记录表把三种结构在同样扰动下的四项指标填进去结构扰动系数max|Δp|谐振峰偏移超调量变化稳态值变化直接型a4 ±0.1%待测待测待测待测级联型第二节 a2 ±0.1%待测待测待测待测并联型第一支路 a1 ±0.1%待测待测待测待测表里的待测由你自己的脚本填不要抄别人的数。填完之后那条结论自然就出来了扰动量相同时直接型的极点位移通常最大、谐振峰偏移最明显级联和并联的响应基本锁在被扰动的那一对零极点附近。原因是直接型的高阶多项式根对系数极其敏感系数只在第六位小数上动一下根可能已经跑出可观的距离级联和并联把高阶问题拆成了若干二阶问题每个二阶节的根只依赖自己的三个系数扰动被局部化了。这句话写进报告时最好配一张四宫格的pzmap对照图未扰动、直接型扰动后、级联型扰动后、并联型扰动后。四张图的坐标范围必须一致否则视觉上会骗人。4. 关键模块实现时域响应、截断频谱与滤波器落地4.1 冲激响应与阶跃响应的双轨求解冲激响应是系统函数反变换的直接结果用符号法拿表达式用数值法拿曲线两边对一下。syms s t Hs (2.8*s 1.99) / ((s^2 0.2*s 0.03) * (s 1) * (s 2)); ht ilaplace(Hs); % 冲激响应表达式默认 t 0 fplot(ht, [0 60]); % 只画 t 0 的部分别乘 heavisideilaplace返回的表达式本身已经隐含了因果性不需要再乘heaviside(t)。乘了之后符号引擎在 t0 处会遇到跳变画出来的图可能只剩半截或者出现一条竖线属于自己给自己找麻烦。[0 60]这个区间是按时间常数选的最靠近虚轴的极点实部是 −0.1时间常数 10 秒画 6 个时间常数足够看到衰减到稳态。阶跃响应同理符号法用ilaplace(Hs/s)数值法用stepb [2.8 1.99]; a [1 3.2 2.63 0.49 0.06]; % 由 conv 展开得到与上面 Hs 的分母一致 sys tf(b, a); step(sys, 60); % 数值阶跃响应与符号结果对照a这五个系数必须和符号表达式里乘开的结果一致[1 3.2 2.63 0.49 0.06]就是 (s²0.2s0.03)(s1)(s2) 展开的系数。数值法和符号法画在同一条曲线上如果形状差得远八成是系数对不上先回去查第 2.3 节。阶跃响应的稳态值可以用终值定理快速核对lim(t→∞)y(t)lim(s→0)s·(Hs/s)H(0)33.2。图上最后稳定在 33 附近说明整条链路没问题。这个核对花十秒钟能省掉后面半小时的排查。4.2 截断长度对频谱的影响与频域卷积定理的验证拓展部分的核心是让人亲眼看到截断就是乘矩形窗。理论频谱是四根冲激ω±1、±2一旦截断频谱就变成理论谱和窗谱的卷积除以 2π冲激被抹开成 Sa 形状这就是泄漏。先构造信号并做数值积分近似r 0.02; % 时间步长对应采样率 50 Hz够覆盖到 2 rad/s T 10; % 截取时长秒 t -T:r:0; % 截取区间长度 T f cos(t) cos(2*t); % 待分析信号 N 200; W 2*pi; % 频率扫描半宽 k -N:N; w k * W / N; % 频率轴共 2N1 个点 F r * f * exp(-1j*w.*t).; % 数值积分1×L 乘 L×(2N1) mag abs(F); ph angle(F) * 180/pi;第四行的矩阵维度是这个脚本最容易翻车的地方。f是 1×L 的行向量exp(-1j*w.*t)是 (2N1)×L 的矩阵因为w.是列向量、t是行向量转置之后变成 L×(2N1)再和f相乘才得到 1×(2N1) 的频谱。如果写成f*exp(-1j*t.*w)维度对不上直接报错如果是老版本容忍了某种广播得到的数也是错的。r是积分权重等价于 dt少了它幅度会差 50 倍。截取长度做三组T5、T10、T20。注意最低频率分量是 ω1周期 2π≈6.28 秒T5 已经小于一个周期属于不满足约束的截取——正好留着当反面教材用。三组结果的差别会非常直观T20 时窗谱主瓣宽度约 4π/T≈0.63ω1 和 ω2 两根谱线分得很开峰值位置准T5 时主瓣宽度约 2.5两根谱线糊成一团峰值位置还会偏移相位谱上出现明显的线性斜坡。这个斜坡来自截取区间不以原点对称相当于给信号乘了一个 e^{-jωT/2} 的相位因子幅度谱不受影响相位谱整体倾斜。做报告时把这一点写清楚比只贴三张图强得多。对称截取比如 −10 到 10相位谱就没有这个斜坡可以拿来做对照证明斜坡确实来自时间平移而不是计算错误。4.3 陷波型滤波器的零极点设计与卷积求输出滤波器设计部分给的系统函数是 H(s)(s²−0.2s1.97)/(s²0.2s1.01)。先把零极点算准这里非常容易记反。b [1 -0.2 1.97]; a [1 0.2 1.01]; z roots(b); % 零点0.1 ± j1.4在右半平面 p roots(a); % 极点-0.1 ± j1.0在左半平面零点实部是 0.1落在右半平面说明这是个非最小相位系统极点实部 −0.1贴着虚轴。零点虚部 1.4、极点虚部 1.0两点在频率轴上错开形成的形状是ω≈1.0 处因极点靠近虚轴而隆起ω≈1.4 处因零点靠近虚轴而凹陷ω 很大时分子分母同阶幅频趋于 1。把这条幅频算出来看数w 0:0.005:4; H freqs(b, a, w); mag abs(H); ph angle(H) * 180/pi; H1 polyval(b, 1j*1) / polyval(a, 1j*1); % ω1 处的复增益 H2 polyval(b, 1j*2) / polyval(a, 1j*2); % ω2 处的复增益 abs(H1), abs(H2)freqs直接给连续系统的频率响应输入是分子分母系数和频率向量。polyval那两行是手算校验把 sjω 代进去求复数值abs拿幅度、angle拿相位。跑出来 |H(j1)|≈4.95|H(j2)|≈0.69比值约 7.2 倍。也就是说输入 f(t)cos(t)cos(2t) 经过这个系统之后ω1 的分量被放大近 5 倍ω2 的分量被压到 0.69两个分量的幅度比从 1:1 变成约 7:1视觉上基本只剩低频那一个。用lsim拿输出波形sys tf(b, a); t 0:0.01:60; f cos(t) cos(2*t); y lsim(sys, f, t); idx t 40; % 跳过瞬态只取稳态段 y steady y(idx);lsim对任意输入求响应本质是数值卷积比符号法求ilaplace(Hs*Fs)稳得多。t 40这个截断点是按瞬态时长估的极点实部 −0.1时间常数 10 秒四倍时间常数后瞬态衰减到 2% 以下可以认为进入稳态。稳态段的 FFT 应该只在 ω1 和 ω2 附近有峰且 ω1 的峰高明显大于 ω2这就把滤波效果从定性描述变成了可测量的数。如果想让 ω2 被压得更狠把零点挪到 ±j2 就行此时 H(s)(s²4)/((s0.5)²1) 这一类形式会在 ω2 处形成真正的零点理论上完全滤除该分量。这是报告里可以加的一段改进讨论能体现设计意识而不只是照抄参数。4.4 音频信号的采集、回放与频谱分析应用部分把前面所有连续时间的直觉换成离散处理坑最多。rec audiorecorder(8000, 16, 1); % 8 kHz 采样、16 bit 量化、单声道 recordblocking(rec, 3); % 录 3 秒阻塞式录完才继续 y getaudiodata(rec); % 取回 [-1,1] 区间的双精度序列 audiowrite(speech.wav, y, 8000); % 存盘采样率必须和录制时一致 [y, fs] audioread(speech.wav); % 读回fs 应为 8000 sound(y, fs); % 回放 N 2^nextpow2(length(y)); % 取 2 的幂次FFT 最快 Y fft(y, N) / length(y); % 归一化幅度才有物理意义 fax (0:N-1) * fs / N; % 频率轴单位 Hz plot(fax(1:N/2), 20*log10(abs(Y(1:N/2)) eps));audiorecorder三个参数分别是采样率、位深、通道数8 kHz 对语音够用奈奎斯特上限 4 kHz覆盖语音主要能量但对音乐明显不够会听出闷。recordblocking是阻塞调用程序会停在那里等 3 秒别以为是卡死了。fft之后除以length(y)是归一化不做这一步幅度值随信号长度变化换一段录音图就完全不一样没法比较。频率轴用(0:N-1)*fs/N生成只画前一半是因为实信号的频谱共轭对称后一半是冗余的。20*log10换成 dB 是为了让弱分量也看得见加eps防止对零取对数出-Inf。看频谱时关注三件事基频在哪语音通常在 100~300 Hz、共振峰在哪元音的前三个峰大致在 500 Hz、1500 Hz、2500 Hz 附近、高频能量衰减得多快。这三点写进报告比画出了频谱图有意义得多。如果频谱上出现明显的高频折叠或者一片糊先查采样率设置和 FFT 点数别怀疑信号本身。5. 避坑与排查这套实验里最容易翻车的五个地方5.1 分母系数和因式形式对不上稳定性结论直接失效现象roots(a)算出来的极点和手算的因式根不一致pzmap上有极点跑到右半平面但系统明明是稳定的。原因分母系数是手工从因式乘开抄下来的抄错一位或者用了一个和因式不匹配的旧数组。多项式系数对根的影响是全局的常数项错 0.01根可能移动 0.1。解决一律用conv展开不用手抄a_check conv([1 0.2 0.03], conv([1 1], [1 2])); % 应得 [1 3.2 2.63 0.49 0.06] max(abs(a - a_check)) % 差值应在 1e-12 量级conv嵌套两次就能处理四阶再多就换成poly配合根向量。这个自检花两行代码是整个实验最划算的保险。5.2 二次因子判别式为零拿到的不是共轭复极点现象题目要求一对靠近虚轴的共轭复数极点但pzmap上显示的是实轴上的两个重合点幅频响应也没有明显的谐振峰。原因s²0.2s0.01 的判别式恰好等于零两个根都是 −0.1退化成二重实极点。这种情况下系统仍然稳定但和题目要求的共轭复极点不是一回事后面解释谐振现象时会对不上。解决把常数项往上调让判别式小于零。取 0.03 得到 −0.1±j0.1414实部仍是 −0.1满足靠近虚轴取 0.0101 得到的虚部太小峰不明显不好看。改完记得同步更新分母系数并重跑稳定判据。报告里把这次修正写进去反而是个加分项说明你验算过。5.3 相频响应画出一条直线或者报错问题出在符号变量没定义现象angle(Fw)报 Undefined function or variable Fw或者对含heaviside的表达式求fourier后画出来是一条奇怪的水平线。原因前一段代码里变量名是Hw这一段写成了Fw或者fourier处理含阶跃函数的符号表达式时返回值不含你预期的频率信息。解决变量名统一或者在每段开头clear后重新定义更重要的是对稳定因果系统跳过fourier直接用subs(Hs, s, 1j*w)做变量替换。这样得到的表达式干净abs和angle都能正常作用。如果一定要用fourier先确认表达式里没有heaviside或者把冲激响应用数值法算出来再做fft。5.4 截断频谱的幅度差了 50 倍或者矩阵维度报错现象频谱图的纵轴数值明显不合理或者执行到F r*f*exp(-j*t*w)时报 Matrix dimensions must agree。原因前者是漏了积分步长r把数值积分当成了求和后者是行向量和列向量的方向没对齐t和w的位置决定了结果矩阵的形状。解决固定写法F r * f * exp(-1j*w.*t).;其中f是行向量w.是列向量t是行向量中间那个矩阵是 (2N1)×L转置后和f相乘得到 1×(2N1)。如果实在记不住直接换成fft加fftshift再手动构造频率轴虽然是对离散序列做变换但至少维度不会错只要注意补零和分辨率的关系就行。5.5 滤波器零极点记反整段滤波结论站不住现象报告里写零点 0.1±j、极点 −0.1±1.4j但画出来的幅频曲线在 ω1.4 处是峰而不是谷和文字描述矛盾。原因直接套用了别人的结论没有自己roots一遍。实际算出来零点是 0.1±j1.4极点是 −0.1±j1.0虚部数值刚好对调了结论就全反了。解决任何零极点结论都从roots(b)和roots(a)出发先算再说。另外注意零点实部是正的位于右半平面这意味着系统是非最小相位的——幅频响应和把零点镜像到左半平面的最小相位系统完全一样但相位响应不同。这一点如果能在报告里指出来说明你真的看懂了零极点图而不只是画了个图。6. 验收与答辩把图、数、结论对上号的检查方法实验报告交上去之前我习惯做一遍三对齐检查图上的特征点、代码里的数、结论里的话三者必须能一一对应。这套检查不花什么时间但能挡掉大部分追问。第一步查系统函数。把roots(a)的输出、pzmap上点的位置、报告里写的极点坐标放在一起三处必须完全一致。尤其注意共轭复极点的虚部0.1414 和 0.14 在图上分不出来在数值上差 1%如果后面要讨论谐振峰频率这个差别会被放大。第二步查稳定性。max(real(p))的符号、文字结论里的稳定/不稳定、阶跃响应是否收敛到有限值这三者必须相互印证。三者只要有一个对不上先停下来查别往下写。第三步查滤波效果。|H(j1)| 和 |H(j2)| 的数值、输出波形里两个分量的幅度比、频谱图上两个峰的高度比这三者应该大致吻合。输出波形的幅度比可以这样估稳态段里找相邻的极大值低频分量周期 2π、高频分量周期 π包络的起伏深度直接反映比值。如果算出来是 7 倍但图上看着只有 2 倍八成是稳态段截早了瞬态还没衰减完。答辩时高频出现的几个问题答案其实都在前面几节里。问为什么级联结构对参数不敏感答二阶基本节的根只依赖本节三个系数扰动被局部化问截断为什么会让谱线变宽答时域乘矩形窗等于频域卷积窗谱主瓣宽度和截取长度成反比问这个滤波器为什么在 ω1.4 处凹陷答零点 0.1±j1.4 靠近虚轴距离只有 0.1形成谷值。这三个问题答得出来基本就稳了。还有一个容易被忽略的点所有图的坐标范围和坐标轴标签。四宫格对照图的横纵坐标必须统一否则视觉对比毫无意义频率轴要标单位rad/s 还是 Hz时间轴要标秒dB 和线性幅度不能混着画。这些细节不加分但错了会扣分。我把这套检查做成一张表每次交报告前逐项打勾检查项判据通过标准分母系数conv展开 vs 硬编码差值 1e-12极点位置roots(a)vs 报告文字完全一致稳定性max(real(p))vs 结论符号一致直流增益H(0) vs 阶跃稳态值误差 1%频响峰值freqs扫描 vs 理论估算频率偏差 5%滤波比|H(j1)|/|H(j2)| vs 波形幅度比同量级坐标轴所有对照图范围一致、有单位这张表最初是因为一次翻车才建起来的当时四张pzmap的坐标范围各不相同直接型扰动后的极点看着几乎没动实际上已经跑出画面外了结论写反了还没察觉。从那以后我每次画对照图都强制走一遍统一坐标范围的流程宁可多写两行xlim/ylim也不相信自己的眼睛。希望这套流程帮到你。本文还有配套的精品资源点击获取