www.pudn.com > GNSS_SDR.zip > tracking.asv, change:2010-02-05,size:16002b


function [trackResults, channel]= tracking(fid, channel, settings) 
% Performs code and carrier tracking for all channels. 
% 跟踪各个通道码和载波 
%[trackResults, channel] = tracking(fid, channel, settings) 
% 
%   Inputs: 
%       fid             - file identifier of the signal record. 
%                       - 文件指针 
%       channel         - PRN, carrier frequencies and code phases of all 
%                       satellites to be tracked (prepared by preRum.m from 
%                       acquisition results). 
%                       - 根据捕获结果PRN,载波频率,码相位建立的通道 
%       settings        - receiver settings. 
%                       - 接收机设置 
%   Outputs: 
%       trackResults    - tracking results (structure array). Contains 
%                       in-phase prompt outputs and absolute starting  
%                       positions of spreading codes, together with other 
%                       observation data from the tracking loops. All are 
%                       saved every millisecond. 
%                       - 跟踪结果(结构体)。包括I通道输出,码偏移以及其他环路 
%                       信息,输出周期1ms 
 
%-------------------------------------------------------------------------- 
%                           SoftGNSS v3.0 
%  
% Copyright (C) Dennis M. Akos 
% Written by Darius Plausinaitis and Dennis M. Akos 
% Based on code by DMAkos Oct-1999 
% 注释翻译:苗剑峰 
%-------------------------------------------------------------------------- 
%This program is free software; you can redistribute it and/or 
%modify it under the terms of the GNU General Public License 
%as published by the Free Software Foundation; either version 2 
%of the License, or (at your option) any later version. 
% 
%This program is distributed in the hope that it will be useful, 
%but WITHOUT ANY WARRANTY; without even the implied warranty of 
%MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 
%GNU General Public License for more details. 
% 
%You should have received a copy of the GNU General Public License 
%along with this program; if not, write to the Free Software 
%Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301, 
%USA. 
%程序为开源程序,请在开源协议规定范围内做改动 
%-------------------------------------------------------------------------- 
 
%CVS record: 
%$Id: tracking.m,v 1.14.2.32 2007/01/30 09:45:12 dpl Exp $ 
 
%% Initialize result structure 初始化跟踪结果结构体========================= 
 
% Channel status 
% 通道状态 
trackResults.status         = '-'; % No tracked signal, or lost lock 无信号或者失锁 
 
% The absolute sample in the record of the C/A code start: 
% C/A码起始点  
trackResults.absoluteSample = zeros(1, settings.msToProcess); 
 
% Freq of the C/A code: 
% C/A码频率 
trackResults.codeFreq       = inf(1, settings.msToProcess); 
 
% Frequency of the tracked carrier wave: 
% 载波频率 
trackResults.carrFreq       = inf(1, settings.msToProcess); 
 
% Outputs from the correlators (In-phase): 
% I通道相关输出 
trackResults.I_P            = zeros(1, settings.msToProcess); 
trackResults.I_E            = zeros(1, settings.msToProcess); 
trackResults.I_L            = zeros(1, settings.msToProcess); 
 
% Outputs from the correlators (Quadrature-phase): 
% Q通道相关输出 
trackResults.Q_E            = zeros(1, settings.msToProcess); 
trackResults.Q_P            = zeros(1, settings.msToProcess); 
trackResults.Q_L            = zeros(1, settings.msToProcess); 
 
trackResults.SN		        = zeros(1, floor(settings.msToProcess/1000)); 
 
% Loop discriminators 
% 检相器输出 
trackResults.dllDiscr       = inf(1, settings.msToProcess); 
trackResults.dllDiscrFilt   = inf(1, settings.msToProcess); 
trackResults.pllDiscr       = inf(1, settings.msToProcess); 
trackResults.pllDiscrFilt   = inf(1, settings.msToProcess); 
 
%--- Copy initial settings for all channels ------------------------------- 
%--- 复制到各个通道 -------------------------------------------------------- 
trackResults = repmat(trackResults, 1, settings.numberOfChannels); 
 
%% Initialize tracking variables 初始化跟踪变量 ============================= 
 
codePeriods = settings.msToProcess; % For GPS one C/A code is one ms 跟踪时间 
 
%--- DLL variables 码环变量---------------------------------------------------- 
% Define early-late offset (in chips) 
% 定义超前滞后偏移量 
earlyLateSpc = settings.dllCorrelatorSpacing; 
 
% Summation interval 
% 积分累加时间 
PDIcode = 0.001; 
 
% Calculate filter coefficient values 
% 计算滤波器参数 
[tau1code, tau2code] = calcLoopCoef(settings.dllNoiseBandwidth, ... 
                                    settings.dllDampingRatio, ... 
                                    1.0); 
 
