Skip to content

Commit

Permalink
Version 1.1
Browse files Browse the repository at this point in the history 
Fixes
- IF Inversion
- Fixed size() function warning in NAVDecoding.m (and similar)
- BDS B2a postNavigation.m satellite health check
- BDS B2a BCAVDecoding.m and satpost.m had a absense of a group delay differential correction factor
- Corrected Tracking.m progress bars
- Fixed GAL E5a ephemeris initializatrion in eph_struct_init.m
- README.md update
  • Loading branch information
@gnsscusdr
gnsscusdr committed Jul 1, 2024
1 parent 3204c93 commit 17eaa44
Show file tree
Hide file tree
Showing 98 changed files with 1,032 additions and 1,708 deletions.
Binary file added .DS_Store
View file
Binary file not shown.
Binary file added BDS/.DS_Store
View file
Binary file not shown.
Binary file added BDS/B1C/.DS_Store
View file
Binary file not shown.
2 changes: 1 addition & 1 deletion BDS/B1C/include/BCNAV1decoding.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -97,7 +97,7 @@
crcDet = comm.CRCDetector([24 23 18 17 14 11 10 7 6 5 4 3 1 0]);

%% B1C data decoding and ephemeris extract ==========================
for i = 1:size(index) % For each occurrence
for i = 1:height(index) % For each occurrence

% Take the B-CNAV1 symbols with one-frame length from data-channel
% prompt correlation values and change them to "1" and "0"
Expand Down
124 changes: 62 additions & 62 deletions BDS/B1C/include/NB_tracking.m
View file
Original file line number Diff line number Diff line change
Expand Up @@ -144,12 +144,12 @@

%% Start processing channels ==============================================
for channelNr = 1:settings.numberOfChannels

% Only process if PRN is non zero (acquisition was successful)
if (channel(channelNr).PRN ~= 0)
% Save additional information - each channel's tracked PRN
trackResults(channelNr).PRN = channel(channelNr).PRN;

% Move the starting point of processing. Can be used to start the
% signal processing at any point in the data record (e.g. for long
% records). In addition skip through that data file to start at the
Expand All @@ -158,53 +158,53 @@
fseek(fid, ...
dataAdaptCoeff*(settings.skipNumberOfBytes + channel(channelNr).codePhase-1), ...
'bof');

% Get a vector with the B1C data-channel code sampled 1x/chip
B1CData = generateDataBOC11(settings,channel(channelNr).PRN);
% Then make it possible to do early and late versions
B1CData = [B1CData(codeLength*2) B1CData B1CData(1)]; %#ok<AGROW>

% Get a vector with the pilot BOC(1,1) spreading waveform
pilotBOC11 = generatePilotBOC11(settings,channel(channelNr).PRN);
% Then make it possible to do early and late versions
pilotBOC11 = [pilotBOC11(codeLength*2) pilotBOC11 pilotBOC11(1)]; %#ok<AGROW>

%--- Perform various initializations ------------------------------

% define initial code frequency basis of NCO,
codeFreq = channel(channelNr).codeFreq;

% Define residual code phase (in chips)
remCodePhase = 0.0;

% Define carrier frequency which is used over whole tracking period
carrFreq = channel(channelNr).acquiredFreq;
carrFreqBasis = channel(channelNr).acquiredFreq;

% Define residual carrier phase
remCarrPhase = 0.0;

%code tracking loop parameters
oldCodeNco = 0.0;
oldCodeError = 0.0;

% Carrier/Costas loop parameters
d2CarrError = 0.0;
dCarrError = 0.0;

% For C/No computation
CNoValue = zeros(1,3);
tempCNoValue = zeros(1,3);

%=== Process the number of specified code periods =================
for loopCnt = 1:NumToProcess

%% GUI update -------------------------------------------------------------
% The GUI is updated every 200 ms. This way Matlab GUI is still
% responsive enough. At the same time Matlab is not occupied
% all the time with GUI task.
if (rem(loopCnt, 20) == 0)

Ln = newline;
trackingStatus=['Tracking: Ch ', int2str(channelNr), ...
' of ', int2str(TrackedNr),Ln ...
Expand All @@ -213,7 +213,7 @@
' of ', int2str(NumToProcess*10), ' msec',Ln...
'Data C/No: ',int2str(CNoValue(1)),' (dB-Hz);',...
' Pilot C/No: ',int2str(CNoValue(2)),' (dB-Hz)'];

