MATLAB直扩通信仿真：PN码、BPSK载波与误码率曲线避坑

简介：面向通信原理课程学习与扩频通信入门的MATLAB仿真设计文档，围绕直接序列扩频系统的建模、解扩与抗干扰性能评估展开，适合电子信息、通信工程专业学生及需要完成课程设计、实验报告或仿真练习的读者参考。压缩包文件总数仅1个，为doc文档，约320KB，集中呈现实验原理、直扩系统组成、m序列与GOLD、WALSH、OVSF码对比、扩频通信优势，以及2kbps信源、扩频因子64、BPSK调制和高斯信道下的仿真参数与结果分析。正文还给出部分程序代码与误码率随信噪比变化的观察思路，可帮助读者快速定位系统框图、伪随机码选择、调制解调流程和性能曲线等关键内容，理解扩频、解扩、处理增益与抗干扰之间的关系，并据此调整参数复现实验、整理图表、撰写分析结论。已有142人学习，适合作为短篇参考资料。

1. MATLAB 直扩仿真跑不通？先看这份文档把 PN 码、载波和误码率串成什么

一份基于 MATLAB 的直接序列扩频（DSSS）通信系统仿真文档，通常夹在课程设计或毕设资料包里：前半部分讲直扩系统组成、m 序列与 GOLD 码、抗干扰原理，后半部分直接甩出一段close all; clear all;的程序和各种无标题坐标轴截图。它要解决的事其实很具体——把「处理增益」「码片速率」「解扩后信噪比」这些课本词，落到一条能扫出误码率瀑布曲线的代码上。适合通信课程设计、毕设选题、以及想快速验证抗干扰性能的工程师。但它不是拿过来按 F5 就能出图的东西：里面缺函数、参数前后打架、载波频率甚至高于采样率。这篇就按我实际重搭一遍的顺序讲。

2. 扩频链路搭建：从 2 kbps 信源到 64 倍码片率的 BPSK 载波

2.1 先把采样率、码片率和载波频率的账算平

直扩仿真翻车最多的不是算法，是时间轴。这份文档里Rb=0.25、Tb=1/Rb算出 4 秒，可正文写的是「信源速率 2 kbps」；fs=4000一边出现，另一边又用fc=4e4做载波。三个数互相对不上，跑出来的东西一定是噪声模样。我一般会先列一张参数自洽表，把所有频率定死，后面所有代码都从这张表推。

这张表的约束关系只有三条，记住就不会乱：fs = Rb × SF × sps；fc < fs/2，实际工程里一般取fc ≥ 2 × Rc让载波甩开基带；Npb必须是整数，否则积分判决按比特切窗时最后一位会对不齐。文档里fc=4e4配fs=4000，载波是采样率的十倍，混叠之后解调出来只剩噪声——这是那份代码最致命的一处，第 5 章还会细说。

2.2 信源与 PN 码扩频：统一用 ±1 表示省掉一半麻烦

信源生成我建议直接映射成 ±1，而不是保留 0/1 再到处sign()。扩频本质是模二加，在 ±1 域上等价于乘法，扩频、解扩、判决三处都能用同一套乘法，代码短且不容易在取反的地方错一位。

rng(2024); % 固定随机种子，保证每次复现一致 Nb = 200; % 仿真比特数 bits = randi([0 1], 1, Nb); % 原始 0/1 比特流 sym = 2*bits - 1; % 映射为 -1 / +1，后续全部在此域运算 pn = mseq(6, [1 6], [1 1 1 1 1 1]); % 6 级 m 序列，周期 2^6-1 = 63 pn = [pn, -1]; % 补一位凑成 64，与扩频因子对齐 sps = 8; % 每码片 8 个采样点 chip_wave = repelem(pn, sps); % 把码片序列展开成采样级波形，长度 512 spread = kron(sym, chip_wave); % 每比特复制一遍 512 点扩频波形

mseq是自己补的函数，第 4 章给实现，这里只要知道它输出 ±1 且长度 63。repelem负责把离散码片变成带采样点的方波，kron则是把每个比特的符号乘到整段扩频波形上——kron比写 for 循环快，也不会在索引边界上出错。