%--- PLL variables 载波环变量--------------------------------------------------- 
% Summation interval 
% 积分累加时间 
PDIcarr = 0.001; 
 
% Calculate filter coefficient values 
% 计算滤波器参数 
[tau1carr, tau2carr] = calcLoopCoef(settings.pllNoiseBandwidth, ... 
                                    settings.pllDampingRatio, ... 
                                    0.25); 
hwb = waitbar(0,'跟踪进行中...'); 
 
%% Start processing channels 开始跟踪通道计算 ================================== 
for channelNr = 1:settings.numberOfChannels 
     
    % Only process if PRN is non zero (acquisition was successful) 
    % 只跟踪捕获成功的卫星 
    if (channel(channelNr).PRN ~= 0)  % channel:preRun 的输出参数<include> 
        % Save additional information - each channel's tracked PRN 
        % 保存卫星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 
        % appropriate sample (corresponding to code phase). Assumes sample 
        % type is schar (or 1 byte per sample)  
        % 跳过开头,并定位到输入信号C/A码起始点 
        fseek(fid, ... 
              settings.skipNumberOfBytes + channel(channelNr).codePhase-1, ... 
              'bof'); 
 
 
        % Get a vector with the C/A code sampled 1x/chip 
        % 构建本地C/A码向量 
        caCode = generateCAcode(channel(channelNr).PRN); 
        % Then make it possible to do early and late versions 
        % 前后各扩展一个码片,保证可以计算偏移 
        caCode = [caCode(1023) caCode caCode(1)]; 
 
        %--- Perform various initializations 初始化 ------------------------------ 
 
        % define initial code frequency basis of NCO 
        % 初始码频率 
        codeFreq      = settings.codeFreqBasis; 
        % 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 
        % 载波环路参数 
        oldCarrNco   = 0.0; 
        oldCarrError = 0.0; 
		 
		% CN0 estimate parameter 
		% CN0估计参数 
		m = 0; 
		R = 0; 
		 
		 
        %=== Process the number of specified code periods 开始跟踪 ========= 
        for loopCnt =  1:codePeriods 
             
%% GUI update 界面更新------------------------------------------------------ 
            % The GUI is updated every 50ms. This way Matlab GUI is still 
            % responsive enough. At the same time Matlab is not occupied 
            % all the time with GUI task. 
            % 50ms更新一次界面,保证Matlab响应速度 
            if (rem(loopCnt, 50) == 0) 
                try 
                    waitbar(loopCnt/codePeriods, ... 
                            hwb, ... 
                            ['跟踪通道 ', int2str(channelNr), ... 
                            ' / ', int2str(settings.numberOfChannels), ... 
                            '; 卫星PRN:', int2str(channel(channelNr).PRN), ... 
                            '; 已完成 ',int2str(loopCnt), ... 
                            ' / ', int2str(codePeriods), ' ms']);                        
                catch 
                    % The progress bar was closed. It is used as a signal 
                    % to stop, "cancel" processing. Exit. 
                    % 跟踪进度条取消时退出跟踪 
                    disp('跟踪被取消...'); 
                    return 
                end 
			end 
			 
	 
%% Read next block of data 读入一段数据 ------------------------------------          
            % Find the size of a "block" or code period in whole samples 
            % Update the phasestep based on code freq (variable) and 
            % sampling frequency (fixed) 
            % 根据码频率(变量)和采样率(固定)计算数据段长度 
            codePhaseStep = codeFreq / settings.samplingFreq; % 每次采样的码片数 
             
            blksize = ceil((settings.codeLength-remCodePhase) / codePhaseStep);% 一周期CA所需的采样数 
             
            % Read in the appropriate number of samples to process this 
            % interation  
            % 读入数据 
            [rawSignal, samplesRead] = fread(fid, ... 
                                             blksize, settings.dataType); 
            rawSignal = rawSignal';  %transpose vector 转置 
             
            % If did not read in enough samples, then could be out of  
            % data - better exit  
            % 无法读入足够数据,退出 
            if (samplesRead ~= blksize) 
                disp('无法读入足够数据,退出!') 
                fclose(fid); 
                return 
            end 
 
%% Set up all the code phase tracking information 建立所有码相位跟踪信息 ---- 
            % Define index into early code vector 
            % 超前码 
            tcode       = (remCodePhase-earlyLateSpc) : ... 
                          codePhaseStep : ... 
                          ((blksize-1)*codePhaseStep+remCodePhase-earlyLateSpc); 
            tcode2      = ceil(tcode) + 1; 
            earlyCode   = caCode(tcode2); 
             
            % Define index into late code vector 
            % 滞后码 
            tcode       = (remCodePhase+earlyLateSpc) : ... 
                          codePhaseStep : ... 
                          ((blksize-1)*codePhaseStep+remCodePhase+earlyLateSpc); 
            tcode2      = ceil(tcode) + 1; 
            lateCode    = caCode(tcode2); 
             
            % Define index into prompt code vector 
            % 当前码 
            tcode       = remCodePhase : ... 
                          codePhaseStep : ... 
                          ((blksize-1)*codePhaseStep+remCodePhase); 
            tcode2      = ceil(tcode) + 1; 
            promptCode  = caCode(tcode2); 
