GPS信号生成：C/A码序列生成【MATLAB实现】

GPS C/A码序列生成【MATLAB实现】

在本文中，将简要介绍GPS C/A码及其生成原理，并且用MATLAB代码实现。

GPS信号与C/A码

GPS的信号主要有三类：载波（carrier）、测距码（也可以说是伪随机噪声码，Pseudo Random Noise，PRN）和导航电文（Navigation message）。GPS卫星在实现定位时就使用以上三种信号进行相应的伪距（Pseudorange）测量或载波相位（Carrier phase）测量。
本文介绍实现伪距测量的一种信号：C/A码（Coarse/Acquisition Code）是用于进行粗略测距和捕获P码的粗码，也称捕获码。周期T为1毫秒，一个周期含有码元即码长N为2^10 - 1=1023，每个码元持续的时间即码元周期t=1 ms/1023=0.977517微秒，相应的码元宽度为293.05米。C/A码是一种公开的明码，可供全球用户免费使用。但C/A码在L1载波上，所以无法精确地消除电离层延迟，测距精度一般为±（2~3）米。C/A码和P码（Precise code）同属PRN码。但是C/A码为民用，P码通常为军用，所以我们在学习和实验时基本上使用的是C/A码。
C/A码有1023个元素（chips），这里面含有512个1与511个0，他们的排列是看起来仿佛随机的，但是却是完全可决定的，因此叫做伪随机噪声序列。

GPS C/A码序列生成原理

GPS C/A码的生成通常需要经过几个步骤：如下图，窄带宽基带信号（baseband signal）经 PRN 码调制后，基带信号的频谱明显变宽。也就是说，信号从narrow bandwidth变成了wide bandwidth。基带信号经过扩频处理后，使用二进制相移键控（BPSK）调制方法对扩频信号进行调制。最后，BPSK 扩频信号与 在L1频段1575.42 MHz 的载波信号一起传输。

GPS的L1频段的码率是50 bps（bit per second）。在这50 bps中，信号的频谱宽度是-50 Hz到50 Hz。另外，C/A码的码率是1.023 MHz，也就是1 ms内可以得到1023 chips的C/A码。
C/A码的生成遵循Gold code原理，这里我放了1967年的论文原文。在C/A码生成时，需要G1和G2两个线性反馈移位寄存器（Linear feedback shift register，LFSR）。每个寄存器都有10位，所以按照Gold code原理生成的序列长度是2^10-1=1023 chips。寄存器G1和G2按照异或（Exclusive or，XOR）的方式生成C/A码序列，如下图所示。一般来说，移位寄存器是一组一位存储单元（存储单元）。在时钟脉冲作用下，每个寄存器的内容将向右移动一位。输出是最后一个单元（最左边的单元或寄存器）的内容。

需要注意的是，G1和G2分别有自己的反馈准则：

对于 G1，在每个时钟脉冲下，寄存器位置 3 和 10 的进行异或运算，并将结果移至寄存器位置 1。其他寄存器向右移动一个位置，最后寄存器位置 10 的内容成为 G1 的输出G1。对于 G2，每个时钟脉冲都会对寄存器位置 2、3、6、7、8、9、10 的进行异或运算，并将结果移至寄存器位置 1。同样，所有其他寄存器的内容也会向右移动一位。然后，选择两个selector，进行异或运算，再offset移位运算，成为G2i。最后，将两个寄存器的输出G1 和 G2i 进行异或，得到每个epoch脉冲上产生的 PRN 码位。
每个PRN对应的selector和offset的选择参考GPS Interface Specification IS-GPS-200 的表，如下图：

像这样重复运算1023次之后，我们就得到1 ms的C/A码序列。

MATLAB代码实现

附上生成C/A码序列的MATLAB代码，输入PRN号，得到相应的C/A码。

function code = CA_generation(prn)
% Generate the C/A code for a given PRN (1-32)
% Input:
%   prn - PRN number (1-32)
% Output:
%   code - 1023-length C/A code sequence% Initial states for Linear feedback shift register (LFSR) of length G1 and G2
G1 = ones(1, 10); 
G2 = ones(1, 10);% PRN taps (satellite-specific tap positions)
prn_taps_table = [2, 6;   % PRN13, 7;   % PRN24, 8;   % PRN35, 9;   % PRN41, 9;   % PRN52, 10;  % PRN61, 8;   % PRN72, 9;   % PRN83, 10;  % PRN92, 3;   % PRN103, 4;   % PRN115, 6;   % PRN126, 7;   % PRN137, 8;   % PRN148, 9;   % PRN159, 10;  % PRN161, 4;   % PRN172, 5;   % PRN183, 6;   % PRN194, 7;   % PRN205, 8;   % PRN216, 9;   % PRN221, 3;   % PRN234, 6;   % PRN245, 7;   % PRN256, 8;   % PRN267, 9;   % PRN278, 10;  % PRN281, 6;   % PRN292, 7;   % PRN303, 8;   % PRN314, 9;   % PRN32
];offset_table=[5;6;7;8;17;18;139;140;141;251;252;254;255;256;257;258;469;470;471;472;473;474;509;512;513;514;515;516;859;860;861;862;
];% Validate PRN input
if prn < 1 || prn > 32error('PRN must be between 1 and 32.');
end% Get the taps for the given PRN
selector = prn_taps_table(prn, :);
select1 = selector(1);
select2 = selector(2);% Get the offset
offset=offset_table(prn,:);% Generate the C/A code sequence
code = zeros(1, 1023);for i = 1:1023% Generate the code bit (modulo-2 addition of selected G2 taps)G1_output = G1(10);G2_output = xor(G2(select1), G2(select2)); % Adjust G2 output by offset% Shift G1feedback1 = xor(G1(3), G1(10));G1 = [feedback1, G1(1:end-1)];% Shift G2feedback2 = xor(G2(2), xor(G2(3), xor(G2(6), xor(G2(8), xor(G2(9), G2(10))))));G2 = [feedback2, G2(1:end-1)];% Apply offset to G2G2_shifted = circshift(G2_output, [0, offset]);code(i) = xor(G1_output, G2_shifted); % Generate C/A code
end% Convert 0 to -1 for numerical representation
code = 1 - 2 * code; % Change logic value 0 to +1 and 1 to -1
end

Github小项目在此，实现GPS信号的生成和获取，喜欢的朋友可以给个小星星。