数字信号处理课程设计实战：MATLAB信号生成、频谱分析与FIR滤波器设计

简介中南大学数字信号处理课程设计任务书是一份面向数字信号处理课程学生的完整实践指导文档适用于需要完成课程设计、掌握MATLAB编程实现及DFT/滤波器设计的学习者。任务书系统梳理了三大内容连续信号采样与谱分析、周期方波信号滤波、音乐信号单/多重回声及混响效果设计并对采样定理、频率分辨率、混叠与截断效应、巴特沃斯与FIR/IIR滤波器原理等关键知识点展开论述同时还包含GUI界面设计要求与完整设计报告规范。下载包内共包含1个doc文件压缩包大小318KB便于直接查看与二次编辑。已有397人学习下载。任务书结构清晰涵盖课程设计目的、题目内容、设计要求、设计原理与思路、程序测试输出结果、总结及源代码清单能为学习者提供从理论分析到程序实现、再到报告撰写的全流程参考。1. 数字信号处理课程设计的转型从任务书到可验证的工程问题拿到中南大学数字信号处理课程设计任务书的第一反应多数人先翻最后几页看要交什么代码报告答辩演示实际上任务书的真正价值不在这些交付物而在它强制把《数字信号处理》教材里的公式串成一条完整的工程链路——信号建模、加噪、滤波、频谱分析、性能验证。课本上每个知识点单独看都懂一旦要求把卷积、FFT、窗函数、FIR滤波器放进同一个程序里协同工作问题就暴露了采样率该定多少滤波器阶数怎么选频谱图横轴怎么标噪声加多大才算“合理”。这篇文章以基于MATLAB的实现为主线顺带给出Python的替代点按“先搭骨架、再填模块、最后验证”的顺序把课程设计中真正会卡住的地方逐个拆开。2. 数字信号处理课程设计的开发环境与系统框图先把任务书拆成可执行的模块2.1 选型MATLAB还是Python取决于任务书里有没有“交互界面”课程设计题目的表述方式直接决定工具链。如果任务书写的是“设计一个信号分析与处理系统实现对信号的时域和频域分析”MATLAB是稳妥选择信号处理工具箱函数齐全、绘图交互成本低、答辩时修改参数实时看波形变化很直观。如果任务书明确要求提交可运行的跨平台程序或允许自由选择语言Python配合numpy、scipy、matplotlib也能完整覆盖只是滤波器和窗函数需要习惯scipy.signal的接口风格。我一般建议用目录结构把任务书里的功能点映射成文件而不是写一个几百行的脚本从头跑到底。一个能应对多数题目的结构如下dsp_project/ ├── main.m # 入口参数配置、调用各模块、输出结果图 ├── config_params.m # 所有可调参数集中在一个文件 ├── signal_gen.m # 信号生成正弦波组合、加噪、异常信号 ├── time_analysis.m # 时域波形绘制、卷积、相关计算 ├── freq_analysis.m # 频域FFT、频谱绘制、窗函数对比 ├── fir_filter_design.m # FIR低通/带通滤波器设计与滤波 └── verification.m # 指标验证SNR、MSE、幅值误差这种拆分的核心逻辑任务书评审时老师会问“采样率为什么取这个值”“滤波器阶数怎么定”参数集中放在config_params.m里答辩时改一个数字重新跑一遍就能现场展示参数对结果的影响。模块化还让每个函数可以单独测试不需要为了验证一个FFT写法把整个系统跑起来。2.2 最小骨架先把环境跑通再填业务逻辑不要一上来就写滤波器。先写一个最小程序验证开发环境可用同时把参数体系建立起来。下面的代码是main.m的初始版本只做两件事生成一段信号画出它的时域波形和频谱。%% main.m 最小骨架 clear; close all; clc; run(config_params.m); % 加载采样率、频率、点数等参数 fs 8000; % 采样率 8 kHz N 1024; % 分析点数 f1 50; % 基波频率 50 Hz f2 650; % 谐波/干扰频率 650 Hz t (0:N-1) / fs; % 时间序列单位秒 x sin(2*pi*f1*t) 0.5*sin(2*pi*f2*t); % 双频叠加信号 figure(Name, 基础信号); subplot(2,1,1); plot(t, x); title(时域波形); xlabel(时间 (s)); ylabel(幅值); grid on; X fft(x); f_axis (0:N-1) * fs / N; % 频率轴单位 Hz subplot(2,1,2); plot(f_axis, abs(X) / N); % 幅值谱归一化 title(幅值谱); xlabel(频率 (Hz)); ylabel(幅值); xlim([0 fs/2]); grid on;这段代码里有几个参数是后续所有模块都要依赖的fs 8000决定了频率分析的上限是 4 kHz奈奎斯特频率N 1024决定了频率分辨率是df fs / N ≈ 7.8125 Hz。f1 取 50 Hz 是工频f2 取 650 Hz 是为了模拟一个远离基波的干扰分量。fft之后用abs(X) / N归一化是因为MATLAB的FFT结果没有做幅度补偿不除以 N 的话幅值谱的峰值会是实际幅值的 N 倍。验证环境跑通的标准是频谱图上 50 Hz 和 650 Hz 处各有一个峰值幅值分别接近 1 和 0.5其他地方接近 0。如果 650 Hz 的峰值不明显先检查是不是f_axis的计算出了问题常见错误是忘记乘fs/N导致横轴单位变成“每帧的样点序号”而非 Hz。提示config_params.m 里不要写死变量后再在 main.m 里重新赋值。让 config_params.m 成为唯一参数来源否则改参数时会漏改某一处排查起来非常被动。2.3 参数表采样率、点数、频率分辨率的三方约束参数之间不是独立的。fs和N的关系可以用一张表梳理清楚这张表在写报告时也是理论依据。参数组合频率分辨率 df最大可分析频率适用场景fs8000, N10247.8125 Hz4000 Hz语音频段通用分析fs8000, N40961.9531 Hz4000 Hz需要分辨 2 Hz 以内的相邻频率fs16000, N102415.625 Hz8000 Hz高频成分丰富但分辨率要求不高fs16000, N40963.9063 Hz8000 Hz兼顾高频范围和分辨率提高频率分辨率有两个方向增加 N或者降低 fs。增加 N 意味着处理时延变大、内存开销增大降低 fs 会压缩可分析的频率范围。课程设计里通常先定 fs因为信号源和任务书里给的参考频率是固定的然后再根据“需要分辨的两个频率之间最小间隔”反推 N。例如要区分 495 Hz 和 505 Hz 两个分量df 至少要小于 10 Hz取 df 7.8125 Hz 的 fs8000、N1024 组合就够用。3. 时域分析模块信号生成、卷积与相关的实现和验证3.1 用“信号模型”承接任务书的第一条要求任务书里最常见的信号描述是“带有噪声的周期性信号”或“含有谐波失真的正弦信号”。这些描述落到代码里需要拆成两部分确定性分量和随机性分量。确定性分量用正弦波的线性叠加表示随机性分量用randn生成的高斯白噪声。下面的代码生成一个“50 Hz 基波 150 Hz 三次谐波 可调强度噪声”的测试信号作为后续所有模块的输入。%% signal_gen.m 可靠复现的信号生成 function [x, t] signal_gen(fs, N, f_base, harmonics, noise_amp) % 输入: % fs 采样率单位 Hz % N 采样点数 % f_base 基波频率单位 Hz % harmonics 结构体数组: 每个元素含 freq 和 amp 字段 % noise_amp 噪声幅值标准差0 表示无噪声 % 输出: % x 生成的信号序列 (1xN) % t 时间轴 (1xN) t (0:N-1) / fs; x zeros(1, N); % 叠加各次谐波 for k 1:length(harmonics) x x harmonics(k).amp * sin(2*pi*harmonics(k).freq * t); end % 添加高斯白噪声 if noise_amp 0 x x noise_amp * randn(1, N); end end这个函数刻意把噪声幅度作为独立参数暴露出来而不是写死在内部。这样做的意义是验证“滤波器的去噪效果”时可以用同一个确定性信号分别生成不含噪和含噪两个版本对比滤波后的输出与原始无噪信号的误差。谐波用结构体数组传入是因为任务书可能要求分析“基波 2~5 次谐波”这种数据结构不需要改动函数体就能覆盖。调用方式harmonics struct(freq, {50, 150, 250}, amp, {1.0, 0.3, 0.1}); [x_clean, t] signal_gen(fs, N, 50, harmonics, 0); [x_noisy, ~] signal_gen(fs, N, 50, harmonics, 0.2);noise_amp 0.2的含义是噪声标准差为 0.2对应的信噪比大约是20*log10(1/0.2) ≈ 14 dB属于比较恶劣的条件。如果任务书里没给噪声强度用 0.1~0.3 区间内的值都比较合理太小体现不出滤波的作用太大会让滤波后残余噪声仍然明显。3.2 卷积实现与循环卷积的坑卷积是数字信号处理课程设计里必然会涉及的知识点但任务书很少让你直接调用conv完事通常是要你用 DFT 实现卷积来验证线性卷积与循环卷积的关系。这里有一个高频出错的点直接用fft做乘法再ifft回来得到的是循环卷积不是线性卷积。%% 验证循环卷积与线性卷积的关系 x [1, 2, 3, 4, 5]; h [1, 1, 1]; linear_conv conv(x, h); % 线性卷积长度为 53-17 Nfft length(x) length(h) - 1; % 关键必须补零到线性卷积长度 X fft(x, Nfft); H fft(h, Nfft); circular_conv ifft(X .* H); % 补零后的循环卷积 线性卷积 disp(线性卷积:, linear_conv); disp(循环卷积(补零后):, real(circular_conv));运行这个脚本会发现两次结果完全一致。但如果把Nfft直接取成length(x)而不补零结果会变成循环卷积中间几个点被首尾叠加污染和conv的结果对不上。任务书里如果要你“用 FFT 加速卷积”这个补零条件Nfft length(x) length(h) - 1是必须写在报告推导里的。3.3 相关计算怎么用互相关在噪声里找周期相关分析的典型任务是“从带噪信号中检测已知频率分量是否存在”。自相关函数对周期信号在滞后等于整数倍周期处会出现峰值而噪声的自相关在非零滞后处迅速衰减。这一条性质可以用来估计信号的基波周期。%% 用自相关估计基波周期 [x_noisy, t] signal_gen(fs, N, 50, harmonics, 0.3); [acf, lags] xcorr(x_noisy, coeff); % 归一化自相关 positive_lags lags 0; % 只看正滞后 [~, max_idx] findpeaks(acf(positive_lags), MinPeakHeight, 0.3); period_samples lags(positive_lags); detected_period period_samples(max_idx(1)); % 第一个显著峰 fprintf(检测到基波周期: %.2f ms (理论值 %.2f ms)\n, ... 1000*detected_period/fs, 1000/50);xcorr的coeff选项把自相关峰值归一化到 1MinPeakHeight0.3过滤掉噪声引起的伪峰。如果检测出的周期对应频率和 50 Hz 偏差超过 5%检查两个地方一是 N 是否太短导致自相关估计方差过大通常 N 要包含至少 10 个信号周期二是噪声幅度是否超过 0.5过强的噪声会让自相关峰被淹没。4. 频域分析模块FFT 的实用参数、窗函数对比与 FIR 滤波器设计4.1 画对频谱图单边谱、幅度归一化和频率轴刻度任务书里“画出信号的频谱图”是最常见的要求但频谱图画法有好几种直接决定图上读出的幅值是否正确。对于实信号双边频谱在正负频率处各有一半能量单边谱把正频率部分的幅值乘以 2 以还原真实幅度。下面的代码实现了标准单边谱绘制%% 单边谱绘制函数 function plot_spectrum(x, fs, title_str) N length(x); X fft(x); half_N floor(N/2) 1; spec abs(X(1:half_N)) / N; spec(2:end-1) spec(2:end-1) * 2; % 除直流外乘以2恢复幅度 freq (0:half_N-1) * fs / N; plot(freq, spec); title(title_str); xlabel(频率 (Hz)); ylabel(幅值); grid on; end这段代码的关键是spec(2:end-1) spec(2:end-1) * 2。直流分量下标 1和奈奎斯特频率分量如果 N 是偶数下标端点是 fs/2在双边谱里没有对称的镜像所以不乘 2。如果把所有点都乘 250 Hz 处的幅值会变成 2和真实幅值 1 对不上答辩时老师一问就露馅。4.2 频谱泄漏的根源与窗函数参数表直接用fft分析有限长信号时如果信号频率不是频率分辨率的整数倍能量会泄漏到相邻频点。这个问题在课程设计里的经典场景是f1 50 Hzfs 8000N 1000 时df 8 Hz50 不是 8 的整数倍频谱峰值会偏离且幅值偏低。解决办法是加窗但不同窗的代价不同。窗函数主瓣宽度近似最大旁瓣衰减幅值精度适用场景矩形窗2df-13 dB最高但泄漏严重频率间隔大且幅值相近汉宁窗4df-31 dB居中通用分析课程设计首选汉明窗4df-41 dB略低旁瓣要求更高时布莱克曼窗6df-57 dB较低强干扰近邻弱信号加窗后的幅值会变化因为窗函数等效于在时域乘一个衰减因子。用汉宁窗时幅值需要乘以 2 补偿周期窗的相干增益约 0.5。如果任务书要求对比“加窗前后的频谱差异”代码里要同时画出两种结果并标注频谱分辨率的变化。4.3 FIR 滤波器设计指标怎么定、参数怎么调滤波任务是课程设计的重头戏要求通常是“设计一个 FIR 低通滤波器滤除 650 Hz 的干扰保留 50 Hz 基波”。最常用的设计函数是fir1但它的截止频率参数需要归一化到奈奎斯特频率且阶数的选取直接影响过渡带宽度和计算量。%% FIR 低通滤波器设计与滤波 fs 8000; f_pass 100; % 通带截止频率 100 Hz f_stop 300; % 阻带起始频率 300 Hz fl f_pass / (fs/2); % 归一化通带边界 0.025 fh f_stop / (fs/2); % 归一化阻带起始 0.075 N_fir 64; % 滤波器阶数 h fir1(N_fir, [fl, fh], low); % 设计 64 阶低通 FIR freqz(h, 1, 1024, fs); % 查看幅频/相频响应 x_filtered filter(h, 1, x_noisy); % 零相位滤波需要与 filtfilt 对比fir1的第二个参数如果是一个双元素向量表示通带与阻带的边界中间是过渡带low指定低通。阶数N_fir 64对应滤波器长度为 65过渡带宽度近似为4 / N_fir * fs / 2 250 Hz从 100 Hz 到 300 Hz 的 200 Hz 过渡带刚好可以被容纳。如果任务书要求更陡的过渡带阶数要相应提高但代价是群延迟变大时域波形起点会出现更多的预振铃。filter和filtfilt的差别是filter输出有相位延迟严格意义上这是因果滤波的正常行为filtfilt做零相位滤波先正向后反向没有相位失真但引入了非因果性。课程设计里如果重点是“滤除干扰”而不关心实时性建议用filtfilt因为滤波器引入的相位偏移会干扰后续的“滤波后信号与原始信号对比”环节。x_filtered filtfilt(h, 1, x_noisy); % 零相位适合离线分析滤波效果的直观验证是对比滤波前后的频谱650 Hz 处的峰值应明显下降50 Hz 幅值应基本不变。如果 650 Hz 处残余峰仍高于基波的 20%说明过渡带设置过宽或阻带衰减不够优先调整N_fir和f_stop的值。注意不要用filter滤波后直接测量输出的信噪比而不考虑群延迟那会把对齐误差也算进噪声里。离线分析场景优先filtfilt。5. 集成验收把分散的数字信号处理模块连成一个可评审的完整系统5.1 量化指标SNR、MSE 和幅值误差缺一不可任务书评审的标准不是“程序能跑”而是“结果能说明问题”。三个指标能覆盖大多数题目要求滤波前后的信噪比SNR改善、滤波输出与无噪参考信号的均方误差MSE、基波幅值测量误差。它们比“从图上看滤得挺干净”更有说服力。%% verification.m 性能指标计算 function [snr_in, snr_out, mse_val, amp_err] verify_system(x_noisy, x_clean, x_filtered, fs, f_base) % 输入 SNR以基波为信号 signal_power mean(x_clean.^2); noise_in mean((x_noisy - x_clean).^2); snr_in 10*log10(signal_power / noise_in); % 输出 SNR residual x_filtered - x_clean; noise_out mean(residual.^2); snr_out 10*log10(signal_power / noise_out); % MSE mse_val mean((x_filtered - x_clean).^2); % 基波幅值测量 N length(x_filtered); X fft(x_filtered); [~, idx] min(abs((0:N-1)*fs/N - f_base)); % 找基波对应频点 amp_measured abs(X(idx)) / N * 2; amp_err abs(amp_measured - 1.0); % 假设基波理论幅值为 1 end实测一组典型参数fs8000, N4096, 噪声标准差 0.2滤波前 SNR ≈ 14 dB滤波后 SNR 通常在 30 dB 以上MSE 在 0.01 量级幅值误差取决于窗函数类型和基波频率是否落在频点整数倍上。这组数字可以当作自测基准如果滤波后 SNR 提升不到 10 dB先怀疑滤波器的通带设得太窄把基波削了如果 MSE 很大但 SNR 提升正常大概率是滤波输出和 x_clean 之间存在相位偏差换成filtfilt再测。5.2 批量验证技巧用一个脚本把所有题目参数跑完课程设计答辩时经常会被问到“参数稍微改一下结果会怎样”现场手改参数容易出错且浪费时间。一个高效做法是把所有待测参数组合写成一个循环输出对比表%% batch_test.m 批量参数验证 noise_levels [0.05, 0.1, 0.2, 0.4]; filter_orders [32, 64, 128]; results []; for ni 1:length(noise_levels) for fi 1:length(filter_orders) % 重新生成信号、重新设计滤波器并滤波 % ... results(end1, :) [noise_levels(ni), filter_orders(fi), snr_out, mse_val]; end end disp(array2table(results, VariableNames, ... {NoiseAmp, FilterOrder, SNR_out_dB, MSE}));批量结果的价值在于能看出规律噪声增大时 SNR 改善量是否稳定、滤波器阶数提高后 MSE 是下降还是出现平台期。这些规律本身就是报告里“参数影响分析”章节的素材不需要额外编写分析程序。5.3 答辩演示的呈现顺序与一个信号对齐技巧演示时的顺序建议是先展示含噪信号的时域和频域特征说明“能看出有干扰但频域才能定位干扰频率”再展示滤波器幅频响应最后展示滤波后的频谱和指标表。这个逻辑链正好对应任务书“分析—设计—验证”的评审点。信号对齐是演示中容易翻车的地方如果把滤波前后波形直接叠加绘图由于filtfilt不引入延迟两条曲线的高频细节在视觉上会有相位差容易让评审误以为滤波引入了畸变。一个有效的对齐技巧是计算滤波输出与输入在前 20 个采样点的互相关系数如果最大相关点出现在滞后 0 处说明对齐正确否则用alignsignals重新对齐后再绘图。MATLAB 的alignsignals函数可以自动估计两个信号之间的延迟并对其[x_align, y_align] alignsignals(x_filtered, x_noisy);这个函数用于对比展示不改变滤波输出本身只是让波形在图上重叠得更整齐避免因为视觉上的细节错位干扰评审对滤波效果的判断。课程设计的收尾工作往往就是把这种细节打磨到位。本文还有配套的精品资源点击获取