rinex_nav_to_eph.m

Simple GPS Simulation with RINEX Data

function rinex_nav_to_eph(rinex_nav_file, eph_file)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% function rinex_nav_to_eph(rinex_nav_file,eph_file)
%
%	This function reads in a standard RINEX ephemeris file and converts		%
%	it to the standard binary eph format for ephemeris data.					%
%																				%
% Input parameters:															%
%  rinex_nav_file:  The name of the file containing the RINEx data				%
%	eph_file:  the desired name of the binary output eph file					%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
format long
check = 0;EOH = 0;
fin = fopen(rinex_nav_file,'rt');
prn_vec=zeros(32,1);
max_vec=zeros(32,1);
min_vec=ones(32,1)*999999;
if fin<0
   disp(sprintf('Error opening RINEX nav file %s',rinex_nav_file))
else
   fout = fopen('deletethis.eph','wb');
   
   %% section to avoid the header information %%
   while EOH == 0
      nextline = fgetl(fin);
      headline = sscanf(nextline,'%s');
      check = strmatch('ENDOFHEADER',headline,'exact');
      if length(check)>0
         if check == 1
            EOH = 1;
         end
      end   
   end
   
   while feof(fin) == 0
      nextline_1 = fgetl(fin);
      %line_1 = sscanf(nextline_1,'%e');
      line_1 = get_data_from_line(nextline_1);
      nextline_2 = fgetl(fin);
      %line_2 = sscanf(nextline_2,'%e');
      line_2 = get_data_from_line(nextline_2);
      nextline_3 = fgetl(fin);
      %line_3 = sscanf(nextline_3,'%e');
      line_3 = get_data_from_line(nextline_3);
      nextline_4 = fgetl(fin);
      %line_4 = sscanf(nextline_4,'%e');
      line_4 = get_data_from_line(nextline_4);
      nextline_5 = fgetl(fin);
      %line_5 = sscanf(nextline_5,'%e');
      line_5 = get_data_from_line(nextline_5);
      nextline_6 = fgetl(fin);
      %line_6 = sscanf(nextline_6,'%e');
      line_6 = get_data_from_line(nextline_6);
      nextline_7 = fgetl(fin);
      %line_7 = sscanf(nextline_7,'%e');
      line_7 = get_data_from_line(nextline_7);
      nextline_8 = fgetl(fin);
      %line_8 = sscanf(nextline_8,'%e');
      line_8 = get_data_from_line(nextline_8);
      
      %% Calculate toc  %%
      days = fix(line_4(1)/86400);
      toc = days*86400+line_1(5)*3600+line_1(6)*60+line_1(7);
      %%  store data in an ephemeris array %%
      eph_array = zeros(26,1);
      eph_array(1) = line_1(1);		% sv or prn  ****
      eph_array(2) = line_6(3);		% week       ****
      eph_array(3) = 0;				% zcount
      eph_array(4) = line_7(3);		% tgd
      eph_array(5) = line_7(4);		% aodc (iodc)
      eph_array(6) = toc;				% toc
      eph_array(7) = line_1(10);		% af2 sv clock drift rate
      eph_array(8) = line_1(9);		% af1 sv clock drift
      eph_array(9) = line_1(8);		% af0 sv clock bias
      eph_array(10) = line_2(1);		% aode (iode)
      eph_array(11) = line_2(2);		% crs
      eph_array(12) = line_2(3);		% delta n
      eph_array(13) = line_2(4);		% m0
      eph_array(14) = line_3(1);		% cuc
      eph_array(15) = line_3(2);		% e
      eph_array(16) = line_3(3);		% cus
      eph_array(17) = line_3(4);		% sqrta
      eph_array(18) = line_4(1);		% toe		*******
      eph_array(19) = line_4(2);		% cic
      eph_array(20) = line_4(3);		% Omega0
      eph_array(21) = line_4(4);		% cis
      eph_array(22) = line_5(1);		% i0
      eph_array(23) = line_5(2);		% crc
      eph_array(24) = line_5(3);		% w
      eph_array(25) = line_5(4);		% omegadot
      eph_array(26) = line_6(1);		% idot
      fwrite(fout,eph_array,'double');
      prn=eph_array(1);
      toe=eph_array(18);
      prn_vec(prn)=prn_vec(prn)+1;
      min_vec(prn)=min([min_vec(prn) toe]);
      max_vec(prn)=max([max_vec(prn) toe]);
   end
   
   fclose(fin);
   fclose(fout);
   disp(sprintf('PRN  # of eph records    t0e time range'))
   disp(sprintf('---  ----------------    --------------'))
   for j=1:32
      if prn_vec(j)==0
         disp(sprintf('%2d       %4d                  n/a',j,0))
      else
         disp(sprintf('%2d       %4d            %6d - %6d',j,prn_vec(j),min_vec(j),max_vec(j)))
      end
   end
   
   binary_eph2ascii('deletethis.eph',eph_file)
   delete deletethis.eph
   
