ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

信号分析与处理实验1-4:采样、FFT、卷积与滤波器设计代码实现

信号分析与处理实验1-4:采样、FFT、卷积与滤波器设计代码实现 简介这份资源是南京邮电大学「信号分析与处理实验」课程的完整实验报告面向正在修读数字信号处理、信号与系统相关课程的高校学生以及需要借助 MATLAB 完成实验与课程设计的学习者。报告围绕信号的产生和运算、连续时间信号的频域分析、信号的时域采样和重建、离散傅里叶变换与应用等核心实验展开并延伸至连续时间系统的频域与复频域分析、数字滤波器设计及信号处理综合应用。包内为 1 个 docx 文档压缩包约 1.02MB内容包含各实验的目的、内容、MATLAB 程序代码、波形图与结果分析可直接对照复现单位阶跃、复指数、正弦等序列的生成与卷积、奇偶分解等操作。目前已有 951 人学习下载适合作为实验预习、报告撰写与 DSP 知识点查漏补缺的参考帮助读者直观理解信号时域与频域特性夯实数字信号处理基础。1. 信号分析与处理实验的四个台阶从采样到滤波到底在练什么很多同学拿到“信号分析与处理实验报告实验1-4”这个任务时第一反应是找模板、抄波形、凑结论。但真正决定报告能不能拿高分、答辩能不能扛住提问的不是排版而是你有没有搞懂这四个实验串起来的那条主线连续信号怎么变成离散序列离散序列怎么在频域被看清系统又怎么把不需要的频率成分滤掉。南京邮电大学这套实验通常围绕采样与重构、离散傅里叶变换DFT/FFT、卷积与系统响应、数字滤波器设计展开实验1到实验4正好对应“时域采样—频域分析—系统作用—滤波器落地”的递进。它适合通信、电子、信息工程方向、已经学过信号与系统但动手偏弱的同学。下面按可复现的方式把每个实验的代码、参数和踩坑点讲清楚让你能自己跑出波形而不是对着别人的截图改数据。2. 实验1与实验2采样、重构和DFT/FFT的代码实现2.1 采样定理在代码里怎么验证从连续正弦到离散序列实验1的核心是采样定理。理论说采样频率必须大于信号最高频率的两倍但报告里只写这句话没有分。你要做的是用代码把“欠采样导致混叠”这个现象画出来。常见做法是用一个连续正弦信号取不同采样率生成离散序列再尝试重构。import numpy as np import matplotlib.pyplot as plt # 连续信号参数频率 f0 5 Hz f0 5 t_continuous np.linspace(0, 1, 10000) # 高密度模拟连续时间 x_continuous np.sin(2 * np.pi * f0 * t_continuous) # 采样频率分别取 50 Hz满足和 8 Hz不满足2*f010 fs_list [50, 8] for fs in fs_list: n np.arange(0, int(fs * 1)) # 1秒内的采样点索引 t_sample n / fs x_sample np.sin(2 * np.pi * f0 * t_sample) plt.figure() plt.plot(t_continuous, x_continuous, labelcontinuous) plt.stem(t_sample, x_sample, linefmtr-, markerfmtro, basefmt , labelffs{fs}Hz) plt.title(fSampling with fs{fs} Hz) plt.xlabel(Time (s)) plt.legend() plt.show()这段代码的逻辑是先用高密度点画出“连续”参考波形再用np.arange按采样率生成离散时间点stem画出采样值。当fs50Hz时离散点能还原出5Hz正弦的形状当fs8Hz时离散点看起来像是一个更低频率的信号这就是混叠。参数上f0是信号频率fs是采样频率判断依据是fs 2*f0。报告里要写清楚欠采样后频谱发生搬移高频成分被折叠到低频重构时无法恢复原信号。注意np.stem画出的离散序列在报告里比plot更规范因为采样信号本质是序列不是连续曲线。2.2 DFT与FFT频谱分辨率、补零和泄漏的三个必调参数实验2通常要求对离散序列做DFT或FFT观察频谱。很多人直接调np.fft.fft就完事但报告里如果没解释频率轴怎么来的、补零有什么用、泄漏怎么产生的分数就上不去。下面是一个完整的频谱分析代码。import numpy as np import matplotlib.pyplot as plt fs 1000 # 采样频率 1000 Hz N 200 # 采样点数 n np.arange(N) f0 50 # 信号频率 50 Hz x np.sin(2 * np.pi * f0 * n / fs) # 不加窗直接FFT X np.fft.fft(x) freq np.fft.fftfreq(N, d1/fs) # 频率轴单位Hz # 只取正半轴 half N // 2 plt.figure() plt.plot(freq[:half], np.abs(X[:half]) / N * 2) plt.title(FFT without window, N200) plt.xlabel(Frequency (Hz)) plt.ylabel(Amplitude) plt.grid(True) plt.show() # 补零到 2048 点观察频谱插值效果 N_zero 2048 X_zero np.fft.fft(x, nN_zero) freq_zero np.fft.fftfreq(N_zero, d1/fs) plt.figure() plt.plot(freq_zero[:N_zero//2], np.abs(X_zero[:N_zero//2]) / N * 2) plt.title(FFT with zero-padding to 2048) plt.xlabel(Frequency (Hz)) plt.grid(True) plt.show()逻辑说明np.fft.fftfreq(N, d1/fs)生成频率轴d是采样间隔。幅度归一化用/N*2是为了让单频正弦的峰值接近真实幅值1。补零到2048点并不会提高真实频率分辨率但会让频谱曲线更平滑方便读数。频率分辨率由fs/N决定这里1000/2005Hz所以50Hz刚好落在第10根谱线上没有泄漏。如果信号频率改成55Hz就会看到频谱能量扩散这就是泄漏。参数上N决定分辨率fs决定奈奎斯特上限补零长度只影响显示密度。参数含义典型设置对结果的影响fs采样频率1000 Hz决定不混叠的最高频率N采样点数200决定频率分辨率 fs/N补零长度FFT点数2048不提高分辨率只平滑频谱窗函数减少泄漏汉宁窗主瓣变宽旁瓣降低提示报告里如果要求“分析频谱分辨率”一定要写出Δf fs/N这个公式并说明补零不能改变Δf。3. 实验3卷积、系统响应和时域到频域的对应关系3.1 用代码验证卷积定理时域卷积等于频域相乘实验3通常涉及离散系统的卷积运算。理论课上讲卷积定理但自己动手验证一遍报告里就能写出别人写不出的细节。常见做法是构造一个输入序列和一个系统冲激响应分别做时域卷积和频域相乘比较结果。import numpy as np import matplotlib.pyplot as plt # 输入序列 x[n]长度为8的矩形窗 x np.ones(8) # 系统冲激响应 h[n]指数衰减 h np.array([1, 0.8, 0.6, 0.4, 0.2]) # 时域卷积 y_time np.convolve(x, h) # 频域相乘先补零到相同长度再做FFT L len(x) len(h) - 1 X np.fft.fft(x, nL) H np.fft.fft(h, nL) Y_freq X * H y_freq np.fft.ifft(Y_freq) plt.figure() plt.subplot(2,1,1) plt.stem(y_time, linefmtb-, markerfmtbo, basefmt ) plt.title(Time-domain convolution) plt.subplot(2,1,2) plt.stem(np.real(y_freq), linefmtr-, markerfmtro, basefmt ) plt.title(Frequency-domain multiplication then IFFT) plt.tight_layout() plt.show() print(Max difference:, np.max(np.abs(y_time - np.real(y_freq))))逻辑说明np.convolve直接做时域卷积长度是len(x)len(h)-1。频域方法要求两个序列补零到相同长度L否则X*H是逐点相乘但长度不匹配。np.fft.ifft的结果可能有极小虚部取实部即可。最后打印最大误差通常在1e-15量级说明两种方法等价。参数上L必须大于等于卷积结果长度否则会发生循环卷积混叠。报告里可以写频域相乘对应时域卷积但必须补零到足够长度否则得到的是循环卷积而不是线性卷积。3.2 系统频率响应从差分方程到幅频特性曲线实验3还可能要求分析一个差分方程描述的系统。比如y[n] 0.5*y[n-1] x[n]这是一个一阶IIR系统。你可以用scipy.signal.freqz直接画幅频响应。import numpy as np import matplotlib.pyplot as plt from scipy import signal # 差分方程y[n] - 0.5*y[n-1] x[n] # 系数a [1, -0.5], b [1] b [1] a [1, -0.5] w, h signal.freqz(b, a, worN512) plt.figure() plt.plot(w / np.pi, 20 * np.log10(np.abs(h))) plt.title(Frequency response of y[n]0.5y[n-1]x[n]) plt.xlabel(Normalized frequency (×π rad/sample)) plt.ylabel(Magnitude (dB)) plt.grid(True) plt.show()逻辑说明signal.freqz返回归一化角频率w和复数频率响应h。20*np.log10(np.abs(h))得到幅频特性分贝值。这个系统在低频增益接近1高频衰减是一个低通滤波器。参数上b是前馈系数a是反馈系数worN是频率采样点数。报告里要写出极点位置决定系统稳定性这里极点在z0.5单位圆内系统稳定。差分方程ba系统类型稳定性y[n]x[n]0.5x[n-1][1,0.5][1]FIR恒稳定y[n]0.5y[n-1]x[n][1][1,-0.5]IIR极点0.5稳定y[n]2y[n-1]x[n][1][1,-2]IIR极点2不稳定注意IIR系统的稳定性看极点是否在单位圆内报告里不要只画幅频曲线就结束要补一句极点分析。4. 实验4数字滤波器设计与报告里的参数对照4.1 FIR滤波器设计窗函数法和参数怎么填实验4通常是设计一个数字滤波器比如低通、高通或带通。FIR滤波器常用窗函数法因为设计过程直观参数好写进报告。下面以设计一个低通FIR为例截止频率归一化到0.3。import numpy as np import matplotlib.pyplot as plt from scipy import signal # 设计低通FIR截止频率0.3归一化阶数64 numtaps 65 # 滤波器长度奇数 cutoff 0.3 # 归一化截止频率Nyquist1 b signal.firwin(numtaps, cutoff, windowhamming) # 查看频率响应 w, h signal.freqz(b, worN1024) plt.figure() plt.plot(w / np.pi, 20 * np.log10(np.abs(h))) plt.title(FIR low-pass, cutoff0.3, Hamming window) plt.xlabel(Normalized frequency (×π rad/sample)) plt.ylabel(Magnitude (dB)) plt.grid(True) plt.show() # 打印滤波器系数方便写进报告 print(Filter coefficients:, np.round(b, 4))逻辑说明signal.firwin的numtaps是滤波器阶数加1奇数长度适合设计高通或带通。cutoff是归一化频率1对应奈奎斯特频率。windowhamming控制旁瓣衰减。阶数越高过渡带越窄但计算量越大。报告里要写清楚截止频率、窗函数类型、阶数、过渡带宽度和阻带衰减。参数上numtaps每增加一倍过渡带大约减半但阻带衰减由窗函数决定汉明窗约53dB。4.2 IIR滤波器设计巴特沃斯和切比雪夫怎么选IIR滤波器用signal.butter或signal.cheby1设计阶数低但相位非线性。实验报告里常见要求是设计一个巴特沃斯低通通带截止0.2阻带截止0.4通带波纹和阻带衰减给定。import numpy as np import matplotlib.pyplot as plt from scipy import signal # 巴特沃斯低通通带截止0.2阻带截止0.4 wp 0.2 ws 0.4 gpass 1 # 通带最大衰减 dB gstop 40 # 阻带最小衰减 dB # 计算最小阶数 N, wn signal.buttord(wp, ws, gpass, gstop) print(fButterworth order: {N}, natural frequency: {wn:.4f}) b, a signal.butter(N, wn, btypelow) w, h signal.freqz(b, a, worN1024) plt.figure() plt.plot(w / np.pi, 20 * np.log10(np.abs(h))) plt.title(fButterworth low-pass, order{N}) plt.xlabel(Normalized frequency (×π rad/sample)) plt.ylabel(Magnitude (dB)) plt.grid(True) plt.show()逻辑说明signal.buttord根据通带、阻带指标自动计算最小阶数N和自然频率wn。signal.butter返回分子分母系数。巴特沃斯幅频特性最平坦切比雪夫允许通带波纹换取更陡的过渡带。报告里要对比同样指标下巴特沃斯阶数更高切比雪夫阶数更低但通带有波纹。参数上gpass和gstop是设计指标wp和ws是归一化频率。滤波器类型函数特点适用场景巴特沃斯signal.butter通带最平坦一般低通、高通切比雪夫Isignal.cheby1通带等波纹过渡带要求陡切比雪夫IIsignal.cheby2阻带等波纹阻带衰减要求高椭圆signal.ellip通带阻带都等波纹阶数最低提示报告里如果要求“比较不同滤波器”不要只贴幅频曲线要列出阶数、过渡带宽度和相位非线性程度。5. 实验报告排错与验证波形不对时先查这五个地方5.1 频谱峰值对不上先查频率轴和归一化很多同学跑完FFT发现峰值不在预期频率最常见的原因是频率轴没算对。np.fft.fftfreq返回的是周期数单位是Hz但如果你用np.arange(N)当频率轴那就全错了。另一个坑是幅度归一化np.abs(X)不除以N峰值会随N变化。正确做法是单频信号除以N/2直流分量除以N。报告里写清楚归一化方式答辩时老师一看就知道你懂。5.2 滤波器输出异常检查系数顺序和初始状态用signal.lfilter滤波时如果输出一开始有剧烈震荡通常是初始状态没设置。lfilter默认零初始状态但IIR滤波器从零状态启动会有瞬态。可以用signal.lfilter_zi计算稳态初始条件或者直接丢弃前若干个输出点。另外b和a的顺序不能反b是前馈a是反馈反了系统就完全变了。from scipy import signal import numpy as np b, a signal.butter(4, 0.2) zi signal.lfilter_zi(b, a) x np.random.randn(100) y, _ signal.lfilter(b, a, x, zizi*x[0])这段代码用lfilter_zi让滤波器从稳态开始避免启动瞬态。参数上zi*x[0]是让初始状态匹配第一个输入样本。5.3 报告数据可复现固定随机种子和记录参数如果实验里用了随机信号一定要在代码开头写np.random.seed(42)否则每次跑出来的波形不一样报告里的图对不上。另外所有关键参数都要在报告里列出来采样率、点数、截止频率、阶数、窗函数类型。这些参数就是实验的可复现性保证。老师如果让你现场重跑你能直接改参数跑出一样的图这比任何结论都有说服力。5.4 用理论值交叉验证别只信一条曲线验证实验结果时不要只看波形“长得像”。比如设计低通滤波器要检查通带内增益是否接近0dB阻带内是否低于设计指标。可以用signal.freqz输出具体数值和理论值对比。卷积实验里时域卷积和频域相乘的最大误差应该在1e-10以下。这些数值验证比截图更有说服力也是报告里最容易拿分的地方。5.5 常见报错和对应处理ValueError: Length of input signal must be a multiple of hop length这类错误通常出现在用scipy.signal做频谱图时检查nperseg和信号长度。LinAlgError出现在滤波器设计阶数过高时降低阶数或改用sos格式。np.fft相关错误多半是输入不是一维数组用np.ravel展平。把这些报错和处理方式写进报告的“问题与解决”部分比只写“程序运行成功”更有内容。本文还有配套的精品资源点击获取
返回列表