%             codesumPhase=codesumPhase+tcode(blksize)+codePhaseStep; 
             
            remCodePhase = (tcode(blksize) + codePhaseStep) - 1023.0; 
%% Generate the carrier frequency to mix the signal to baseband 产生载波---- 
            time    = (0:blksize) ./ settings.samplingFreq; 
             
            % Get the argument to sin/cos functions 
            
            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 
            % 产生同相、正交信号 
            carrCos = cos(trigarg(1:blksize)); 
            carrSin = sin(trigarg(1:blksize)); 
 
%% Generate the six standard accumulated values 计算6组累加量--------------- 
            % First mix to baseband 
            % 剥离载波 
            qBasebandSignal = carrCos .* rawSignal; 
            iBasebandSignal = carrSin .* rawSignal; 
 
            % Now get early, late, and prompt values for each 
            % 剥离CA码,得到6路相关结果 
            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); 
			 
%% Calculator CN0 计算载噪比----------------------------------------------- 
 
            m = m + I_P^2 + Q_P^2 ; 
            R = R + (I_P^2 + Q_P^2)^2; 
 
%% Find PLL error and update carrier NCO 更新载波环------------------------- 
 
            % Implement carrier loop discriminator (phase detector) 
            % 载波检相器 
            carrError = atan(Q_P / I_P) / (2.0 * pi); 
 
            % Implement carrier loop filter and generate NCO command 
            % 数字二阶低通滤波器 
            carrNco = oldCarrNco + (tau2carr/tau1carr) * ... 
                (carrError - oldCarrError) + carrError * (PDIcarr/tau1carr); 
            oldCarrNco   = carrNco; 
            oldCarrError = carrError; 
 
            % Modify carrier freq based on NCO command 
            % 更新载波频率 
            carrFreq = carrFreqBasis + carrNco; 
 
 
            trackResults(channelNr).carrFreq(loopCnt) = carrFreq; 
 
%% Find DLL error and update code NCO 更新码环------------------------------ 
 
            % 码环检相器 
            codeError = (sqrt(I_E * I_E + Q_E * Q_E) - sqrt(I_L * I_L + Q_L * Q_L)) / ... 
                (sqrt(I_E * I_E + Q_E * Q_E) + sqrt(I_L * I_L + Q_L * Q_L)); 
             
            % Implement code loop filter and generate NCO command 
            % 数字二阶低通滤波器 
            codeNco = oldCodeNco + (tau2code/tau1code) * ... 
                (codeError - oldCodeError) + codeError * (PDIcode/tau1code); 
            oldCodeNco   = codeNco; 
            oldCodeError = codeError; 
             
            % Modify code freq based on NCO command 
            % 更新码频率 
            codeFreq = settings.codeFreqBasis - codeNco; 
             
            trackResults(channelNr).codeFreq(loopCnt) = codeFreq; 
 
%% Record various measures to show in postprocessing 记录数据--------------- 
            % Record sample number (based on 8bit samples) 
            % 记录采样数据文件指针 
            trackResults(channelNr).absoluteSample(loopCnt) = ftell(fid); 
 
            trackResults(channelNr).dllDiscr(loopCnt)       = codeError*0.5; 
            trackResults(channelNr).dllDiscrFilt(loopCnt)   = codeNco; 
            trackResults(channelNr).pllDiscr(loopCnt)       = carrError; 
            trackResults(channelNr).pllDiscrFilt(loopCnt)   = carrNco; 
            trackResults(channelNr).coderem(loopCnt)        = remCodePhase; 
 
            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; 
			 
			rem_code(loopCnt)=remCodePhase; 
			rem_carr(loopCnt)=remCarrPhase; 
%% Calculator CN0 计算载噪比 ----------------------------------------------- 
 
			if mod(loopCnt,1000) == 0 
				 
				m = m / 1000; 
				R = R / 1000; 
 
				SN = 10 * log10( ( 4*m^2 ) / (R-m^2) -3 ) + 24 ; 
				trackResults(channelNr).SN(loopCnt/1000) = SN; 
 
				m = 0; 
				R = 0;			 
 
			end			 
			 
        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  
 
 
% Close the waitbar 
% 关闭进度条 
close(hwb)