ARTICLE DETAIL

资讯详情

深耕网站视觉设计与运营推广的一线实战洞察。

信号与系统实验:采样率、FFT与可复现仿真

信号与系统实验:采样率、FFT与可复现仿真 简介北京理工大学信号与系统实验报告完整记录了基于MATLAB的信号时域描述与运算实验是信息工程类本科生学习信号与系统课程的实用参考。资源面向初学者系统梳理连续时间信号与离散时间信号的向量表示法、符号对象表示法并逐一讲解信号相加、相乘、微分、积分以及位移、反转、尺度变换等时域运算的MATLAB实现途径。报告结合常用函数与绘图命令给出多组典型连续信号和离散信号的波形绘制实例覆盖heaviside、rectpuls、tripuls、square、sinc、sawtooth等常用信号产生函数便于读者理解可视化操作和编程细节。压缩包内共1个docx文档大小2.23MB已有1461人学习下载。这份实验报告既适合实验课前预习、课堂对照也可作为课后撰写实验报告的完整范本能有效帮助读者掌握信号时域特性的分析方法。1. 信号与系统实验报告从公式表到可复现仿真北理工信号与系统实验报告最值得写的不是那几页手推公式而是三个数字采样率 fs、观察时长 T、FFT 点数 N。这三个参数一旦含糊时域波形、频谱峰值和卷积结果就都没有参照系。也正因如此做这份报告的思路和工程师做仿真验证是一致的先固定实验载体再写实现代码最后用图表和误差指标把现象钉在结论里。本文从连续时间系统实验和离散时间系统实验两条主线出发给出可直接改用的实验流程、代码段和参数表并在最后列出报告里常见又不容易被看出的六个边界问题。新手照着走能复现老手可以直接跳到自己需要的参数说明。2. 信号与系统实验的载体设计采样率、时窗与依赖版本实验报告如果一上来就写表达式和代码后边很容易在数据解释上自相矛盾。先花几行把实验对象说清楚采样频率应当至少高于关心的最高频率两倍这是老老实实做完实验的硬约束同时时间窗长度决定频率分辨率想分辨相近的两个频点必须拉长观察窗。下面先给出实验设计时最常用的参数表再说明怎么在 Python 环境下把这个参数表放进代码里。2.1 采样率、观察时长和点数怎么定常见的参数选取会按实验内容分成几类我一般先在报告里放一张表把每个实验的采样条件一次写清实验内容关心最高频率 f_max推荐采样率 fs观察时长 T备注一阶 RC 冲激响应100 Hz1000 Hz 以上大于 5 倍时间常数尾部不能截断二阶振荡系统频响200 Hz1000 Hz包含 10 个振荡周期频率分辨率等于 1/T离散卷积验证无采样要求整数便于对齐约等于 x 与 h 的点数和full 和 same 边界不同信号加窗谱分析500 Hz2000 Hz0.5 s 以上同时看主瓣和泄漏fs 的下限由奈奎斯特条件给出工程上通常取关心频率的 5 到 10 倍FFT 分辨率是 fs/N所以 N 不一定是 2 的幂能被 fs 和 T 整除即可。报告里写“fs1000、T1”比只写“N1000”信息更完整前者能看出时间轴真实跨度后者看不出有没有补零。2.1.1 频率轴与点数换算的固定写法把频率轴放进代码时可以用下面这一段前后台对应起来import numpy as np fs 1000.0 # 采样率单位 Hz T 1.0 # 观察时长单位 s N int(fs * T) # 总采样点数无补零时等于 fs*T t np.arange(N) / fs df fs / N # 频率分辨率 freq np.fft.fftfreq(N, d1.0 / fs)这段代码里df fs / N是频谱中相邻两根谱线的间隔代表报告里能不能分辨两个邻近频点fftfreq返回的是正负频率两半画单边谱时需要手动截取。不要直接写freq np.arange(N) / N * fs虽然数值等价但少了负频率说明后文做逆变换时容易混。2.2 用 SciPy 搭出最小可复现环境报告不一定要用 MATLAB。Python 这边用 NumPy、SciPy、Matplotlib 足以覆盖课程实验而且把代码和图表放在同一个脚本或 Jupyter 文件里更容易追踪。先把环境固定下来python -m venv .venv source .venv/bin/activate # Windows 下用 .venv\Scripts\activate pip install numpy scipy matplotlib装完后在报告的“实验环境”一节贴两行验证代码import numpy as np import scipy.signal as sig print(numpy, np.__version__) # 输出当前 numpy 版本 print(scipy, sig.__version__) # 输出当前 scipy 版本第一段命令创建独立虚拟环境避免把包装进系统 Python第二段打印版本号是为了让拿到报告的人能判断数值差异是环境造成的还是算法造成的。SciPy 的signal子模块提供lfilter、freqz、fftconvolve等核心函数版本差异主要体现在滤波器边界处理与数值细节上所以版本应当固定。如果学校要求 MATLAB思路完全一致这边函数名和 MATLAB 高度对应方便对照。2.3 输入信号离散步的约定加窗与去直流实验信号写“x(t)sin(2π·50t)”之后还要在代码里明确离散化方式。最常见做法是把连续信号直接按采样间隔切出来A 1.0 # 幅值单位 V f_sig 50.0 # 信号频率单位 Hz x A * np.sin(2 * np.pi * f_sig * t) # x[n] 对应 x(t_n)建议先小幅移去直流分量去除由数值计算和量化引入的常数偏移但画图时保留原始幅值单位。这里的t由np.arange(N)/fs生成因此x[n]其实是对连续信号的等间隔采样n 与时间一一对应。若实验关心频谱泄漏报告里必须写明“对 x 使用矩形窗未再补零”若用了汉宁窗后续峰值的三分贝带宽会改变结论要相应调整。提示报告里不应只有“x np.sin(...)”一句而要写清楚采样率、x 的单位和频率轴的单位。单位缺失是报告中最常见的扣分点。3. 实验核心计算怎么写卷积、FFT 归一化与 lfilter 仿真进入实现后重点不再是函数怎么调用而是每个函数结果对应物理哪个量。卷积输出 y 是关于采样点 n 的序列FFT 输出 X 的单位是幅度谱密度还是幅度谱系数lfilter 的初始条件是什么这些说不清报告结论就会悬空。下面三节分别把这三块写成能直接用的代码。3.1 卷积的三种写法手写、numpy 与 fftconvolve先写手写卷积它的价值在于验证定义而不在于性能def my_conv(x, h): len_out len(x) len(h) - 1 y np.zeros(len_out) for n in range(len_out): for k in range(len(x)): if 0 n - k len(h): y[n] x[k] * h[n - k] return y # 用随机序列与 numpy 实现对照误差应接近 0 x_test np.random.rand(16) h_test np.random.rand(8) np.testing.assert_allclose(my_conv(x_test, h_test), np.convolve(x_test, h_test))循环内n-k是第二个序列的下标边界条件决定了输出长度len(x)len(h)-1用np.testing.assert_allclose做相等性断言比人工看图可靠。这里np.random.rand没有设随机种子后面正式实验建议先固定种子再生成噪声让报告示例可重放。方法复杂度边界控制报告用途手写双重循环O(NM)自己控制但易错展示卷积定义np.convolveO(NM)C 层实现full/same/valid 三种快速验证长度scipy.signal.fftconvolveO(N log N)full 长度按线性卷积处理长序列大数组full/same/valid 的含义分别是输出全长、截断到输入长度、只取完全重叠部分。报告里如果只给一条曲线至少要说明用的哪一种否则后续误差比较没有对齐基准。3.2 连续时间傅里叶变换的数值近似fft 的正确归一化教材里的连续时间傅里叶变换是积分式而np.fft.fft是离散时间序列的变换数值上两者差一个1/fs因子。很多报告只画abs(fft(x))结果纵轴和时域信号幅值对不上。标准写法是X np.fft.fft(x) / fs # 幅度谱密度单位 V/Hz freq np.fft.fftfreq(N, d1.0 / fs) half N // 2 X_pos X[:half] # 正频率部分 freq_pos freq[:half] A_single np.abs(X_pos).copy() A_single[1:] * 2 # 单边谱除直流外能量折半np.fft.fft(x)得到的是带量纲的离散频谱系数除以fs后近似对积分式∫ x(t)e^{-j2πft}dt做矩形法求和。注意这里是除以fs而不是除以 N这是和周期序列 DFT 最大的区别。如果实验对象本来就是离散周期序列则反过来用X/N得到傅里叶级数系数。输出A_single后峰值应接近输入幅值 A。若峰值小很多先检查窗函数再检查归一化。分辨率df1/T不变峰值的微小下降通常来自矩形窗泄漏而不是算法错误。3.3 用 lfilter 仿真 LTI 系统并用 freqz 验证频响离散时间系统实验里最常见的是给差分方程求输出。以下面这个二阶系统为例b np.array([0.2, 0.5, 0.3]) # 分子系数对应 x[n], x[n-1], x[n-2] a np.array([1.0, -0.8, 0.2]) # 分母系数a[0] 必须为 1 y sig.lfilter(b, a, x) # 零初始条件计算零状态响应lfilter的b与a按差分方程组织根在 z 平面上决定系统稳定性。默认初始条件是零因果系统的零状态响应正是这一情况。若要做零输入响应需要手工设置zi参数这里不展开。要验证频响用freqz同时给出幅度和相位w, H sig.freqz(b, a, worN2048, fsfs) # w 的单位是 Hz phase np.unwrap(np.angle(H)) # 相位展开避免 ±π 跳变 group_delay -np.diff(phase) / np.diff(w) w_mid 0.5 * (w[1:] w[:-1]) # 群时延曲线的频率轴中点报告的群时延曲线应该画在w_mid上否则会比真实曲线偏半个采样点。unwrap让相位曲线连续而不是折叠在 ±π 之间线性相位系统尤其需要这一步。4. 信号与系统实验报告的数据组织图表、误差指标与结论对齐实验做完、代码能跑通之后组织数据其实是独立的一道工序。图要在报告中承担信息传递而不是把能画的波形全部堆上去。我的固定做法是一个结论对应一张图一张图对应一个坐标含义再配一个定量误差数字。4.1 报告正文最稳定的三图结构最稳定的结构是左中右三图左图输入x[n]中图系统冲激响应h[n]右图输出y[n]。这张图能直接回答“输入是什么、系统是什么、输出是什么”这三个问题。import matplotlib.pyplot as plt fig, ax plt.subplots(1, 3, figsize(12, 3.2), constrained_layoutTrue) ax[0].plot(np.arange(len(x)), x, lw0.8) ax[0].set_title(输入信号 x[n]) ax[0].set_xlabel(采样点 n) ax[1].plot(np.arange(len(h)), h, lw0.8) ax[1].set_title(系统冲激响应 h[n]) ax[1].set_xlabel(采样点 n) ax[2].plot(np.arange(len(y)), y, lw0.8) ax[2].set_title(零状态响应 y[n]) ax[2].set_xlabel(采样点 n)为什么横轴写采样点 n 而不是时间 t如果实验信号由np.arange(N)产生且 fs1000两者数值相同但单位会在频谱部分引起混淆。时间轴写“t/s”采样点轴写“n”这段措辞本身就是在把代码与物理单位对齐。报告若要黑白打印图里还要用实线、虚线、点划线区分数据颜色不作为唯一区分通道。4.2 用误差指标代替“看上去差不多”实验报告常见的弱句是“两者基本一致”。换成标量就能避免这种暧昧。若y_true是解析解或高精度参考y_sim是本次实验输出e y_true[:len(y_sim)] - y_sim mae np.mean(np.abs(e)) # 平均绝对误差 rmse np.sqrt(np.mean(e ** 2)) # 均方根误差 peak np.max(np.abs(e)) # 峰值误差 rho np.corrcoef(y_true[:len(y_sim)], y_sim)[0, 1] # 相关系数指标含义可以汇总成下表指标表达式对实验报告的意义MAEmean(|e|)平均偏移程度RMSEsqrt(mean(e^2))对个别大误差更敏感峰值误差max(|e|)找出最坏采样点相关系数ρ波形趋势一致但看不出常数偏置相关系数很高不代表系统无差别。比如输出整体放大 1.2 倍相关性仍接近 1但 RMSE 会直接暴露问题。所以报告里把 MAE、RMSE 和相关系数一起贴出不要只挑最好看的数字。4.3 结论部分怎么写把现象与公式、参数对齐结论不是心得体会而是把观测结果映射回原理。我写结论时会给固定句型。例如“y 的峰值位于 n203相对输入峰值延迟 3 个采样点与理论群延迟 3.02 个点一致差异来源于 1 Hz 的频率分辨率量化”。这一句包含现象、理论、边界三个信息。另一个常用句“FFT 单边谱在 50 Hz 处幅度为 0.98 V比输入幅值低 2%与矩形窗对主瓣泄漏的估算相符。”这种句式可以直接套进任何实验内容。结论末尾要带一句“上述数值在 fs1000 Hz、T1 s、矩形窗条件下成立”否则评审者无法判断数字是不是只在本地机器上偶然成立。提示结论里不要写“误差来源有很多”。把误差来源具体到采样率、截断长度或窗函数中的某一个报告质量会立刻上一个台阶。5. 信号与系统实验报告最容易翻车的六个边界问题5.1 采样率过低混叠被当成真实频谱当信号含 100 Hz 分量而 fs 取 90 Hz会看到一个被折叠到 10 Hz 的虚假峰。判定方法很简单把 fs 翻倍重新做一次 FFT若峰值明显移动先怀疑混叠。报告里应给出奈奎斯特频率 fs/2并在代码开头加assert fs 2 * f_sig。5.2 FFT 幅度归一化错误把fft(x)直接除以 N 会得到约等于 A/2 的峰值原因是在连续时间积分式中用的是1/fs而不是1/N。自检思路很直接构造 A1 的正弦检查单边谱峰值是否接近 1。1% 到 3% 的误差通常来自频谱泄漏而不是归一化错误。5.3 时域截断与窗函数选择互相矛盾同样一段正弦矩形窗主瓣最窄但旁瓣高Hamming 窗旁瓣低但主瓣展宽。如果实验目的写“优化频率分辨率”却用 Hamming结论就应该讨论主瓣展宽否则会出现峰值下降归因错误。报告里写清“加窗方式是什么、为什么用这个窗”比反复换窗重跑更重要。5.4 lfilter 初始条件不等于稳态lfilter(b, a, x)默认从零初始状态开始输出前几个点是暂态。若实验关心稳态响应要么丢弃前几百点要么先让系统预热若干采样周期。报告里出现“信号开头有明显突变”不是代码 bug是零初始条件下的正常现象应当写进实验观察而不是删掉。5.5 相位图不 unwrap 导致群时延曲线断裂频率响应相位被限制在 ±π 之间直接对np.angle(H)做差分会在 180 度处出现假跳变。正确写法是np.diff(np.unwrap(np.angle(H)))。对线性相位系统这一小段处理能避免一半的返工。5.6 随机信号实验不固定随机种子产生噪声信号时np.random.randn每次结果不同。报告要让别人复现在实验环境部分写上np.random.seed(42)并把 numpy、scipy 版本号一并贴在附录。当有人问“为什么我们俩跑出来的时域包络不一样”时你应该先对比版本号而不是去改代码参数。本文还有配套的精品资源点击获取
返回列表