MATLAB扩频仿真:m序列/Gold序列/Walsh序列原理与选型
简介面向通信工程研究人员与数字通信系统设计者这份资源聚焦M序列、Gold序列与Walsh序列的生成及性能仿真解决CDMA等场景下扩频码选择与抗干扰能力评估问题。压缩包共2个文件包含1个MATLAB主程序脚本和1个说明文档包体仅4KB脚本可直接运行并支持修改序列长度、信道模型、噪声类型等参数适合课程设计、毕业设计或通信系统预研。当前已有156人学习下载。配套内容围绕三种序列的自相关/互相关特性展开结合信噪比分析和误码率计算直观对比不同序列在多址访问与抗噪声方面的表现其中包含自相关函数与互相关函数的对比结果能清晰观察序列的旁瓣特性与正交性差异。说明文档对仿真流程和结果解读做了梳理可帮助读者快速复现实验并为实际系统设计提供选型与优化参考。1. M序列/Gold序列/Walsh序列数字通信仿真里三种扩频码先搞清楚再动手做CDMA或直扩DSSS链路仿真时扩频码往往是第一个翻车点。很多人图省事用randi随机生成一串数当扩频码结果跑自相关发现旁瓣高得离谱接收端相关积分根本提不出尖锐峰扩频增益被白白吃掉。这套MATLAB数字通信系统仿真资源解决的就是这个问题用m序列、Gold序列和Walsh序列三种经典扩频码把序列生成、相关性分析、DSSS基带链路误码率仿真串成一条完整验证链路。适合三类人通信原理课设要交仿真报告的学生、做CDMA或直扩系统预研的工程师、想快速搞懂扩频码选型逻辑的入门者。上手路径很直接——先跑自相关三条曲线再跑误码率对比几十行脚本就能看出三种码的性格差异。2. 扩频序列的原理与选型为什么CDMA仿真里这三种码缺一不可直扩系统的本质是用带宽换信噪比低速数据乘上一路高速扩频码码速率决定带宽扩展倍数码的统计特性决定系统抗多址干扰和码同步捕获的能力。m序列、Gold序列、Walsh序列是三种最常被拿来对比的候选码m序列胜在结构简单、自相关尖锐Gold序列胜在互相关有界Walsh序列胜在完全正交。但问题在于它们的天赋恰好也是软肋——选错码型后面捕获、解扩、误码率统计会连锁翻车。这章先把生成原理讲透再谈应用场景怎么选。2.1 m序列LFSR移位寄存器与本原多项式的关系m序列的全称是最大长度线性反馈移位寄存器序列。一个n级寄存器在反馈逻辑驱动下最多有2的n次方个状态其中全零状态会卡死所以能遍历其余2的n次方减1个状态的序列就是最大长度序列周期正好是2^n−1。反馈抽头组合能不能跑出最长周期取决于它对应的特征多项式是不是本原多项式。MATLAB通信工具箱可以直接查本原多项式列表例如n5时x^5x^21是本原的n7时x^7x^31是本原的这两个组合足够应付绝大多数课程设计和预研验证。function m mseq(n, taps, init, len) % n : 寄存器级数 % taps : 反馈抽头位置从1开始计数1为最低位 % init : 初始状态向量长度n不能全零 % len : 输出长度默认 2^n - 1 if nargin 4 len 2^n - 1; end reg init(:); m zeros(1, len); for k 1:len m(k) reg(n); % 最高位作为输出码片 fb mod(sum(reg(taps)), 2); % 抽头位置做模二加(异或) reg [fb, reg(1:n-1)]; % 反馈值进入最低位整体右移 end end这段代码有三个关键点。一是输出永远取最高位reg(n)换一个输出位会得到相位不同的同一条序列周期不变但码字起点变了后面解扩时本地码必须和发送端同一个输出位。二是反馈值是所有抽头位置的模二和抽头组合里必须包含第n位否则反馈环断开序列会退化成短周期。三是初始状态不能全零否则LFSR一直输出零生成Gold序列时两条m序列的init也不能相同。生成完成后用m(m0) -1转成双极性码扩频乘法直接用±1相乘比if判断快也更容易读。验证m序列生成是否正确最快的方法是算周期自相关主峰为1其余时延旁瓣恒为-1/N。如果旁瓣不是这个值先怀疑抽头配错或init不是非零状态不要急着往下跑误码率。2.2 Gold序列优选对异或生成互相关从不可控变成有界m序列自相关漂亮互相关却不可控。随便找两条同长m序列互相关峰值可能接近1两个用户同时发射时互相干扰严重这在CDMA场景里不可接受。Gold的解决办法很朴素找一对满足特定条件的m序列业内叫优选对preferred pair。优选对中任意一条序列与另一条循环移位后的序列互相关只取三个固定值最大值被压到一个可预测量的上界。把这两条m序列逐位异或就得到一条Gold序列第二条序列每次循环移位一次再异或就能得到一族周期相同、互相关都有界的Gold码码族大小为2^n1。function g gold_seq(m1, m2, shift) % m1, m2 : 长度相同的0/1 m序列必须是优选对 % shift : 对m2循环移位后再异或产生族内不同码 g xor(m1, circshift(m2, shift)); end调用时有三个注意点。两条m序列必须同长度、同周期否则异或结果不再具有Gold性质这在代码上表现为两个mseq返回向量长度不一致时直接报错或生成错误码。m2的移位量shift取值在0到2^n−1之间不同shift对应码族内不同用户码道。仿真前把0/1转成±1用g*2-1。我实际做预研时习惯先一次性生成五组不同shift的Gold码打印两两之间的互相关峰值确认都落进理论三值集合{-1, -t(n), t(n)-2}内再进主仿真避免后面误码率曲线异常时分不清是码的问题还是链路实现的问题。2.3 Walsh序列Hadamard正交码同步时完美、异步时脆弱Walsh序列不来自LFSR而是来自Hadamard矩阵。矩阵按递归规则构造1阶是[1]把当前矩阵H拼成[H H; H -H]就得到下一阶。2^k阶的Hadamard矩阵任意两行内积为零把每一行当一条码同一矩阵内任意两条码完全正交理论上完全消除多用户间干扰——前提是所有用户码片严格同步到达接收端。定时偏差一旦出现正交性立刻瓦解性能可能比m序列还差。function W walsh_mat(k) % k : 生成 2^k 阶 Hadamard 矩阵 W 1; for i 1:k W [W W; W -W]; % 递归扩展 end end % 取第 r 行作为 Walsh 码 r 2; walsh walsh_mat(3); % 8阶矩阵 walsh_code walsh(r, :); % 第2行长度8使用时有三个坑。第一扩频因子SF必须等于矩阵阶数2^kSF8就取8阶矩阵SF64取64阶不能拿8阶矩阵的码去配64倍扩频的链路。第二矩阵元素本来就是±1不要画蛇添足做0/1到双极性的转换。第三Walsh码的周期自相关特性很差非零时延下旁瓣并不低意味着它不能靠相关峰做码捕获实际系统里Walsh码通常配一个独立的同步码多用m序列先完成定时再用正交码区分用户。三种码的特征对比如下表选型时直接对着场景挑。序列类型生成结构长度自相关特性互相关特性同步要求典型用途m序列n级LFSR2^n−1主峰尖锐旁瓣-1/N不可控中等捕获、测距、同步头Gold序列两条优选m异或2^n−1尖锐旁瓣低三值有界中等CDMA用户码、抗多址Walsh序列Hadamard矩阵行2^k无尖锐峰旁瓣偏高同步时严格为零高前向正交信道化3. MATLAB仿真搭建从序列生成到误码率曲线一份能直接跑的代码流程这套资源里的主仿真脚本和我平时搭DSSS链路的口径基本一致先定参数后生成码再按“扩频—加噪—解扩—判决”四步走。资源里的脚本按模块拆开mseq.m、gold_seq.m、walsh_mat.m、误码率统计函数和主仿真脚本分开存放下载后按顺序调用就行。下面把每一步拆开讲代码可以直接抄进自己的工程换码型只需要换一行。3.1 序列生成模块抽头、初始状态和周期这三个参数先定死生成三种码之前先把寄存器级数n、抽头taps、初始状态init三个参数定死。这三个参数决定码的周期、相位和正交性。我推荐n7起步周期127足够看出相关特性跑误码率也快想验证多用户场景再换n10或n15但蒙特卡洛仿真会明显变慢不要一上来就用大周期。n 7; % 寄存器级数周期 2^7-1127 % 本原多项式 x^7 x^3 1抽头 [7 3] m1 mseq(n, [7 3], [1 0 0 0 0 0 0]); m2 mseq(n, [7 3], [0 1 0 0 0 0 0]); % 换 init得到另一条 m 序列 gold gold_seq(m1, m2, 5); % shift5得到一条 Gold 码 walsh walsh_mat(3); % 8 阶 Hadamard 矩阵 walsh_code8 walsh(2, :); % 取第 2 行做码长度 8 % 转双极性用于扩频乘法 sm 2*m1 - 1; sg 2*gold - 1;m1和m2用的是同一个本原多项式只是初始状态不同得到的m序列在周期和频谱特性上等价但相位不同这样异或出来的Gold序列才有意义。gold的移位量取了5也可以换成其他shift生成族内另外的码字。注意这里m序列长度127Walsh码长度8单用户AWGN场景下做误码率对比不需要强行把长度拉齐——因为Eb/N0定标已经把所有扩频增益归一到比特能量上码长只决定码片速率和占用带宽不改变AWGN下的理论误码率上界。这个结论在3.4节会再验证一次。3.2 DSSS扩频调制与解扩码相位对齐是解扩正确的前提扩频的数学操作很干净一个数据比特乘整条扩频码。用m序列时每个比特对应127个码片码片速率就是比特率乘127。发送端把数据逐比特与码相乘再串起来就得到基带扩频信号。SF 127; % 扩频因子等于码长 Nbits 2000; data 2*randi([0 1], 1, Nbits) - 1; % ±1 数据 tx []; for k 1:Nbits tx [tx, data(k)*sm]; % 每比特乘整条码 end加噪和定标是这里最常出错的地方。码片幅度设为1时一个比特的能量就是SF即EbSF。给定Eb/N0后噪声单边功率谱密度N0SF/EbN0_linear码片上噪声方差σ²N0/2所以sigma直接写成sqrt(SF/(2*EbN0))这个换算关系是后面误码率曲线能否贴住理论值的命根子。EbN0dB 0:2:10; ber zeros(size(EbN0dB)); for i 1:length(EbN0dB) EbN0 10^(EbN0dB(i)/10); % 线性值 sigma sqrt(SF/(2*EbN0)); % 码片噪声标准差 rx tx sigma*randn(1, length(tx)); rx_mat reshape(rx, SF, Nbits); % 按比特窗口切片 y sum(rx_mat .* repmat(sm, 1, Nbits), 1) / SF; d (y 0); % 硬判决 ber(i) sum(d ~ (data 0)) / Nbits; end解扩的关键在于reshape窗口与本地码完全对齐。rx_mat的第k列是第k个比特的SF个码片本地码sm必须是发送端同一个码、同一个相位逐列相乘再求和相当于做了一次相关积分。如果收发两端码相位错一位相关积分拿不到主峰误码率直接崩溃到0.5附近。正规流程里解扩之前要先做码捕获确认峰值时延再开始统计误码率这个我在第6章会给出一个完整示例。3.3 蒙特卡洛误码率统计Eb/N0、噪声方差与仿真次数的换算关系上面循环里每个信噪比点固定2000个比特低信噪比够用但到8dB以后错误比特只剩个位数曲线会剧烈抖动甚至出现误码为0的假象。我一般改成错误计数停止准则每个信噪比点累积至少100个错误比特再退出同时设一个最大比特数兜底防止高信噪比下while循环跑不完。这样低信噪比快速收敛、高信噪比平滑稳定。function ber mono_carlo(code, EbN0dB, maxbits) % code : 双极性扩频码 % 每个信噪比点至少统计 100 个错误比特 SF length(code); for i 1:length(EbN0dB) EbN0 10^(EbN0dB(i)/10); sigma sqrt(SF/(2*EbN0)); % 噪声定标 errs 0; bits 0; rng(42); % 固定随机种子方便复现 while errs 100 bits maxbits nbits 1000; data 2*randi([0 1], 1, nbits) - 1; tx kron(data, code); % 向量化扩频等价于逐比特循环 rx tx sigma*randn(1, length(tx)); y reshape(rx, SF, nbits); dec sum(y .* repmat(code, 1, nbits), 1) 0; errs errs sum(dec ~ (data 0)); bits bits nbits; end ber(i) errs / bits; end endkron在这里等价于前面逐比特扩频的循环但MATLAB向量化写法性能高一个量级仿真周期长时差别很明显。rng(42)固定随机种子是为了让每次仿真结果可复现交报告或者排查问题时能确定不是随机数差异造成的曲线漂移。100个错误数是误码率统计的工程惯例想更平滑可以提到200代价是仿真时间翻倍。3.4 运行结果判读先看自相关再看BER别上来就盯曲线主仿真跑完先别急着看三条误码率曲线的重合度按这个顺序自检。第一确认m序列周期自相关旁瓣在-1/127附近Gold互相关峰值落在理论范围内第二把仿真误码率和理论BPSK AWGN曲线画在一起偏差超过0.5dB就说明噪声定标或者解扩有问题第三再看三种码的结果是否彼此重合。eb 10.^(EbN0dB/10); th qfunc(sqrt(2*eb)); % BPSK AWGN 理论误码率 semilogy(EbN0dB, ber, o); hold on; semilogy(EbN0dB, th, -); legend(仿真,理论); grid on; xlabel(Eb/N0 (dB)); ylabel(BER);主仿真参数表如下跑完对着表核对一遍基本能定位八成的问题。参数本次取值说明寄存器级数 n7m序列周期127扩频因子 SF127 / 8每比特码片数信噪比范围0~10 dB步进2 dB共6个点每点最小错误数100控制曲线平滑度maxbits 兜底5e6防止高信噪比死循环三条曲线重合是正确的不是bug。单用户AWGN下扩频不带来信噪比增益扩频码的区别要放进多用户、异步、捕获场景才体现。所以看到重合说明链路实现是对的看不到重合反而要回头查代码。4. 三种序列的性能对比自相关、互相关与抗噪声能力差在哪里单用户AWGN场景三种码几乎没差别真实系统是多用户、带定时偏差的区分度全部来自码字本身的统计特性。这章用相关函数把三种码摊开看再上双用户干扰实验结论会非常直观。4.1 自相关与互相关的MATLAB实测峰值旁瓣怎么量化相关函数的量化指标就两个自相关峰值旁瓣、互相关峰值。旁瓣越低捕获越不容易错锁互相关峰值越低多用户干扰越小。function R corr_cyclic(a, b) % 周期互相关结果归一化 % ab 时得到自相关 N length(a); R zeros(1, N); for tau 0:N-1 R(tau1) sum(a .* circshift(b, tau)) / N; end end % 实测 Ra corr_cyclic(sm, sm); side max(abs(Ra(2:end))); % 排除主峰 fprintf(m序列自相关旁瓣: %g\n, side); % 理论 -1/127 Rg corr_cyclic(sg, sm); fprintf(Gold互相关峰值: %g\n, max(abs(Rg))); % 理论上界约17/127对n7的m序列周期自相关旁瓣恒为-1/127约-0.0079Gold互相关峰值被压在约17/127约0.134这就是它的三值有界特性。Walsh码同步时互相关严格为零但自相关在非零时延下会出现多个近等高旁瓣没有尖锐主峰这个特征直接判了它“不能用于捕获”的死刑。提示corr_cyclic的第二参数用circshift循环移位对应的是周期相关。实际系统中如果码片定时抖动超过一个码片非周期相关特性也要测但仿真验证阶段先看周期相关就够了。4.2 单用户AWGN下的BER对比三者理论一致差异在别处把3.3节的mono_carlo分别喂给三种码一行代码跑一组。ber_m mono_carlo(sm, EbN0dB, 5e6); ber_g mono_carlo(sg, EbN0dB, 5e6); ber_w mono_carlo(walsh_code8, EbN0dB, 5e6);跑完三条曲线全部落在理论BPSK曲线上。原因在于解扩是线性相关运算加性白噪声经过相关积分后统计特性不变扩频码只决定信号在时间轴上的形态不改变信噪比。所以单用户AWGN场景下“选哪种码都行”。这里给一句直白结论单用户BER仿真验证的是链路实现正确性不是码型优劣码型的价值要到多用户和捕获场景才体现。很多人第一次跑出重合曲线以为出了问题实际上这是正确的、且必须出现的结果。4.3 双用户干扰场景Walsh正交优势与m/Gold的有界互相关把链路扩展成两个用户用户2的信号对用户1构成多址干扰。两个用户功率相等、定时对齐时用户1解扩输出里除了自己的数据和噪声还有用户2的贡献贡献大小由两条码的互相关决定。% 双用户CDMA用户2用Gold码sg Nbits 2000; data2 2*randi([0 1],1,Nbits)-1; tx2 kron(data2, sg); % 用户2扩频信号 rx1 tx1 tx2 noise; % 接收叠加 y1 reshape(rx1, SF, Nbits); dec1 sum(y1 .* repmat(sm,1,Nbits), 1) 0; % 干扰系数 I sum(sm .* sg) / SF; % 用户间干扰归一化系数当两用户定时对齐、功率相等时解扩后用户1的输出里包含data2 * I这一项。Walsh码场景下I严格为零多址干扰被完全消除这是IS-95前向链路用Walsh码的根本原因。m/Gold场景下I是循环互相关值有界但非零用户数增多时干扰功率按用户数累加所以码型选择必须结合用户数和同步质量。异步场景下Walsh的正交性完全破坏I不再为零干扰反而比Gold更猛因为它的自相关旁瓣没有规律可预判。选型结论可以压缩成一张表。场景推荐码型原因同步多用户前向CDMAWalsh互相关严格为零异步多用户 / 随机接入Gold互相关有界、码族大同步头 / 捕获 / 测距m序列自相关主峰尖锐5. 避坑指南扩频序列仿真里最常见的五个翻车点这套仿真代码我拆过好几版也帮人排查过不少报错。下面五个问题是学生和入门工程师最容易撞上的统一按现象、原因、解决记录每一条都是真实踩过的坑。5.1 现象Gold序列生成出来全是零原因有两个可能一是两条m序列的初始状态相同异或后每一位都是0二是mseq函数内部初始状态设为全零LFSR进入死锁输出恒为0。前者是使用问题后者是参数问题。解决生成后立刻检查sum(gold)和unique(gold)全零立刻能看出来。确认m1和m2的init不同且都不是全零向量换一组shift重试如果还是全零单独打印m1和m2的前32位看两条序列是否完全相同。5.2 现象解扩以后误码率比不扩频还高这是最迷惑人的一个坑。原因通常是三选一本地码与发送码相位没对齐相关积分取不到主峰reshape切窗口时SF用错积分窗口错位或者是噪声定标重复既按Eb/N0加了噪声又在判决门限里二次折算噪声功率。解决先跑一个开环验证——在无噪声条件下解扩输出应该是干净的数据±1有噪声时先算接收信号与本地码的循环互相关找到峰值对应的时延再解扩。如果BER曲线整体比理论高但不发散优先怀疑sigma公式里SF的位置写错。5.3 现象Walsh码维数错位正交性直接消失Hadamard矩阵阶数2^k必须等于扩频因子。常见错误是SF64却拿了8阶矩阵的行或者反过来正交性在数学上就不成立多用户干扰实验完全失真。解决用walsh_mat(3)取8阶矩阵时链路里SF同步改成8。取行后做一次正交性自检sum(walsh_code8 .* walsh(3,:))应该等于0整矩阵验证用W*W/8 eye(8)。这两行自检代码一分钟就能跑完能拦住后面所有误判。5.4 现象高信噪比下误码率曲线抖动剧烈或出现零误码原因就是每个信噪比点固定仿真比特数8dB以上实际误码率降到1e-4以下2000个比特里可能一个错误都没有统计上完全失效。解决改用3.3节的错误计数停止准则每个点至少攒100个错误比特或者改用半解析方法只仿真信道衰落系数用qfunc直接算条件误码率再取平均。仿真时间从几分钟拉到几十分钟时优先检查是不是还在用固定比特数方案。5.5 现象m序列周期对不上理论值2^n−1抽头多项式不是本原多项式序列周期会缩短或者taps索引把1基和0基搞反反馈抽头取错位置再或者init全零进入死锁。MATLAB里最隐蔽的是多项式对但taps写反比如x^7x^31对应[7 3]错写成[3 7]实际抽头位置变了。解决先用通信工具箱的本原多项式表核对比如gfprimfd(7, min)可以列出7阶本原多项式然后打印输出序列数周期长度最后用周期自相关旁瓣-1/N做校验三项全过再往下走。我一般把这三项检查写成一个自检脚本每次换多项式都跑一遍。6. 进阶验证用滑动相关捕获器把三种序列的同步性能一次测清楚前面的BER仿真默认码相位已经对齐但真实接收机第一个任务是“不知道码从哪里开始”。这一步可以用滑动相关捕获器验证三种码的同步能力也是判断码型选得对不对的最直接实验。做法是发送端对扩频信号做一个人为时延接收端用本地码在滑动窗口里算相关能量峰值对应的位置就是估计时延。delay 50; % 发送端人为时延 rx circshift(tx, delay); % 接收信号 code sm; % 换 sg / walsh_code8 对比 SF length(code); metric zeros(1, SF); for tau 0:SF-1 seg rx(1:SF); % 取一个比特窗口 metric(tau1) abs(sum(seg .* circshift(code, tau))); end [~, est] max(metric); fprintf(估计时延%d理论%d\n, est-1, delay);对m序列和Gold序列相关谱会出现唯一尖锐主峰估计时延与理论值完全一致这是它们能承担同步头职责的根本原因。Walsh码则会出现多个近等高的旁瓣峰最大峰值的位置随机漂移捕获器容易锁错时延——这就是为什么纯Walsh系统必须额外配一路同步码。把三种码分别跑一遍这个实验比看十页理论推导都直观。注意滑动相关捕获的窗口必须正好等于一个扩频周期窗口小了相关积分能量不足窗口大了会把相邻比特数据混进来。从那次以后我每次跑扩频仿真都强制先走一遍“相关性三连”自相关旁瓣、互相关峰值、捕获谱线确认码本身没问题再进误码率统计。这个习惯帮我堵住了大量本来要浪费一整天的排查。这套资源里的脚本从序列生成到捕获验证都齐下载后按章节顺序跑每步输出的都是可复现的数值希望帮到你。本文还有配套的精品资源点击获取