processrinex.m

利用gps多天线载波相位技术

%%========================================
%%     Toolbox for attitude determination
%%     Zhen Dai
%%     [e-mail]
%%     ZESS, University of Siegen, Germany
%%     Last Modified  : 1.Sep.2008
%%========================================
%% Functions:
%%      1. RINEX observation file reading 
%%      2. RINEX navigation file reading
%%      3. Synchronize the observation data from all antennas
%% Input parameters:
%%      ------------------------------------------------------------------------
%%                          Variables from the GUI 
%%          1. vFileObsNames  
%%                  the RINEX observation file names
%%          2. total_tnt                     
%%                  total number of valid baseline input
%%          3. filenav
%%                  name of the RINEX navigation file
%%          4. smoothing_interval
%%                  smoothing interval (could be 0)
%%          5. total_epoch
%%                  how many epochs to be processed (-999=all data)
%%          6. maskangle
%%                  elevation mask angle                  
%%          7. mBaseline -> magnitude of baselines 
%%                  magnitude of baselines 
%%      -----------------------------------------------------------------------
%% Output: (To the data file DataSatRecord)
%%      mMeasurement->Satellite IDs, Code, Phase of all antennas
%%      mCommonEpoch-> common epochs of the data from all antennas
%%      mEpochStr->common epochs in the form of STRING
%%      mBaseline -> magnitude of baselines  
%%      smoothing_interval ->  smoothing interval (could be 0)
%%      maskangle -> elevation mask angle  

function ProcessRINEX()

totalGPSsatellite=32;
load DataGUI
%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%%          Read each RINEX observation file
%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
for i=1:1:total_ant,
    filename=vFileObsNames{i}; 
    datafilename=LoadRinexOBS(filename,i,total_epoch);
end

%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%%                  Analyse the GPS observations
%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%% Check the common epochs and minimal common interval
clear vFitstEpoch vLastEpoch vInterval
for i=1:1:total_ant,
    clear mSatID vTime interval mCode mPhase
    filename=sprintf('DataRecord%1d',i);
    load (filename);
    totalepoch=length(vTime);
    vFitstEpoch(i)=vTime(1);  
    vLastEpoch(i)=vTime(totalepoch);
    vInterval(i)=interval;
end

%% The first and last common epochs
first_epoch_common=max(vFitstEpoch);
last_epoch_common=min(vLastEpoch);
%% Find the least common multiple of the sampling rates
%% multiplying 100 is to solve the high data rate like 20Hz
interval_common=round(vInterval(1)*100);
for i=2:1:total_ant;
    interval_common=lcm(interval_common,round(vInterval(i)*100));
end
interval_common=interval_common/100;
if abs(interval_common-round(interval_common))<1e-3,
    interval_common=round(interval_common);
end
%% No data overlap
if (first_epoch_common>=last_epoch_common),
    error('No overlapped observation');
end
%% Too large common sampling rate
if interval_common>(last_epoch_common-first_epoch_common)
    error('Difficult to implement synchronization');
end

%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%% Re-arrange the observation data at synchronized epochs
%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%% We only concern the synchronized data 
total_common_epoch=...
    round(last_epoch_common-first_epoch_common)/interval_common;
%% Initialize the matrix recording the observation
clear mMeasurement mCommonEpoch mEpochStr
mMeasurement=zeros(total_ant,total_common_epoch,5,totalGPSsatellite);
mCommonEpoch=zeros(total_common_epoch,1);
%% Now assign the matrix mMeasurement
%% including the satellite IDs, pseudorange, phase data
%% of the synchronized epochs
for antid=1:1:total_ant,       
    clear mSatID vTime interval mCode mPhase
    filename=sprintf('DataRecord%1d',antid);
    load (filename);
    total_epoch=length(vTime);
    num_epoch=1;
    %% pick up the common epochs from the original measurement
    for j=1:1:total_epoch,
        if abs(vTime(j)-(first_epoch_common+(num_epoch-1)*interval_common))<interval_common,
            vSatID=find(mCode(:,j)~=0);
            %% satellite IDs
            mMeasurement(antid,num_epoch,1,vSatID)=1;
            %% pseudorange
            mMeasurement(antid,num_epoch,2,vSatID)=mCode(vSatID,j);
            %% phase data
            mMeasurement(antid,num_epoch,3,vSatID)=mPhase(vSatID,j); 
            %% common (synchronized) epochs
            if antid==1,
                mCommonEpoch(num_epoch)=vTime(j);
                mEpochStr{num_epoch}=vTimeStr{j};
            end 
            num_epoch=num_epoch+1;
            if num_epoch>total_common_epoch, break; end
        end
    end %% end for epochs 
end %% loop for antenna 

%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%%               Analyse the RINEX navigation file
%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 
%% Result will be saved into a data file "RinexNav"
LoadRinexNav(filenav); 

%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
%%                Save the data into a data file
%% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
text{1}=sprintf('Analysis of Rinex files is finished.');
text{2}=sprintf('-----------------------------------------------------------------------');
text{3}=sprintf('Common epoch=%d',length(mCommonEpoch));
text{4}=sprintf('Number of antennas(valid baselines)=%d',total_ant);
text{5}=sprintf('Smoothing interval=%d (epochs)',smoothing_interval);
text{6}=sprintf('Elevation mask angle=%d (degree)',maskangle);
text{7}=sprintf('Data rate=%3.2f (s)',interval_common);
text{8}=sprintf('-----------------------------------------------------------------------');
text{9}=sprintf('');
text{10}=sprintf('Data processing starts in 5 seconds, please check the progress.' );

%% Write the output file "Results.txt"
file_output=fopen('Results.txt','wt');
vClock=clock;
fprintf(file_output,'%s %d:%d:%2.0f\n',date,vClock(4),vClock(5),vClock(6));
for i=1:1:total_ant,
    filename=vFileObsNames{i}; 
    fprintf(file_output,'Observation data of #%d antenna : %s\n',i,filename);
end
fprintf(file_output,'Ephemerides data : %s\n', filenav);
fprintf(file_output,'%s\n', text{3});
fprintf(file_output,'%s\n', text{4});
fprintf(file_output,'%s\n', text{5});
fprintf(file_output,'%s\n', text{6});
fprintf(file_output,'\n');
fclose(file_output);

msg1=msgbox(text); 

%% Save the data to a data file for further use
save('DataSatRecord','mMeasurement','mCommonEpoch','mEpochStr','maskangle','smoothing_interval','mBaseline','num_baseline','-v6')