提示：m 序列周期是 2^n−1，永远是奇数。扩频因子取 64 这种 2 的幂时，一定要显式补齐或截断，否则kron出来的长度对不上fs/Rb，后面切窗就全错。

补位方式有两种：把最后一位复制一遍（我上面这样），或者干脆把扩频因子定成 63，让Rc = 63 × Rb。前者对处理增益几乎无影响，后者参数更干净。我倾向后者，因为一旦 SF 是 63 这种非整幂，参数表里每个数都能互相整除，排查时少一类玄学问题。

2.3 BPSK 载波调制：文档里那段 I/Q 拆分其实是 QPSK 的骨架

文档里有一段很迷惑的代码：先reshape(d,2,N/2)拆出 I、Q 两路，再用I.*cos(...) - Q.*sin(...)合成，但坐标轴标题写的是「BPSK 信号」，后面解调又只做了一路相干相乘。这是典型的拿 QPSK 模板改 BPSK 改了一半。真正的 BPSK 只需要一路：

t = (0:length(spread)-1) / fs; % 时间轴，单位秒 carrier = cos(2*pi*fc*t); % 载波，必须满足 fc < fs/2 tx = spread .* carrier; % BPSK 调制：符号直接乘载波

这段的逻辑很直白：spread是 ±1 的基带波形，乘以载波后，+1 对应 0 相位，−1 对应 π 相位，这就是 BPSK。幅度保持 1，功率归一化后是 0.5 W，awgn函数用'measured'模式时会自动测量，所以不需要额外缩放。

如果要复现文档里那套 I/Q 结构，就得把fc处的正交两路都用上，那已经是 QPSK 的范畴，扩频因子、判决门限、误码率公式全都要重推。课设里图省事直接照抄那段，结果 BER 曲线跟理论差一个数量级，多半就是这里：I/Q 只解了一路，相当于丢掉一半信号能量。

3. 解扩与判决：相干解调后那三级滤波到底在滤什么

3.1 相干解调的同步假设与载波相乘

扩频通信的接收端要先解扩再解调，还是先解调再解扩？直扩系统里两种都有，这份文档走的是「先载波相干解调下变频到基带，再乘本地 PN 码解扩」的路线，工程上叫基带解扩。它成立的前提是收发载波同频同相——仿真里可以直接用同一个carrier变量，现实中要靠载波恢复环。

snr_dB = -12; % 输入信噪比，注意是负值，理由见第 6 章 rx = awgn(tx, snr_dB, 'measured'); baseband = rx .* carrier * 2; % 相干解调，乘 2 补偿 cos 相乘的 1/2 因子

.* carrier * 2这一步把射频搬回基带，同时产生一个 2fc 的高频分量，靠后面的低通滤掉。乘 2 是因为cos·cos = 0.5 + 0.5cos(2ωt)，不补这个系数，解调后幅度只剩一半，判决门限虽然还是 0，但跟噪声的相对关系会看着别扭。

注意：仿真里收发共用carrier是理想假设。真实链路里载波相位差几十度，BPSK 的判决就会开始翻转。要验证相位敏感度，可以手动给接收载波加一个phase_err再乘，看 BER 在多少度时开始崩，这个实验比看曲线有意思。

3.2 低通滤波的截止频率怎么定

低通的任务只有一个：把 2fc 分量干掉，同时尽量不伤基带方波。截止频率取4 × Rb是比较稳的选择——比信息带宽宽一点，方波的边沿保留得住，2fc 又远在通带外。

fc_lp = 4 * Rb; % 截止频率 800 Hz wc = fc_lp / (fs/2); % fir1 要求归一化到奈奎斯特频率 b = fir1(64, wc); % 64 阶线性相位 FIR y = filtfilt(b, 1, baseband); % 零相位滤波，避免群延迟让判决窗偏移

三个参数各有用处。阶数 64 决定了过渡带宽度，阶数太低 2fc 压不干净，太高则运算慢、边界瞬态长；归一化截止wc必须除以fs/2，写成fc_lp/fs是最常见的低级错误，滤出来会是一片平滑的直线；filtfilt做正反两遍滤波，相位抵消，代价是边缘要丢一段，所以后面切窗时我从头部对齐开始取。

