简介：资源聚焦数字通信中的BPSK调制与采样率变换，面向信号处理、通信原理方向的学生及MATLAB开发者，提供可直接运行的上采样与下采样实现。BPSK通过载波相位传递0/1信息，上采样在序列中插零提高采样率，下采样则按整数因子抽取以降低数据率，三者是构建数字收发信机时反复用到的操作。压缩包内共2个文件，均为.m脚本，总大小仅918B，包含UpSampling.m与DownSampling.m两个核心函数，分别实现插值补零和抽取降速，二者常配合低通滤波器使用，用于防止混叠、匹配符号速率或构造多速率仿真链路。已有296名学习者浏览或下载，适合对照教材理解插值、抽取、抗混叠滤波等概念，并快速嵌入到BPSK调制解调实验中。借助这两个函数，可观察信号频谱在采样率变化前后的形态差异，也可进一步扩展为多速率通信系统的基础模块，对课程设计或科研入门均有直接帮助。
1. 从BPSK的相位跳变说起：为什么要动采样率
如果你在MATLAB里做BPSK链路仿真，一定遇到过这样的问题：符号速率是1 MHz，但后级DAC或射频前端希望输入12.8 MHz采样率的数据，这时候就需要在基带符号序列里插入零值、做脉冲成形，再通过下采样把数据率降到解码器能处理的节奏。这个“插零”和“丢点”的过程就是上采样和下采样。值得注意的反直觉结论是：插零并不会增加信息量，只是把频谱周期变密；而抽取也绝不是每隔 M 点扔一个样本，而是先滤波再丢弃。围绕 UpSampling.m 和 DownSampling.m 两个自定义MATLAB函数，把BPSK基带仿真里采样率转换的数学背景、参数设置和验证技巧一次说透。适合做物理层算法验证、软件无线电原型和通信链路仿真的工程师。
2. 上采样/下采样的数学本质与频谱解释
2.1 上采样就是插零，不是插值
上采样的第一步本质上是序列扩展，也就是在相邻样本间插入 L-1 个零。若原始序列为 x[n]，插零后的序列在 m 为 L 的整数倍时等于 x[m/L]，其余位置为 0。在MATLAB里，这个操作一行就能完成，但它的频谱变化常常被忽略。插零后的序列频谱等于原始频谱在数字频率轴上被压缩 L 倍，也就是 X_u(e^{jw}) = X(e^{jwL})。这意味着原始频谱在 0 到 π/L 之外的成分会被压缩到该区间内，并以 2π/L 为周期复制出镜像。如果直接把插零结果送给DAC，镜像会变成高频噪声，所以完整的上采样必须紧跟一个低通滤波器，滤掉这些镜像。
这也顺便澄清了摘要里“上采样也称为插值”的歧义：单独插零只是 upsample，插零加上抗镜像滤波组合起来才叫插值。很多资料把两者混为一谈，是工程里最常见的概念坑。看下面这段MATLAB代码，能直观看到插零后的频谱复制效果：
fs = 1000; % 原始采样率 1000 Hz f0 = 50; % 信号频率 50 Hz N = 256; % 样本数 n = 0:N-1; x = cos(2*pi*f0/fs*n); % 单频信号 L = 4; % 上采样倍数 xu = zeros(1, L*N); xu(1:L:end) = x; % 插零 figure; subplot(2,1,1); freqz(x, 1, 1024, fs); % 原始频谱 title('Original'); subplot(2,1,2); freqz(xu, 1, 2048, fs*L); % 插零后，采样率变成 fs*L title('After zero-insertion');freqz 的第四参数会按新采样率 fs*L 标出频率轴，你会明显看到原始频谱被复制成 4 份。插零本身不改变傅里叶变换的幅值，但采样率升高后，数字频率被重新归一化，频谱被压缩，镜像随之出现。这就是为什么插值滤波器需要带 L 倍增益：插零让能量分布被稀释，滤波后如果不乘 L，输出幅度会只有原始信号的 1/L。很多人在BPSK仿真时发现星座点缩到 ±0.5 而不是 ±1，多半就是漏了这个增益。
2.2 下采样前必须过防混叠滤波器
抽取器的输出是 y[m] = x[Mm]，M 为抽取倍数。直接丢掉中间样本会让原频谱中高于 π/M 的分量折叠回低频，形成混叠。从时域看，抽取相当于把数字频率轴拉伸 M 倍，原本不重叠的频谱周期叠在一起。为了防止混叠，必须先经过一个截止频率为 π/M 的低通滤波器，把不需要的高频分量滤掉。
这个顺序不能反过来，因为抽取是不可逆的，折叠一旦发生，就无法判断某个低频成分是原始的，还是从高频折叠下来的。工程上防混叠滤波器通常用FIR实现，MATLAB里直接用 fir1 设计，截止频率给 1/M，这是相对原奈奎斯特频率的归一化值。下面的例子构造了一个 30 Hz 与 280 Hz 混合信号，采样率 1000 Hz，抽取 2 倍后新奈奎斯特频率是 250 Hz，280 Hz 分量如果不提前滤除，会折叠到 220 Hz 附近：
fs = 1000; t = 0:0.001:1; x = cos(2*pi*30*t) + 0.5*cos(2*pi*280*t); M = 2; order = 32; h = fir1(order, 1/M); % 截止 250Hz，归一化 0.5 xf = filter(h, 1, x); y_direct = x(1:M:end); % 不滤波直接抽取 y_filtered = xf(1:M:end); % 先滤波再抽取 figure; subplot(2,1,1); pwelch(y_direct); title('Direct decimation'); subplot(2,1,2); pwelch(y_filtered); title('Filtered then decimation');直接抽取的功率谱在 220 Hz 附近会出现一个伪峰，那就是 280 Hz 折叠下来的；而先滤波再抽取的频谱里只有 30 Hz 成分。需要注意的是，防混叠滤波器不需要补偿增益，因为抽取只是减少了样本密度，没有改变样本值本身。相对于上采样乘 L，这个差异就是很多人写 DownSampling.m 时容易出错的地方。
2.3 用频谱图理解采样率转换
| 操作 | 时域表现 | 频谱变化 | 必须配套 |
|---|---|---|---|
| 上采样插零 | 样本间插入 L-1 个零 | 频谱压缩 L 倍，生成镜像 | 低通抗镜像，截止 π/L，增益 L |
| 下采样抽取 | 每 M 点取 1 点 | 频谱展宽 M 倍，高频折叠 | 低通抗混叠，截止 π/M，增益 1 |
这张表是我调链路时最常看的。先判断你正在做哪一步，再决定滤波器截止和增益，很多BER异常都是因为顺序或增益弄反了。MATLAB 内置的resample会自动完成抗镜像/抗混叠滤波，但使用自定义的 UpSampling.m、DownSampling.m 时，这些细节就必须自己管，这也是为什么值得把这两个函数拆开写清楚。
3. BPSK基带信号里的UpSampling.m与DownSampling.m实现
3.1 先搭一个BPSK发射机
先把BPSK基带发射端搭起来。这里用根升余弦脉冲成形，每个符号采样 4 个点，相当于对符号序列做了一个 4 倍上采样：
% 发射参数 M = 4; % 每符号采样数，等价于上采样倍数 nBits = 5000; bits = randi([0 1], nBits, 1); sym = 2*bits - 1; % BPSK: 0 -> -1, 1 -> +1 hTx = rcosdesign(0.35, 6, M, 'sqrt'); % 根升余弦脉冲，span=6 up = upsample(sym, M); % 插零到符号采样率 tx = filter(hTx, 1, up); % 脉冲成形rcosdesign(rolloff, span, sps, 'sqrt')返回根升余弦FIR系数，rolloff 取 0.35 时，带宽是符号速率乘以 1.35。upsample是MATLAB内置上采样，只做插零。我们的 UpSampling.m 在仿真链路里要替代它，并额外加上抗镜像滤波控制。如果没有 Communications Toolbox，也可以用 fir1 按升余弦频率响应近似设计，但这里用rcosdesign更直观，参数也更容易对齐。
3.2 UpSampling.m：插零 + 抗镜像滤波
把 UpSampling.m 实现为一个完整的上采样函数，同时支持只插零和插零后滤波两种模式：
function [y, h] = UpSampling(x, L, type, h) % UpSampling - 上采样（插值） % 输入: % x : 输入序列，行向量或列向量 % L : 上采样倍数，正整数 % type: 'zero' 只插零; 'fir' 插零+低通滤波（默认） % h : 可选低通滤波器系数，默认自动设计 % 输出: % y : 上采样后的序列 % h : 实际使用的滤波器系数 if nargin < 3, type = 'fir'; end if nargin < 4, h = []; end if L == 1 y = x; h = 1; return; end % 第一步：插零 y = zeros(1, length(x)*L); y(1:L:end) = x(:).'; % 如果只要原始插零效果，直接返回 if strcmpi(type, 'zero') return; end % 第二步：设计抗镜像低通 if isempty(h) Ntap = min(64, max(8, 16*L)); % 经验阶数 h = L * fir1(Ntap, 1/L); % 截止 1/L，增益补偿 L end % 滤波 y = filter(h, 1, y); end这个函数把“插零”和“滤波”分成了两步。type='zero' 时，它等价于内置upsample；type='fir' 时才完成完整的插值。滤波器增益乘 L 是为了恢复插零损失的能量，这也是和 DownSampling.m 最大的区别。经验阶数用16*L估算，基本能保证镜像抑制在 40dB 以上；如果信道后级对频谱纯度要求高，再加大到 32*L。
3.3 DownSampling.m：防混叠滤波 + 抽取
对应的 DownSampling.m 实现如下，支持可选的抽取相位参数，方便后面做定时偏移搜索：
function [y, h] = DownSampling(x, M, phase, h) % DownSampling - 下采样（抽取） % 输入: % x : 输入序列 % M : 抽取倍数 % phase : 抽取相位，1~M，默认1 % h : 可选防混叠滤波器，默认 fir1(N, 1/M) % 输出: % y : 抽取后的序列 if nargin < 3, phase = 1; end if nargin < 4, h = []; end if M == 1 y = x; h = 1; return; end if isempty(h) Ntap = min(64, max(8, 16*M)); h = fir1(Ntap, 1/M); % 截止 1/M，不需要增益补偿 end xf = filter(h, 1, x); y = xf(phase:M:end); endphase参数表示从滤波后序列的第几个样本开始抽取，取值范围 1 到 M。默认 phase=1 就是最简单的固定相位。为什么不需要乘 M？因为抽取只是从原样本中选点，没有改变样本的绝对幅度，而上采样插零引入了大量零值，能量被稀释，所以必须乘 L。可以把这个区别记成一句话：上采样要补能量，下采样要防混叠。
在BPSK接收端，DownSampling.m 通常放在匹配滤波之后：
rx = awgn(tx, 10, 'measured'); % 加噪声 hRx = fliplr(hTx); % 匹配滤波 rf = filter(hRx, 1, rx); rxDown = DownSampling(rf, M); % 降到符号率 % 这里 rxDown 已经是每符号 1 个样本 % 由于滤波器延迟，符号中心不一定落在 rxDown(1)，需要定时相位搜索这段代码先加噪声，再用发射滤波器的翻转系数做匹配滤波，最后抽取到符号率。单纯这样还不能直接算误码率，因为抽出来的序列可能有固定相位偏移，所以我把定时搜索放在第 4 章专门讨论。
3.4 参数怎么调：倍数、阶数、截止频率
| 参数 | 含义 | 典型范围 | 调参方向 |
|---|---|---|---|
| L / M | 上/下采样倍数 | 2 / 3 / 4 / 8 | 满足系统级采样率匹配 |
| Ntap | 滤波器阶数 | 8 ~ 64 | 阶数越高过渡带越窄，但延迟越大 |
| 截止频率 | 归一化到Nyquist | 上采样 1/L，下采样 1/M | 不能超过该值，否则混叠或镜像残留 |
| 增益 | 滤波器幅度 | 上采样乘 L，下采样不乘 | 不乘 L 时星座点幅度不对 |
| rolloff | 成形滚降 | 0.2 ~ 0.5 | 越小带宽窄，但实现滤波器更难 |
我一般先设 L/M=2，滤波器阶数 32，跑通链路后看星座图和误码率，再逐步提高倍数。高倍数场景如 L=16 时，16*L 的滤波器阶数会超过 64 被钳制，镜像抑制可能不够，此时需要手动传入更高阶的 h。不要迷信“有滤波器就行”，镜像抑制不足时 BER 曲线会出现在高信噪比段不收敛的情况。
4. 采样率转换的工程坑：从误码率到滤波器延迟
4.1 奈奎斯特与滚降系数之间的妥协
BPSK基带信号的单边带宽约等于符号速率的一半，实际使用根升余弦成形后，双边带宽变为 (1+β)Rs，β 是滚降系数。上采样倍数 L 如果很小，比如 L=2，镜像就会落在离主谱很近的位置。以符号速率 Rs 为参考，镜像中心在 L*Rs 处，但镜像的滚降边沿会和主谱边沿靠近，抗镜像滤波器的过渡带被压缩得极窄，必须用高阶滤波器才能分离。
一个比较实用的经验是：L 至少取 4，让镜像中心落在主谱带宽外足够远的位置。仿真场景下镜像是“频域假象”，但到了真实DAC或射频前端，镜像会混入带外发射，导致频谱模板超标。如果你在做 AD9361 这类射频前端的数据接口设计，通常希望基带采样率至少是符号率的 4 到 8 倍，就是为了给抗镜像滤波留出过渡带空间。
4.2 群延迟和符号定时偏移
FIR滤波器只要系数对称，群延迟就是常数 (N-1)/2 个样本。在 UpSampling.m 和 DownSampling.m 级联后，总延迟可能不是 L 或 M 的整数倍，符号采样点就会整体偏移。最常见的现象是直接取rxDown(1:M:end)作为符号序列，结果星座图旋转或 BER 接近 0.5。
匹配滤波器和抽取器都引入延迟，但抽取器只选固定相位，所以真正要解决的是：在 M 个可能的抽取相位里，找出离符号中心最近的那一个。下面这段代码用BPSK判决结果做参考，扫描相位并计算欧式距离，选择误差最小者：
rxFilt = filter(hRx, 1, rx); % 匹配滤波后的序列，每符号 M 点 evmAll = zeros(1, M); for k = 1:M probes = rxFilt(k:M:end); % 第 k 个抽取相位 probes = probes(:); ref = sign(probes); % BPSK 参考: +1 / -1 evmAll(k) = mean(abs(probes - ref).^2); end [~, bestOffset] = min(evmAll);这段代码不需要知道滤波器群延迟是多少，直接幅度搜索。在无噪声或中高信噪比下，sign(probes)作为参考足够可靠。放在低信噪比下会引入判决误差，但还是能给出一个可用的初值，真正的系统里应该在帧同步头之后做定时恢复。
4.3 多相结构替代直接插值/抽取
直接在高采样率上用长FIR滤波，计算量很大。上采样 L 倍时，插零后大部分样本为零，长滤波器的高样本率卷积存在明显浪费。多相结构把原型低通按相位分解成 L 个子滤波器，每个输入样本只触发一个子滤波器，计算量从 O(N*L) 降到 O(N)。下面的代码示意了上采样多相输出的等价写法：
function y = polyphaseUp(x, L, h) N = length(h); % 将滤波器系数按相位排列成 L 行 hp = reshape([h zeros(1, ceil(N/L)*L - N)], L, []); outLen = length(x) * L; y = zeros(1, outLen); for k = 1:L hk = hp(k, :); conv_part = conv(x, hk); y(k:L:end) = y(k:L:end) + conv_part(1:ceil(outLen/L)); end end这段代码不是性能最优的多相实现，但能把结构表达清楚：输入序列分别与 L 个相位滤波器卷积，再交错叠加到输出。实际工程中，MATLAB 内置的resample(x, L, M)就是有理数采样率转换的多相实现，内部会自动设计滤波器并做延迟补偿。替换自定义函数时，一定要重新验证采样点对齐，因为 resample 的相位约定和手写版本可能不同。
4.4 验证方法：EVM、误码率、频谱观测
链路调不通时，我按三个维度做验证：
- 频谱观测：用 pwelch 看插零后频谱镜像是否被抑制，下采样后是否还有折叠分量。如果镜像或混叠残留低于 40dB，在要求较高的链路里会有可测的EVM恶化。
- 星座图与 EVM：BPSK 只有两个点，但 EVM 能反映幅度和相位误差。
comm.EVM可以直接用，也可以按 4.2 节的欧式距离手算。 - 误码率曲线：理论 BPSK 误码率是 Q(sqrt(2*EbN0))。如果仿真曲线在高信噪比出现平层，优先查采样相位偏移；如果一开始就差 3dB，查滤波器增益或截止频率。
| 检查项 | 工具/方法 | 期望结果 |
|---|---|---|
| 镜像抑制 | pwelch(UpSampling输出) | 主谱与镜像间隔 > 40dB |
| 混叠抑制 | pwelch(DownSampling输出) | 无折叠分量 |
| 相位对齐 | 扫描M个抽取相位 | 最佳相位EVM最低 |
| 星座点幅度 | plot(real, imag) | 各点接近 ±1，散点小 |
5. 在MATLAB里做采样率转换模块：延迟补偿与冲激响应验证
5.1 用冲激响应检查增益和延迟
任何把 UpSampling.m 和 DownSampling.m 组合进新链路的时候，我都会先做一次冲激响应测试，再跑随机比特。冲激输入能一次性暴露滤波器增益、群延迟和抽取相位三个问题，比用随机数据看波形高效得多。
L = 4; M = 2; x = [1 zeros(1, 127)]; % 单位冲激 [y1, hUp] = UpSampling(x, L, 'fir'); % 上采样 [y2, hDown] = DownSampling(y1, M); % 下采样 % 上采样域的总群延迟近似值 delay_up = (length(hUp)-1)/2; delay_down = (length(hDown)-1)/2; total_delay = delay_up + delay_down; % 实际冲激中心 [~, actual_center] = max(abs(y2)); expected_center = floor(total_delay / M); disp(['expected: ', num2str(expected_center), ... ', actual: ', num2str(actual_center)]);理想情况下actual_center和expected_center应该一致，或者至少相差小于 1。如果差异很大，说明 DownSampling 的抽取相位选得不对，或者 UpSampling 里滤波器的延迟没有按设想对齐。更稳妥的做法是直接用finddelay(y2, x_expected)估计整体延迟，然后用circshift做序列对齐。
另一个实用检查是能量对比：输入冲激经过线性滤波和抽取后，输出能量应与输入能量保持同等量级。如果输出能量明显变小，优先检查 UpSampling.m 里有没有乘 L；如果输出能量正常但星座点发散，检查防混叠滤波器截止是否过低，把带内信号也削掉了。
把冲激测试、偏移搜索和BER统计放进同一个回归脚本，每次修改滤波器阶数或采样倍数都先跑一遍。例如在脚本入口保留这样一段：
% 回归入口: 修改L/M或滤波器后，先跑冲激再跑BER % runTest_srate_conversion('L', 4, 'M', 2, 'EbN0', 10);跑通这个回归脚本后，再去修改发射机里的成形滤波器参数，你会更容易判断性能变化到底来自脉冲成形，还是来自采样率转换。等到从纯仿真过渡到 AD9361 实际链路时，采样率转换这一块的数值行为和预期一致，剩下的就是射频前端校准和定时同步的对接了。
本文还有配套的精品资源，点击获取