end


function data=get_data_from_line(line)

line=strrep(line,'D','E');
data=sscanf(line,'%e');

function binary_eph2ascii(binary_eph_file_name, ascii_eph_file_name)
% function binary_eph2ascii(binary_eph_file_name, ascii_eph_file_name)
%
% Converts from old binary ephemeris file format to new ascii file format

rtd=180/pi;

% Open up the binary ephemeris file
f=fopen(binary_eph_file_name, 'rb', 'ieee-le');

% Open up the ascii ephemeris file
fout = fopen(ascii_eph_file_name, 'w');

% Loop through each record
while ~feof(f)

    % Read in a new raw record
    [raw_ephem, num_read] = fread(f, [26,1], 'double');

    if num_read == 26

        eph.prn      = raw_ephem(1);
        eph.week     = raw_ephem(2);
        eph.t0e      = raw_ephem(18);
        eph.sqrt_a   = raw_ephem(17);
        eph.e        = raw_ephem(15);
        eph.M0       = raw_ephem(13);
        eph.i0       = raw_ephem(22);
        eph.Omega0   = raw_ephem(20);
        eph.omega    = raw_ephem(24);
        eph.idot     = raw_ephem(26);
        eph.Omegadot = raw_ephem(25);
        eph.delta_n  = raw_ephem(12);
        eph.Cuc      = raw_ephem(14);
        eph.Cus      = raw_ephem(16);
        eph.Crc      = raw_ephem(23);
        eph.Crs      = raw_ephem(11);
        eph.Cic      = raw_ephem(19);
        eph.Cis      = raw_ephem(21);
        eph.toc      = raw_ephem(6);
        eph.af0      = raw_ephem(9);
        eph.af1      = raw_ephem(8);
        eph.af2      = raw_ephem(7);
        eph.tgd      = raw_ephem(4);
        eph.valid    = 1;

        fprintf(fout,'********  Ephemeris for PRN-%2d  ********\n',eph.prn);
        fprintf(fout,'PRN:                       %d\n',eph.prn);
        fprintf(fout,'Week:                      %2d\n',eph.week);
        fprintf(fout,'Time of Ephemeris (s):     %2d\n',eph.t0e);
        fprintf(fout,'Square Root of A(m^1/2):   %.8f\n',eph.sqrt_a);
        fprintf(fout,'Eccentricity:              %.15f\n',eph.e);
        fprintf(fout,'Mean anomaly (deg):        %.15f\n',eph.M0*rtd);
        fprintf(fout,'Inclination (deg):         %.15f\n',eph.i0*rtd);
        fprintf(fout,'Long of asc. node (deg):   %.15f\n',eph.Omega0*rtd);
        fprintf(fout,'Argument of perigee (deg): %.15f\n',eph.omega*rtd);
        fprintf(fout,'d-inclination/dt (deg/s):  %.15e\n',eph.idot*rtd);
        fprintf(fout,'d-Omega/dt (deg/s):        %.15e\n',eph.Omegadot*rtd);
        fprintf(fout,'DeltaN (mean motion corr): %.14e\n',eph.delta_n);
        fprintf(fout,'Cuc:                       %.14e\n',eph.Cuc);
        fprintf(fout,'Cus:                       %.14e\n',eph.Cus);
        fprintf(fout,'Crc:                       %.14e\n',eph.Crc);
        fprintf(fout,'Crs:                       %.14e\n',eph.Crs);
        fprintf(fout,'Cic:                       %.14e\n',eph.Cic);
        fprintf(fout,'Cis:                       %.14e\n',eph.Cis);
        fprintf(fout,'Clock ref. time (s):       %2d\n',eph.toc);
        fprintf(fout,'af0:                       %.14e\n',eph.af0);
        fprintf(fout,'af1:                       %.14e\n',eph.af1);
        fprintf(fout,'af2:                       %.14e\n',eph.af2);
        fprintf(fout,'Group delay (s):           %.14e\n\n',eph.tgd);

    end
    
end

fclose(f);
fclose(fout);