文档里用了一个自定义的lpf(f, rf, fc/10000)，把截止频率写成fc/10000这种量纲不明的形式，换一组参数就失效。自己重写时把截止频率和采样率的关系写明白，比什么都强。

3.3 积分清洗与硬判决的判决窗

解扩就是再乘一次本地 PN 码，然后按比特周期做积分清洗。这一步等价于匹配滤波，把 512 个采样点的能量累加起来，噪声按 √N 增长、信号按 N 增长，信噪比改善 10·log10(512) ≈ 27 dB。

despread = y .* chip_wave; % 与发端完全相同的扩频波形相乘 Npb = fs / Rb; % 每比特 512 点 despread = despread(1 : Nb*Npb); % 截齐，丢掉滤波边缘瞬态 blk = reshape(despread, Npb, Nb); % 每列是一个比特的采样点 metric = sum(blk, 1); % 按列求和，积分清洗 rx_bits = metric > 0; % 硬判决 ber = mean(rx_bits ~= bits); % 误码率

reshape按列切是最省事的写法，前提是长度严格等于Nb*Npb，多一个点少一个点都会报维度错，所以前面那行截齐不能省。metric > 0的判决门限是 0，因为信号是 ±1 对称的；如果前面某处引入了直流偏置，门限就要改成mean(metric)，这也是排查「误码率恒定 0.5」时第一个要看的地方。

第 3 章这三步串起来，链路就闭环了。剩下的问题是文档里那几个被调用的函数根本不存在——mseq、goldseq、T2F、lpf全是黑的，不补上，代码连语法都过不去。

4. 文档里没给的四个函数：mseq、goldseq、T2F、lpf 自己补

4.1 m 序列：抽头、初相与周期

m 序列由线性反馈移位寄存器产生，输出长度 2^n−1，抽头位置对应本原多项式。文档里stage=3, ptap1=[1 3]，三级寄存器只有 7 位长，拿来当扩频码明显太短，处理增益不到 9 dB，曲线根本压不下来。我一般用 6 级或 7 级，抽头按常用的本原多项式取。

function seq = mseq(stage, taps, init) % stage : 寄存器级数，序列周期为 2^stage - 1 % taps : 反馈抽头位置（1-based，按本原多项式取） % init : 初相，不能全 0，否则退化为全 0 序列 N = 2^stage - 1; reg = init(:).'; % 行向量，末位为输出端 seq = zeros(1, N); for n = 1:N seq(n) = reg(end); % 每拍从末位输出 fb = mod(sum(reg(taps)), 2); % 抽头异或（模 2 求和） reg = [fb, reg(1:end-1)]; % 右移并注入反馈位 end seq = 2*seq - 1; % 0/1 映射到 -1/+1 end

6 级取taps = [1 6]，对应 x^6 + x + 1，周期 63。初相千万别写全 0，全 0 状态下反馈恒为 0，输出一串 0 映射成 −1，看着像扩频码其实是直流。判断方法很简单：sum(seq)的绝对值应该等于 1（m 序列比 0 多一个），如果不等于 1，说明初相或抽头取错了。

4.2 GOLD 序列与优选对

文档里提到 GOLD 码是「两个等长 m 序列优选对模二加」，这句话对，但落到代码上最麻烦的是「优选对」三个字——不是任意两个 m 序列都能配。工程上通常查表拿现成的抽头对，6 级常用的两组是[1 6]和[2 5]。

function code = goldseq(stage, taps1, taps2, init, user) % 两个 m 序列模 2 加得到 GOLD 码，user 为循环移位量，用于产生不同地址码 m1 = mseq(stage, taps1, init); m2 = mseq(stage, taps2, init); base = (m1 .* m2); % ±1 域上的乘法等价于模 2 加 code = circshift(base, [0, user-1]); % 循环移位产生同族不同码字 end

