plot_sim.m

来自「Simple GPS Simulation with RINEX Data」· M 代码 · 共 145 行

function plot_sim(offset,xlims,PRN)

% offset = str2num(offset);
% PRN = str2num(PRN);
% xlims = str2num(xlims);
global s0;

[s0 t] = create_code_samples_for_demo(PRN,0);
[s1 t] = create_code_samples_for_demo(PRN,offset);
S2 = s0.*s1;
S2_mean = mean(S2);
s0f=fft(s0);
R=ifft(s0f.*conj(s0f));
R=fftshift(R)/length(R);

H = figure(1);
set(H,'name','Code Shift Demo');
set(H,'units','normalized');
set(H,'position',[.35 .1 .5 .8]);
subplot(411);
title1 = sprintf('PRN %1.0f Code (40 out of 1023 chips)',PRN);
plot(t,s0);title(title1);
set(gca,'xlim',[-20 20]);
set(gca,'ylim',[-1.1 1.1]);
subplot(412);
title2 = sprintf('Offset PRN %1.0f Code (Offset = %3.1f chips)',PRN,offset);
plot(t,s1);title(title2);
set(gca,'xlim',[-20 20]);
set(gca,'ylim',[-1.1 1.1]);
subplot(413);
plot(t,S2);title('Correlator Output');
set(gca,'xlim',[-20 20]);
set(gca,'ylim',[-1.1 1.1]);
subplot(414);
plot(t,R,offset,S2_mean,'o');title('Auto-Correlation Curve');
xlabel('Offset (chips)');set(gca,'xlim',xlims);
set(gca,'ylim',[-0.2 1.1]);

%------------------------------------------
%------------------------------------------

function [code_samples, time_vec] = create_code_samples_for_demo(prn, offset)
% function code_samples = create_code_samples_for_demo(prn, offset)

global ca_code_sequence;
fc=1.023E6;  % Assumes constant sample rate!
Tc=1/fc;

sample_interval=1/100 * Tc;
ca_code_sequence=ca_code(prn);
time_vec=[-.0005:sample_interval:0.0005-sample_interval];
start_time = offset*Tc;

% Next, generate code sequence
start_time=start_time-eps;
dt=time_vec-start_time;
code_index=mod(ceil(dt*fc),length(ca_code_sequence));

% Correct for 1023 mod 1023 which returns 0 but should be 1023
code_index(find(code_index==0))=length(ca_code_sequence);
code_samples=ca_code_sequence(code_index);

time_vec = time_vec * fc;

%------------------------------------------
%------------------------------------------

function code = ca_code(prn)
% function code = ca_code(prn)
% This function returns the C/A code sequence (1023 bytes long,
% represented as 1023x1 vector of integers with values of +/- 1

% Error messages

if nargin~=1
   error('insufficient number of input argumnets')
end
if prn<0 | prn > 37
   error('invalid prn: must be between 1 and 37')
end

g2shift_vector=[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 ...
   863 950 947 948 950];

g2shift=g2shift_vector(prn);


%
% Generate G1 code
%

% load shift regeister

for i=1:10,
   reg(i) = -1;
end

% Generate code

for i=1:1023,
   g1(i) = reg(10);
   save1 = reg(3)*reg(10);
   for j=9:-1:1,
      reg(j+1)=reg(j);
   end
   reg(1)=save1;
end

%
% Generate G2 code
%

% load shift regeister

for i=1:10,
   reg(i) = -1;
end

% Generate code

for i=1:1023,
   g2(i) = reg(10);
   save2 = reg(2)*reg(3)*reg(6)*reg(8)*reg(9)*reg(10);
   for j=9:-1:1,
      reg(j+1)=reg(j);
   end
   reg(1)=save2;
end

% Shift the G2 code

for i=1:1023,
   k = i+g2shift;
   if k > 1023,
     k = k-1023;
   end
   g2tmp(k) = g2(i);
end

g2 = g2tmp;

% Form the C/A code by multiplying G1 and G2

code = g1.*g2;