try
waitbar(loopCnt/NumToProcess, ...
hwb, ...
Expand All @@ -225,39 +225,39 @@
return
end
end

%% Read next block of data ------------------------------------------------
% Record sample number (based on 8bit samples)
trackResults(channelNr).absoluteSample(loopCnt) =(ftell(fid))/dataAdaptCoeff;

% Update the phasestep based on code freq (variable) and
% sampling frequency (fixed)
codePhaseStep = codeFreq / settings.samplingFreq;

% Find the size of a "block" or code period in whole samples
blksize = ceil((codeLength-remCodePhase) / codePhaseStep);

% Read in the appropriate number of samples to process this
% interation
[rawSignal, samplesRead] = fread(fid, ...
dataAdaptCoeff*blksize, settings.dataType);

rawSignal = rawSignal';

if (dataAdaptCoeff==2)
rawSignal1=rawSignal(1:2:end);
rawSignal2=rawSignal(2:2:end);
rawSignal = rawSignal1 + 1i .* rawSignal2; % transpose vector
end

% If did not read in enough samples, then could be out of
% data - better exit
if (samplesRead ~= dataAdaptCoeff*blksize)
disp('Not able to read the specified number of samples for tracking, exiting!')
delete(hwb);
return
end

%% Set up all the code phase tracking information -------------------------
% Save remCodePhase for current correlation
trackResults(channelNr).remCodePhase(loopCnt) = remCodePhase;
Expand All @@ -267,147 +267,147 @@
((blksize-1) * codePhaseStep + remCodePhase- earlyLateSpc)*2;
tcode2 = ceil(tcode) + 1;
earlyCode = B1CData(tcode2);

% For pilot channel signal tracking
p11_earlyCode = pilotBOC11(tcode2);

% Define index into late code vector
tcode = (remCodePhase+earlyLateSpc)*2 : ...
codePhaseStep*2 : ...
((blksize-1)*codePhaseStep+remCodePhase+earlyLateSpc)*2;
tcode2 = ceil(tcode) + 1;
lateCode = B1CData(tcode2);

% For pilot channel signal tracking
p11_lateCode = pilotBOC11(tcode2);

% Define index into prompt code vector
tcode = remCodePhase*2 : ...
codePhaseStep*2 : ...
((blksize-1)*codePhaseStep+remCodePhase)*2;
tcode2 = ceil(tcode) + 1;
promptCode = B1CData(tcode2);

% For pilot channel signal tracking
p11_promptCode = pilotBOC11(tcode2);

remCodePhase = tcode(blksize)/2 + codePhaseStep - codeLength;

%% Generate the carrier frequency to mix the signal to baseband -----------

% Save remCarrPhase for current correlation
trackResults(channelNr).remCarrPhase(loopCnt) = remCarrPhase;

% Get the argument to sin/cos functions
time = (0:blksize) ./ settings.samplingFreq;
trigarg = ((carrFreq * 2.0 * pi) .* time) + remCarrPhase;
remCarrPhase = rem(trigarg(blksize+1), (2 * pi));

% Finally compute the signal to mix the collected data to
% bandband
carrsig = exp(1i .* trigarg(1:blksize));
carrsig = exp(-1i .* trigarg(1:blksize));

%% Generate the six standard accumulated values ---------------------------
% First mix to baseband
qBasebandSignal = real(carrsig .* rawSignal);
iBasebandSignal = imag(carrsig .* rawSignal);
iBasebandSignal = real(carrsig .* rawSignal);
qBasebandSignal = imag(carrsig .* rawSignal);

% Now get early, late, and prompt values for each
I_E = sum(earlyCode .* iBasebandSignal);
Q_E = sum(earlyCode .* qBasebandSignal);
I_P = sum(promptCode .* iBasebandSignal);
Q_P = sum(promptCode .* qBasebandSignal);
I_L = sum(lateCode .* iBasebandSignal);
Q_L = sum(lateCode .* qBasebandSignal);

% Correlation values for pilot BOC(1,1) spreading waveform
p11_I_E = sum(p11_earlyCode .* iBasebandSignal);
p11_Q_E = sum(p11_earlyCode .* qBasebandSignal);
p11_I_P = sum(p11_promptCode .* iBasebandSignal);
p11_Q_P = sum(p11_promptCode .* qBasebandSignal);
p11_I_L = sum(p11_lateCode .* iBasebandSignal);
p11_Q_L = sum(p11_lateCode .* qBasebandSignal);

%% Find PLL error and update carrier NCO ----------------------------------

% Implement carrier loop discriminator (phase detector)
carrError = atan(Q_P / I_P) / (2.0 * pi);

% Combined code tracking error estimation using data and pilot
% chaannel signals
% Pilot channel carrier phase is pi/2 rad ahead of the
% data channel carrier phase.atan is not affectede by
% the secondary code modulation
p11_carrError = atan(-p11_I_P/p11_Q_P)/ (2.0 * pi);

% Composite carrier tracking error
carrError = (carrError *11 + p11_carrError*29)/40;

% Implement carrier loop filter and generate NCO command
d2CarrError = d2CarrError + carrError * pf3;
dCarrError = d2CarrError + carrError * pf2 + dCarrError;
carrNco = dCarrError + carrError * pf1;

% Save carrier frequency for current correlation
trackResults(channelNr).carrFreq(loopCnt) = carrFreq;
% Modify carrier freq based on NCO command
carrFreq = carrFreqBasis + carrNco;

%% Find DLL error and update code NCO -------------------------------------
codeError = (sqrt(I_E ^2 + Q_E ^2) - sqrt(I_L ^2 + Q_L ^2)) / ...
(sqrt(I_E ^2 + Q_E ^2) + sqrt(I_L ^2 + Q_L ^2)) * (1-earlyLateSpc);

% Combined code tracking error estimation using data and pilot
% chaannel signals
p11_codeError = (sqrt(p11_I_E ^2 + p11_Q_E ^2) - sqrt(p11_I_L ^2 + p11_Q_L ^2)) / ...
(sqrt(p11_I_E ^2 + p11_Q_E ^2) + sqrt(p11_I_L ^2 + p11_Q_L ^2))* (1-earlyLateSpc);
% Composite code tracking error
codeError = (codeError*11 + p11_codeError*29)/40;

% Implement code loop filter and generate NCO command
codeNco = oldCodeNco + (tau2code/tau1code) * ...
(codeError - oldCodeError) + codeError * (PDIcode/tau1code);
oldCodeNco = codeNco;
oldCodeError = codeError;

% Save code frequency for current correlation
trackResults(channelNr).codeFreq(loopCnt) = codeFreq;
% Modify code freq based on NCO command
codeFreq = channel(channelNr).codeFreq - codeNco;

%% Record various measures to show in postprocessing ----------------------


trackResults(channelNr).dllDiscr(loopCnt) = codeError;
trackResults(channelNr).dllDiscrFilt(loopCnt) = codeNco;
trackResults(channelNr).pllDiscr(loopCnt) = carrError;
trackResults(channelNr).pllDiscrFilt(loopCnt) = carrNco;

trackResults(channelNr).I_E(loopCnt) = I_E;
trackResults(channelNr).I_P(loopCnt) = I_P;
trackResults(channelNr).I_L(loopCnt) = I_L;
trackResults(channelNr).Q_E(loopCnt) = Q_E;
trackResults(channelNr).Q_P(loopCnt) = Q_P;
trackResults(channelNr).Q_L(loopCnt) = Q_L;

trackResults(channelNr).Pilot_I_P(loopCnt) = p11_I_P ;
trackResults(channelNr).Pilot_Q_P(loopCnt) = p11_Q_P;
%% CNo calculation --------------------------------------------------------

if (rem(loopCnt,settings.CNoInterval)==0)
% Computation of CNo and PLL detector output
[CNoValue, PllDetector]= ...
Calc_CNo_PLD(trackResults(channelNr),settings,loopCnt);

CNoCnt = loopCnt/settings.CNoInterval;

% Save C/No for data channel: a o.5-0.5 filter is used to
% smooth the results
trackResults(channelNr).DataCNo(CNoCnt) = ...
CNoValue(1) * 0.5 + tempCNoValue(1) * 0.5;
% Save PLL lock detector output for data channel
trackResults(channelNr).DataPLD(CNoCnt) = PllDetector(1);

% Save C/No and PLL lock detector output for pilot channel
trackResults(channelNr).PilotCNo(CNoCnt) = ...
CNoValue(2) * 0.5 + tempCNoValue(2) * 0.5;
Expand All @@ -416,14 +416,14 @@
trackResults(channelNr).PilotPLD(CNoCnt) = PllDetector(2);
end
tempCNoValue = CNoValue;

end % for loopCnt

% If we got so far, this means that the tracking was successful
% Now we only copy status, but it can be update by a lock detector
% if implemented
trackResults(channelNr).status = channel(channelNr).status;

end % if a PRN is assigned
end % for channelNr

Expand Down
Loading

0 comments on commit 17eaa44

Please sign in to comment.