注意在 ±1 域上「模二加」变成乘法，这一步很多人卡住：mod(a+b,2)出来还是 0/1，跟 m 序列的 ±1 输出混在一起就乱了。统一在 ±1 域做乘法，结果自然也是 ±1。

至于 WALSH 和 OVSF，文档里只是原理性介绍，没有对应的仿真段，硬要补就得引入 Hadamard 矩阵生成和多速率正交码树，篇幅会膨胀。课设级别的仿真用 m 序列或 GOLD 码足够说明问题，把这两种吃透比四种都浅尝一遍更值。

4.3 T2F 与 lpf：能自己写就别依赖工具箱

文档里[f,rf]=T2F(td,jk_code)这种调用，T2F显然是时域转频域的自定义函数。实现很短，但有两个细节容易错：频率轴要fftshift对中，幅度要除以点数。

function [f, X] = T2F(t, x) % 时域转频域，返回双边谱，频率轴以 0 为中心 N = length(x); dt = t(2) - t(1); X = fftshift(fft(x)) / N; % 除以 N 使幅度与原始信号对应 f = (-N/2 : N/2-1) / (N*dt); % 频率分辨率 1/(N*dt) end

lpf就更好办，直接用fir1加filtfilt替代，不必自己造轮子：

function y = lpf(x, fc, fs, order) % 零相位低通滤波，fc 为截止频率（Hz），order 为 FIR 阶数 wc = fc / (fs/2); % 归一化到奈奎斯特频率 b = fir1(order, wc); y = filtfilt(b, 1, x); end

提示：fir1和filtfilt都在信号处理工具箱里。如果目标机器上没装，用conv(x, ones(1,K)/K, 'same')做移动平均也能凑合，代价是通带不够平、过渡带很宽，误码率曲线会比理论值差 1~2 dB，写报告时如实说明就行。

四个函数补完，整条链路才算真正能跑。我建议把这几个函数单独存成.m文件放在同一目录，主脚本用clear开头但不clear all，避免把函数也清掉。

5. 避坑与排查：这份代码最常翻车的五个地方

5.1 载波频率比采样率还高，解调出来全是噪声

现象：程序不报错，但图 5「扩频解调信号」和后面所有波形都是一片毛刺，误码率恒定在 0.5 附近。

原因：文档里fc=4e4而fs=4000，载波频率是采样率的 10 倍，离散化之后发生严重混叠，cos(2*pi*fc*t)采样出来的序列已经不是那个频率的正弦了，相干相乘当然解不出东西。

解决：把约束fc < fs/2写进参数表，且留足余量。本文取fc = 2 × Rc = 25600、fs = 102400，载波只有采样率的四分之一，波形上每个载波周期 4 个采样点，肉眼能看出正弦形状，这是最直观的自检方式——画出来的载波如果像噪声，后面的都不用看了。

5.2 Rb 的数值和文字描述对不上，时间轴全错

现象：坐标轴标着t(ms)，但波形拉伸得离谱，一个比特占了半张图；或者反过来，100 个比特挤成一条线。

原因：文档里Rb=0.25、Tb=1/Rb=4，单位是秒；正文却写「信源速率 2 kbps」。两个数差 8000 倍，linspace出来的时间轴自然对不上，axis([0,200,...])这类硬编码范围更是雪上加霜。

解决：所有时间一律用秒，只在画图时乘 1000 标成 ms。参数表里把Rb、Tb、fs、Npb四个数并排写出来，检查Npb = fs/Rb是不是整数，不是整数就调sps或Rb，直到整除为止。这一步花五分钟，能省后面两小时。

5.3 m 序列长度 63 配扩频因子 64，kron直接报维度错

现象：Error using kron, Vectors must have the same length或者更隐蔽的——不报错但误码率异常高。

原因：6 级 m 序列周期是 2^6−1=63，扩频因子写 64，repelem出来 504 点，reshape按 512 切必然错位，错一位之后所有比特的判决窗都偏了。

解决：要么补一位（pn = [pn, -1]），要么把 SF 改成 63。补位之后要重新验算Npb = 63 × 8 = 504，参数表跟着改。我一般选 SF=63，因为所有数都能整除，kron、reshape、sum三个地方都不用特殊处理。

