Matlab中北斗B1I、B1C、B2a信号仿真：码生成、调制与验证

简介基于Matlab的北斗三号信号模拟程序包面向通信、卫星导航方向的本科生与研究生用于学习B1I、B1C、B2a三种公开服务信号的生成方法与频谱特性可服务于课程设计、毕业设计或科研预研。压缩包内共6个M文件大小仅5KB代码轻量且模块划分清晰包含BDS主程序、B1C/B1I/B2a数据与导频分支生成脚本以及legendre序列生成工具便于按需修改参数、替换调制方式并二次开发。目前已有164人学习下载。通过运行Matlab仿真可获得三种北斗信号的时域波形与码序列设计理解BOC、QMBOC等调制方式的实现原理与信号差异同时还能借代码注释和运行结果快速定位参数设置问题免去从零搭建的麻烦。整套代码适合希望尽快入门北斗信号级仿真的中初级学习者也可作为课堂教学的补充演示材料。1. 为什么在 Matlab 里模拟 B1I、B1C 和 B2a 三种北斗信号北斗三代的民用信号分成三条典型链路B1I 从北斗二号沿用至今工作在 1561.098 MHzB1C 落在 1575.42 MHz和 GPS L1 同频B2a 在 1176.45 MHz与 Galileo E5a、GPS L5 共频。这三者恰好覆盖了接收机研发里最常遇到的三种情况兼容老协议、跨系统互操作、宽带高精度测距。射频前端工程师要判断这三个频点能不能共用一个采样通道基带算法工程师要算清 BOC 主瓣内的功率损失捕获工程师则需要一段没有多径、没有干扰的干净波形来调门限和积分时长。Matlab 的向量化能力适合把码表、副载波、载波合成按层级拆开处理不必一上来就上 FPGA。这种仿真工程通常不大但坑集中在码表来源、副载波相位和采样率匹配上下面按从理论到落地的顺序把关键步骤和参数逐一拆开。2. 三种北斗信号的码结构差异与主码生成实现2.1 码率、码长、周期和调制方式对照写任何代码之前先把参数表立起来后面所有上采样倍数、积分周期、文件存储估算都从这张表推出来。三个信号的主码周期差异非常大直接决定捕获模块的搜索策略。信号中心频率 (MHz)码速率 (Mcps)码长 (码片)主码周期调制方式B1I1561.0982.04620461 msBPSK-R(2)D1 电文 50 bpsB1C1575.421.0231023010 ms数据分量 BOC(1,1)导频分量 QMBOC(6,1,4/33)B2a1176.4510.23102301 ms数据/导频均为 BPSK-R(10)B-CNAV2 50 bps读这张表有一个容易忽略的结论B1C 的码速率只有 B1I 的一半但码周期是 B1I 的十倍所以捕获时的多普勒搜索网格要比 B1I 密得多。再看存储量B2a 一毫秒就要落 10230 个码片按每码片 6 个采样计算100 毫秒就是 6.138M 个复数点。仿真时长、通道数和数据精度没有定下来之前不要先写调制代码否则后面会因为内存不够反复返工。2.2 B1I 主码2046 码片的 Gold 码与抽头表B1I 的测距码由两个 11 级移位寄存器按 Gold 码思路构造主码长度是 2046 而不是 GPS L1 的 1023原因是北斗规定在 1 ms 内播发 2046 码片速率为 2.046 Mcps。工程实现里最常见的是按 ICD 中每颗卫星的相位抽头表逐星替换生成参数。function code bds_b1i_prn(prn) % B1I 主码生成演示结构按 ICD 思路实现 % 两个 11 级移位寄存器初始状态全 1 g1 ones(1, 11); g2 ones(1, 11); nChip 2046; code zeros(1, nChip); % 演示用相位抽头真实工程必须按 ICD 表逐卫星替换 g2Phase g2DelayTable(prn); % 反馈抽头示例正式使用以 ICD 信号空间接口文件为准 g1TapIdx [1 7 8 9 10 11]; g2TapIdx [1 6 8 9 10 11]; for k 1:nChip code(k) xor(g1(g1TapIdx(1)), g2(g2Phase)); fb1 mod(sum(g1(g1TapIdx)), 2); fb2 mod(sum(g2(g2TapIdx)), 2); g1 [fb1, g1(1:10)]; g2 [fb2, g2(1:10)]; end code 1 - 2*code; % 0/1 映射为 1/-1 end移位寄存器用行向量表示反馈位取对应抽头位置的异或和每循环一次整体右移一位循环 2046 次得到完整主码。g2Phase决定本卫星取自 G2 寄存器的哪一级不同卫星用不同相位保证码间互相关性足够低。如果抽头位置填错自相关主峰会变钝捕获时会出现多个近等高的旁瓣这是调试中最难定位的一类问题。很多现成工程干脆不写抽头表而是把 ICD 附录里的码序列文本解析成int8矩阵启动时一次性读入内存。对演示工程来说这种做法更稳妥。2.3 B1C 与 B2a 的 10230 码片主码查表生成更可靠B1C 和 B2a 的主码长度都是 10230文档里统称为分层码。理论上可以用两个十级移位寄存器实时递推但国标和 ICD 都给定了每颗卫星的明确码序列工程上直接在初始化阶段读取mat或文本码表即可避免长周期递推出现相位漂移。function code b1c_primary(prn) % 读取 B1C 主码prn 范围按 ICD 取 1~63 % 返回 1x10230 的 1/-1 向量 persistent cache if isempty(cache) S load(bds_b1c_primary.mat); % 63x10230 int8 表 cache double(S.codes); end code cache(prn, :); endB2a 的读表函数结构完全相同只需要把数据源换成bds_b2a_primary.mat。用persistent缓存码表可以避免每次调用都产生磁盘 IO否则在逐毫秒循环的捕获仿真里文件读取时间会超过信号生成本身。查表法的内存代价很低63 颗卫星乘 10230 码片约 0.63 MB完全可接受。2.4 多帧扩频序列的边界与数据位对齐信号仿真通常要连续生成 20 到 100 毫秒的发射波形最容易出问题的是数据位与码边界的对齐关系。B1I 的 D1 电文速率 50 bps一个数据比特占 20 个码周期码序列的重复必须和导航位符号保持同步不能把导航位按逐码片独立随机来做。一个正确示例numFrames 20; % 20 个 B1I 主码周期对应 1 个数据比特 chipSeq repmat(code, 1, numFrames); bitSeq repmat(2*randi([0 1])-1, 1, length(chipSeq)); b1iChips bitSeq .* chipSeq;repmat先重复单个主码周期bitSeq用常数常量扩到相同长度乘起来得到码片级发射序列。如果数据位翻转落在任意码片中间后面的相干积分峰值会掉 3 dB 左右表现就是捕获图上某颗星明显比相邻星弱。数据位必须按 20 ms 粗粒度对齐这是仿真发射端和真实接收机最容易被同时忽略的细节。3. Matlab 模拟链路BPSK-R、BOC(1,1) 和 QMBOC 调制合成3.1 三条支路的公共处理流程信号模拟顺序按“码表 → 数据/导频分组 → 副载波调制 → 带宽整形 → 载波混频 → 多信号合成”来组织。建议一开始只做基带复包络先验证码对齐再把载波项加进去。每条支路独立生成和验证后再合路可以减少三路信号混在一起时难以定位的相位问题。3.2 B1I 基带信号码片与数据符号直接相乘fs 4.092e6; % 基带采样率取码速率 2 倍 spc 2; % 每码片 2 个采样点 code bds_b1i_prn(1); navBit 2*randi([0 1]) - 1; chip_vec navBit * code; bb_b1i reshape(repmat(chip_vec, spc, 1), [], 1).;repmat按列展开后再整形就把每个码片复制成spc个采样点。基带采样率 4.092 MHz 对 B1I 主瓣 4.09 MHz 带宽刚好覆盖做码相关的功能验证足够。真正带射频前端效应时需要增加滤波器和采样率h_b1i fir1(64, 2.3e6/(fs/2)); bb_b1i_f filter(h_b1i, 1, bb_b1i);通带边缘取 2.3 MHz 是为了保留主瓣的绝大部分能量同时抑制带外噪声放大。滤波器阶数偏低时过渡带会拖尾后续加多普勒和干扰信号时频域形状不够真实。3.3 B1C 的 BOC(1,1) 与 QMBOC 副载波实现B1C 比 B1I 多一层副载波数据分量用 BOC(1,1)导频分量用 QMBOC(6,1,4/33)。BOC 的直观含义是把普通码片和一个更高频的方波相乘每个码片被切成 2n 个小段正负交替的方波在频谱上形成双峰结构。先在码片级生成副载波序列function sc bocSubcarrier(n, spc) % n: BOC(n,m) 中的 n 值 % spc: 每个码片的采样数必须能被 2n 整除 segLen spc / (2*n); seg [ones(1, segLen), -ones(1, segLen)]; sc repmat(seg, 1, n*2); end对于 BOC(1,1)取spc4时副载波序列是[1 1 -1 -1]对 BOC(6,1)取spc12副载波会切成 12 小段。为了让两个分量共用同一套上采样参数通常直接选能同时被 2 和 12 整除的spc比如 12 或 24。接下来合成 B1C 数据与导频通道spc 12; b1cCode b1c_primary(4); dataSeq 2*randi([0 1], 1, 10230) - 1; pilotSeq 2*randi([0 1], 1, 10230) - 1; sc11 bocSubcarrier(1, spc); sc61 bocSubcarrier(6, spc); sigData dataSeq .* sc11; sigPilot pilotSeq .* (sqrt(4/33)*sc11 sqrt(29/33)*sc61); bb_b1c sigData sigPilot;这里sqrt(4/33)和sqrt(29/33)对应 QMBOC 中 BOC(1,1) 与 BOC(6,1) 的功率分配比例两个系数的平方和为 1保证导频通道总功率归一。不同版本的 ICD 对这个比例的写法略有差异动手前要确认参照的是哪一版直接沿用旧代码的等权系数会让 B1C 频谱主瓣形状失真。3.4 B2a 信号与三路中频合成B2a 使用 BPSK-R(10)没有副载波码速率 10.23 Mcps码长 10230。基带主瓣宽度约 20.46 MHz采样率必须覆盖到这个宽度否则相关峰会被明显拉宽。b2aCode b2a_primary(6); b2aData 2*randi([0 1], 1, 10230) - 1; bb_b2a b2aCode .* b2aData; b2aUp reshape(repmat(bb_b2a, 6, 1), [], 1).; % 上采样到 61.38 MHz三路信号可以在同一个中频坐标系下合成。取公共采样率 61.38 MHz各支路先按各自的码速率上采到该采样率再乘不同的载波偏移后在时域相加。B1I 与 B1C 的频率差是 14.322 MHz在 61.38 MHz 采样率下可以有效排布而不混叠。偏移量加中频带宽只要不超过采样率一半合路结果就能直接作为软件接收机的输入文件。3.5 采样率、存储与典型参数表信号码速率每码片采样基带采样率100 ms 复数 single 存储B1I2.046 Mcps816.368 MHz约 13 MBB1C1.023 Mcps1212.276 MHz约 10 MBB2a10.23 Mcps661.38 MHz约 49 MB存储按complex single估算即 I/Q 各占 4 字节。如果仿真时长超过 200 毫秒应在生成过程中分段落盘避免内存峰值过高。上表里 B1C 的每码片采样数取 12是为了把 BOC(6,1) 的方波沿保留完整B1I 取 8 是滤波和存储的折中B2a 取 6 可以保证基带带宽覆盖主瓣。各参数可按实际接收机前端带宽调整但不建议低于表中的最小值。4. 用频谱、相关峰和相位连续性验证 B1I、B1C、B2a 仿真正确性4.1 用 pwelch 观察调制特征频谱是最快的正确性检查。B1C 生成后执行[pxx, f] pwelch(bb_b1c, hann(4096), 2048, 4096, 12.276e6); plot(f/1e6, 10*log10(pxx));B1C 频谱应出现两个对称主瓣两瓣间隔约 2.046 MHz中心频点有明显凹陷这是 BOC(1,1) 的指纹。如果看到的是单峰基本可以确定副载波没有乘进信号链。B1I 是单峰主瓣宽度约 4.09 MHzB2a 频谱更宽主瓣接近 20.46 MHz。三路信号放同一张图对比时先把峰值归一化到 0 dB再比较零点和凹点位置比直接看纵轴绝对值可靠。4.2 用 xcorr 检查码相位把发射基带信号与本地复制码做互相关可以验证码相位延迟与发射端设置一致[corrOut, lags] xcorr(bb_b1c, b1cCodeUp, coeff); [~, idx] max(abs(corrOut)); delayChip abs(lags(idx)) / spc;本地复制码必须用与发射端相同的spc上采样。如果峰位置偏移量等于你预设的延迟采样数说明码链路正确如果峰宽明显超标问题多出在滤波器阶数不足或采样率过低。测试时先关掉数据位翻转单独验证码对齐再加导航位观察相关峰衰减这样能快速区分码表问题和调制问题。4.3 数据位翻转带来的相关损失相关窗口跨过相邻数据比特边界时窗口内会有一段码片乘积为负。以 B1I 为例1 个数据比特占 20 ms若相干积分取 10 ms 且恰好跨边界一半码位的乘积反相相关幅度掉 3 dB。捕获实现里看到某颗星峰值显著低于其他星大多不是射频问题而是积分长度与导航位边界没有对齐。对策是发射端严格按符号边界生成数据接收端用非相干累加即对每个码相位的abs(corrOut).^2求和后再取峰。4.4 载波相位连续性与 BOC 极性跳变把发射基带与本地码共轭相乘观察相位随时间变化ph angle(bb_b1c .* conj(b1cCodeUp)); plot(ph);正常情况下相位包络只随载波剩余频率做线性变化。若在码周期边界出现固定幅度的阶跃且中间段连续多半是本地 BOC 副载波起始沿与发射端差了半个码片。这种极性跳变会引发码环抖动表现是伪距测量结果在小数部分来回跳。修法是回到bocSubcarrier里确认每个码片第一段是1还是-1不同版本标准对该参考点的定义不完全一致必须与仿真所依据的 ICD 版本保持一致。5. 把 B1I、B1C、B2a 仿真接到接收机链路前的实用参数处理5.1 统一导出为 int16 复数文件后处理常需要把仿真数据喂给 C 或 Python复数落盘格式要统一成 I/Q 交替tx tx / max(abs(tx)) * 32767; fid fopen(bds_three_band_iq.bin, w); iq [real(tx(:)), imag(tx(:))]; fwrite(fid, iq., int16); fclose(fid);先做峰值归一化再转int16避免削顶。iq.转置把矩阵变成 I0 Q0 I1 Q1 的连续排列这是多数软件接收机读取 bin 文件的标准布局。转成int16后数据量比double少四倍对长时间捕获验证很有价值。5.2 快速加入多普勒和多径场景接收机高动态测试不需要修改载波频率设置直接在基带上乘相位斜坡即可fd 5000; txDoppler tx .* exp(1j*2*pi*fd*t);fd从 500 Hz 到 5 kHz 范围内取值可以覆盖常见车载动态。多径场景用延迟副本叠加实现延迟码片数可取 0.5 到 2.0 个码片幅度设为主信号的 0.3 到 0.8 倍。需要注意幅度超过主信号时码环会错误锁定到延迟副本仿真结果会出现伪距台阶。5.3 量化位数与捕获灵敏度的折中接收机射频前端通常输出 2 bit 或 4 bit 量化数据仿真侧可以直接截断验证q round(tx * (2^(bitDepth-1)) / max(abs(tx)));量化位数每减少 1 bit捕获灵敏度大约损失 0.1 到 0.5 dB具体值受输入信号白噪声功率影响。同一份未经量化的仿真数据分别截成 2 bit 和 4 bit 后跑捕获函数对比相关峰值的下降量就能快速确定后续 FPGA 实现的量化位宽。把这段量化整数流直接接到现成的捕获函数里B1I、B1C、B2a 三路信号从码表到接收机前端的完整验证闭环就成立了。本文还有配套的精品资源点击获取