5.4 SNR 扫描范围没扣处理增益，曲线是一条贴地的直线

现象：SNR = 1:1:15扫完，z1全是 0，画出来一条横线贴在坐标轴底部，看不出任何瀑布形状。

原因：扩频因子 64 带来 18 dB 处理增益，输入 SNR 到 1 dB 时，解扩后等效信噪比已经接近 19 dB，误码率早就低到仿不出来了。扫描区间整体偏高了 20 dB 以上。

解决：扫描范围按SNR_in = SNR_eff − Gp反推，想看到 10⁻¹ 到 10⁻⁴ 的过渡，等效 SNR 得在 0~10 dB，那输入就该从 −18 dB 扫到 −8 dB。本文用snr_dB = -20:1:0，曲线正好落在图中间。

5.5plot里的变量名和计算出来的不是同一个

现象：Undefined function or variable 'QPSK'，或者图能出来但曲线明显不对。

原因：文档里BPSK = (I.*cos(...) - Q.*sin(...))/sqrt(2)，后面却写plot(tq, QPSK)，而坐标轴标题又是「BPSK 信号」。这是从 QPSK 模板改过来时留下的残骸，变量名、图标题、注释三处互相矛盾。

解决：命名统一成「信号在链路中的位置」，比如tx_bb（基带）、tx_rf（射频）、rx_bb（解调后）、rx_ds（解扩后），每个名字自带位置信息，改名时不容易漏。MATLAB 里还可以开clearvars而不是clear all，让未定义变量在第一时间暴露出来而不是被上一轮运行的残留值掩盖。

6. 把误码率曲线跑成理论形状：一次参数扫描的完整做法

曲线跑出来只是开始，能不能跟理论对上才是这份仿真有没有价值的判据。BPSK 在 AWGN 下的理论误码率是Q(sqrt(2·Eb/N0))，直扩之后，解扩前的输入 SNR 与解扩后的等效 SNR 之间差一个处理增益，再加上每比特积分带来的过采样增益，实际关系是：

Gp = 10*log10(SF); % 扩频处理增益，SF=64 时约 18.06 dB Gos = 10*log10(sps); % 过采样积分增益，sps=8 时约 9.03 dB snr_eff = snr_dB + Gp + Gos; % 解扩积分后的等效信噪比 ber_th = qfunc(sqrt(2*10.^(snr_eff/10))); % 理论误码率

这里有个容易忽略的点：Gp和Gos不能混着算。扩频增益来自码片相关，过采样增益来自匹配滤波，两者物理来源不同，但在这套「先低通再积分」的实现里会叠加。我第一次对曲线时只加了 18 dB，实测比理论好 9 dB，查了半天以为是噪声功率算错，后来才反应过来是每码片 8 个采样点也被积分进去了。

做对照实验时建议分两组。第一组保留filtfilt低通，直接扫snr_dB = -20:1:0，看仿真曲线与ber_th的贴合度；第二组把低通那行注释掉，只做相干相乘加积分，此时等效增益变成Gp + Gos的完整形式，两条曲线应该几乎重合。如果第二组也贴不上，问题就回到 5.1 和 5.3——载波或码长对不上，跟噪声模型没关系。

再把扩频因子扫一遍，从 7、15、31 到 63，固定snr_dB = -15，看误码率随 SF 怎么掉。理论上 SF 每翻一倍，处理增益 +3 dB，曲线该往左平移 3 dB。实测里 31→63 那一段改善往往不到 3 dB，原因是 63 长度的 m 序列自相关旁瓣比 31 差，相关峰不够尖。这个偏差就是仿真能告诉你、课本不会告诉你的东西，写进报告比抄一段原理有分量得多。

最后给个复查习惯。从那以后我每次拿到这类课程设计文档，都先做一件事：把Rb、SF、sps、fs、fc五个数抄在纸上，逐个验算fs = Rb×SF×sps、fc < fs/2、Npb = fs/Rb是否为整数、SF是否为 2^n−1。四条全过再敲第一行代码，比跑通之后回头找 bug 省事得多。希望帮到你。

本文还有配套的精品资源，点击获取