source: ntrip/trunk/BNC/src/RTCM3/ephemeris.cpp @ 4891

Last change on this file since 4891 was 4891, checked in by mervart, 8 years ago
File size: 37.3 KB
Line 
1#include <math.h>
2#include <sstream>
3#include <iostream>
4#include <iomanip>
5#include <cstring>
6
7#include <newmatio.h>
8
9#include "ephemeris.h"
10#include "bncutils.h"
11#include "timeutils.h"
12#include "bnctime.h"
13#include "bncapp.h"
14
15using namespace std;
16
17// Returns CRC24
18////////////////////////////////////////////////////////////////////////////
19static unsigned long CRC24(long size, const unsigned char *buf) {
20  unsigned long crc = 0;
21  int ii;
22  while (size--) {
23    crc ^= (*buf++) << (16);
24    for(ii = 0; ii < 8; ii++) {
25      crc <<= 1;
26      if (crc & 0x1000000) {
27        crc ^= 0x01864cfb;
28      }
29    }
30  }
31  return crc;
32}
33
34// Set GPS Satellite Position
35////////////////////////////////////////////////////////////////////////////
36void t_ephGPS::set(const gpsephemeris* ee) {
37
38  _receptDateTime = currentDateAndTimeGPS();
39
40  _prn = QString("G%1").arg(ee->satellite, 2, 10, QChar('0'));
41
42  _TOC.set(ee->GPSweek, ee->TOC);
43  _clock_bias      = ee->clock_bias;
44  _clock_drift     = ee->clock_drift;
45  _clock_driftrate = ee->clock_driftrate;
46
47  _IODE     = ee->IODE;
48  _Crs      = ee->Crs;
49  _Delta_n  = ee->Delta_n;
50  _M0       = ee->M0;
51
52  _Cuc      = ee->Cuc;
53  _e        = ee->e;
54  _Cus      = ee->Cus;
55  _sqrt_A   = ee->sqrt_A;
56
57  _TOEsec   = ee->TOE;
58  _Cic      = ee->Cic;
59  _OMEGA0   = ee->OMEGA0;
60  _Cis      = ee->Cis;
61
62  _i0       = ee->i0;
63  _Crc      = ee->Crc;
64  _omega    = ee->omega;
65  _OMEGADOT = ee->OMEGADOT;
66
67  _IDOT     = ee->IDOT;
68  _L2Codes  = 0.0;
69  _TOEweek  = ee->GPSweek;
70  _L2PFlag  = 0.0;
71
72  if (ee->URAindex <= 6) {
73    _ura = ceil(10.0*pow(2.0, 1.0+((double)ee->URAindex)/2.0))/10.0;
74  }
75  else {
76    _ura = ceil(10.0*pow(2.0, ((double)ee->URAindex)/2.0))/10.0;
77  }
78  _health   = ee->SVhealth;
79  _TGD      = ee->TGD;
80  _IODC     = ee->IODC;
81
82  _TOT         = 0.9999e9;
83  _fitInterval = 0.0;
84
85  _ok       = true;
86}
87
88// Compute GPS Satellite Position (virtual)
89////////////////////////////////////////////////////////////////////////////
90void t_ephGPS::position(int GPSweek, double GPSweeks, 
91                        double* xc,
92                        double* vv) const {
93
94
95  static const double omegaEarth = 7292115.1467e-11;
96  static const double gmWGS      = 398.6005e12;
97
98  memset(xc, 0, 4*sizeof(double));
99  memset(vv, 0, 3*sizeof(double));
100
101  double a0 = _sqrt_A * _sqrt_A;
102  if (a0 == 0) {
103    return;
104  }
105
106  double n0 = sqrt(gmWGS/(a0*a0*a0));
107
108  bncTime tt(GPSweek, GPSweeks);
109  double tk = tt - bncTime(int(_TOEweek), _TOEsec);
110
111  double n  = n0 + _Delta_n;
112  double M  = _M0 + n*tk;
113  double E  = M;
114  double E_last;
115  do {
116    E_last = E;
117    E = M + _e*sin(E);
118  } while ( fabs(E-E_last)*a0 > 0.001 );
119  double v      = 2.0*atan( sqrt( (1.0 + _e)/(1.0 - _e) )*tan( E/2 ) );
120  double u0     = v + _omega;
121  double sin2u0 = sin(2*u0);
122  double cos2u0 = cos(2*u0);
123  double r      = a0*(1 - _e*cos(E)) + _Crc*cos2u0 + _Crs*sin2u0;
124  double i      = _i0 + _IDOT*tk + _Cic*cos2u0 + _Cis*sin2u0;
125  double u      = u0 + _Cuc*cos2u0 + _Cus*sin2u0;
126  double xp     = r*cos(u);
127  double yp     = r*sin(u);
128  double OM     = _OMEGA0 + (_OMEGADOT - omegaEarth)*tk - 
129                   omegaEarth*_TOEsec;
130 
131  double sinom = sin(OM);
132  double cosom = cos(OM);
133  double sini  = sin(i);
134  double cosi  = cos(i);
135  xc[0] = xp*cosom - yp*cosi*sinom;
136  xc[1] = xp*sinom + yp*cosi*cosom;
137  xc[2] = yp*sini;                 
138 
139  double tc = tt - _TOC;
140  xc[3] = _clock_bias + _clock_drift*tc + _clock_driftrate*tc*tc;
141
142  // Velocity
143  // --------
144  double tanv2 = tan(v/2);
145  double dEdM  = 1 / (1 - _e*cos(E));
146  double dotv  = sqrt((1.0 + _e)/(1.0 - _e)) / cos(E/2)/cos(E/2) / (1 + tanv2*tanv2) 
147               * dEdM * n;
148  double dotu  = dotv + (-_Cuc*sin2u0 + _Cus*cos2u0)*2*dotv;
149  double dotom = _OMEGADOT - omegaEarth;
150  double doti  = _IDOT + (-_Cic*sin2u0 + _Cis*cos2u0)*2*dotv;
151  double dotr  = a0 * _e*sin(E) * dEdM * n
152                + (-_Crc*sin2u0 + _Crs*cos2u0)*2*dotv;
153  double dotx  = dotr*cos(u) - r*sin(u)*dotu;
154  double doty  = dotr*sin(u) + r*cos(u)*dotu;
155
156  vv[0]  = cosom   *dotx  - cosi*sinom   *doty      // dX / dr
157           - xp*sinom*dotom - yp*cosi*cosom*dotom   // dX / dOMEGA
158                       + yp*sini*sinom*doti;        // dX / di
159
160  vv[1]  = sinom   *dotx  + cosi*cosom   *doty
161           + xp*cosom*dotom - yp*cosi*sinom*dotom
162                          - yp*sini*cosom*doti;
163
164  vv[2]  = sini    *doty  + yp*cosi      *doti;
165
166  // Relativistic Correction
167  // -----------------------
168  xc[3] -= 2.0 * (xc[0]*vv[0] + xc[1]*vv[1] + xc[2]*vv[2]) / t_CST::c / t_CST::c;
169}
170
171// build up RTCM3 for GPS
172////////////////////////////////////////////////////////////////////////////
173#define GPSTOINT(type, value) static_cast<type>(round(value))
174
175#define GPSADDBITS(a, b) {bitbuffer = (bitbuffer<<(a)) \
176                       |(GPSTOINT(long long,b)&((1ULL<<a)-1)); \
177                       numbits += (a); \
178                       while(numbits >= 8) { \
179                       buffer[size++] = bitbuffer>>(numbits-8);numbits -= 8;}}
180
181#define GPSADDBITSFLOAT(a,b,c) {long long i = GPSTOINT(long long,(b)/(c)); \
182                             GPSADDBITS(a,i)};
183
184int t_ephGPS::RTCM3(unsigned char *buffer) {
185
186  unsigned char *startbuffer = buffer;
187  buffer= buffer+3;
188  int size = 0;
189  int numbits = 0;
190  unsigned long long bitbuffer = 0;
191  if (_ura <= 2.40){
192    _ura = 0;
193  }
194  else if (_ura <= 3.40){
195    _ura = 1;
196  }
197  else if (_ura <= 6.85){
198    _ura = 2;
199  }
200  else if (_ura <= 9.65){
201    _ura = 3;
202  }
203  else if (_ura <= 13.65){
204    _ura = 4;
205  }
206  else if (_ura <= 24.00){
207    _ura = 5;
208  }
209  else if (_ura <= 48.00){
210    _ura = 6;
211  }
212  else if (_ura <= 96.00){
213    _ura = 7;
214  }
215  else if (_ura <= 192.00){
216    _ura = 8;
217  }
218  else if (_ura <= 384.00){
219    _ura = 9;
220  }
221  else if (_ura <= 768.00){
222    _ura = 10;
223  }
224  else if (_ura <= 1536.00){
225    _ura = 11;
226  }
227  else if (_ura <= 1536.00){
228    _ura = 12;
229  }
230  else if (_ura <= 2072.00){
231    _ura = 13;
232  }
233  else if (_ura <= 6144.00){
234    _ura = 14;
235  }
236  else{
237    _ura = 15;
238  }
239
240  GPSADDBITS(12, 1019)
241  GPSADDBITS(6,_prn.right((_prn.length()-1)).toInt())
242  GPSADDBITS(10, _TOC.gpsw())
243  GPSADDBITS(4, _ura)
244  GPSADDBITS(2,_L2Codes)
245  GPSADDBITSFLOAT(14, _IDOT, M_PI/static_cast<double>(1<<30)
246  /static_cast<double>(1<<13))
247  GPSADDBITS(8, _IODE)
248  GPSADDBITS(16, static_cast<int>(_TOC.gpssec())>>4)
249  GPSADDBITSFLOAT(8, _clock_driftrate, 1.0/static_cast<double>(1<<30)
250  /static_cast<double>(1<<25))
251  GPSADDBITSFLOAT(16, _clock_drift, 1.0/static_cast<double>(1<<30)
252  /static_cast<double>(1<<13))
253  GPSADDBITSFLOAT(22, _clock_bias, 1.0/static_cast<double>(1<<30)
254  /static_cast<double>(1<<1))
255  GPSADDBITS(10, _IODC)
256  GPSADDBITSFLOAT(16, _Crs, 1.0/static_cast<double>(1<<5))
257  GPSADDBITSFLOAT(16, _Delta_n, M_PI/static_cast<double>(1<<30)
258  /static_cast<double>(1<<13))
259  GPSADDBITSFLOAT(32, _M0, M_PI/static_cast<double>(1<<30)/static_cast<double>(1<<1))
260  GPSADDBITSFLOAT(16, _Cuc, 1.0/static_cast<double>(1<<29))
261  GPSADDBITSFLOAT(32, _e, 1.0/static_cast<double>(1<<30)/static_cast<double>(1<<3))
262  GPSADDBITSFLOAT(16, _Cus, 1.0/static_cast<double>(1<<29))
263  GPSADDBITSFLOAT(32, _sqrt_A, 1.0/static_cast<double>(1<<19))
264  GPSADDBITS(16, static_cast<int>(_TOEsec)>>4)
265  GPSADDBITSFLOAT(16, _Cic, 1.0/static_cast<double>(1<<29))
266  GPSADDBITSFLOAT(32, _OMEGA0, M_PI/static_cast<double>(1<<30)
267  /static_cast<double>(1<<1))
268  GPSADDBITSFLOAT(16, _Cis, 1.0/static_cast<double>(1<<29))
269  GPSADDBITSFLOAT(32, _i0, M_PI/static_cast<double>(1<<30)/static_cast<double>(1<<1))
270  GPSADDBITSFLOAT(16, _Crc, 1.0/static_cast<double>(1<<5))
271  GPSADDBITSFLOAT(32, _omega, M_PI/static_cast<double>(1<<30)
272  /static_cast<double>(1<<1))
273  GPSADDBITSFLOAT(24, _OMEGADOT, M_PI/static_cast<double>(1<<30)
274  /static_cast<double>(1<<13))
275  GPSADDBITSFLOAT(8, _TGD, 1.0/static_cast<double>(1<<30)/static_cast<double>(1<<1))
276  GPSADDBITS(6, _health) 
277  GPSADDBITS(1, _L2PFlag)
278  GPSADDBITS(1, 0) /* GPS fit interval */
279
280  startbuffer[0]=0xD3;
281  startbuffer[1]=(size >> 8);
282  startbuffer[2]=size;
283  unsigned long  i = CRC24(size+3, startbuffer);
284  buffer[size++] = i >> 16;
285  buffer[size++] = i >> 8;
286  buffer[size++] = i;
287  size += 3;
288  return size;
289}
290
291// Derivative of the state vector using a simple force model (static)
292////////////////////////////////////////////////////////////////////////////
293ColumnVector t_ephGlo::glo_deriv(double /* tt */, const ColumnVector& xv,
294                                 double* acc) {
295
296  // State vector components
297  // -----------------------
298  ColumnVector rr = xv.rows(1,3);
299  ColumnVector vv = xv.rows(4,6);
300
301  // Acceleration
302  // ------------
303  static const double GM    = 398.60044e12;
304  static const double AE    = 6378136.0;
305  static const double OMEGA = 7292115.e-11;
306  static const double C20   = -1082.6257e-6;
307
308  double rho = rr.norm_Frobenius();
309  double t1  = -GM/(rho*rho*rho);
310  double t2  = 3.0/2.0 * C20 * (GM*AE*AE) / (rho*rho*rho*rho*rho);
311  double t3  = OMEGA * OMEGA;
312  double t4  = 2.0 * OMEGA;
313  double z2  = rr(3) * rr(3);
314
315  // Vector of derivatives
316  // ---------------------
317  ColumnVector va(6);
318  va(1) = vv(1);
319  va(2) = vv(2);
320  va(3) = vv(3);
321  va(4) = (t1 + t2*(1.0-5.0*z2/(rho*rho)) + t3) * rr(1) + t4*vv(2) + acc[0]; 
322  va(5) = (t1 + t2*(1.0-5.0*z2/(rho*rho)) + t3) * rr(2) - t4*vv(1) + acc[1]; 
323  va(6) = (t1 + t2*(3.0-5.0*z2/(rho*rho))     ) * rr(3)            + acc[2];
324
325  return va;
326}
327
328// Compute Glonass Satellite Position (virtual)
329////////////////////////////////////////////////////////////////////////////
330void t_ephGlo::position(int GPSweek, double GPSweeks, 
331                        double* xc, double* vv) const {
332
333  static const double nominalStep = 10.0;
334
335  memset(xc, 0, 4*sizeof(double));
336  memset(vv, 0, 3*sizeof(double));
337
338  double dtPos = bncTime(GPSweek, GPSweeks) - _tt;
339
340  int nSteps  = int(fabs(dtPos) / nominalStep) + 1;
341  double step = dtPos / nSteps;
342
343  double acc[3];
344  acc[0] = _x_acceleration * 1.e3;
345  acc[1] = _y_acceleration * 1.e3;
346  acc[2] = _z_acceleration * 1.e3;
347  for (int ii = 1; ii <= nSteps; ii++) { 
348    _xv = rungeKutta4(_tt.gpssec(), _xv, step, acc, glo_deriv);
349    _tt = _tt + step;
350  }
351
352  // Position and Velocity
353  // ---------------------
354  xc[0] = _xv(1);
355  xc[1] = _xv(2);
356  xc[2] = _xv(3);
357
358  vv[0] = _xv(4);
359  vv[1] = _xv(5);
360  vv[2] = _xv(6);
361
362  // Clock Correction
363  // ----------------
364  double dtClk = bncTime(GPSweek, GPSweeks) - _TOC;
365  xc[3] = -_tau + _gamma * dtClk;
366}
367
368// IOD of Glonass Ephemeris (virtual)
369////////////////////////////////////////////////////////////////////////////
370int t_ephGlo::IOD() const {
371  bncTime tMoscow = _TOC - _gps_utc + 3 * 3600.0;
372  return int(tMoscow.daysec() / 900);
373}
374
375// Set Glonass Ephemeris
376////////////////////////////////////////////////////////////////////////////
377void t_ephGlo::set(const glonassephemeris* ee) {
378
379  _receptDateTime = currentDateAndTimeGPS();
380
381  _prn = QString("R%1").arg(ee->almanac_number, 2, 10, QChar('0'));
382
383  int ww  = ee->GPSWeek;
384  int tow = ee->GPSTOW; 
385  updatetime(&ww, &tow, ee->tb*1000, 0);  // Moscow -> GPS
386
387  // Check the day once more
388  // -----------------------
389  {
390    const double secPerDay  = 24 * 3600.0;
391    const double secPerWeek = 7 * secPerDay;
392    int ww_old  = ww;
393    int tow_old = tow;
394    int    currentWeek;
395    double currentSec;
396    currentGPSWeeks(currentWeek, currentSec);
397    bncTime currentTime(currentWeek, currentSec);
398    bncTime hTime(ww, (double) tow);
399
400    bool changed = false;
401    if      (hTime - currentTime >  secPerDay/2.0) {
402      changed = true;
403      tow -= int(secPerDay);
404      if (tow < 0) {
405        tow += int(secPerWeek);
406        ww  -= 1;
407      }
408    }
409    else if (hTime - currentTime < -secPerDay/2.0) {
410      changed = true;
411      tow += int(secPerDay);
412      if (tow > secPerWeek) {
413        tow -= int(secPerWeek);
414        ww  += 1;
415      }
416    }
417
418    if (changed && ((bncApp*) qApp)->mode() == bncApp::batchPostProcessing) {
419      bncTime newHTime(ww, (double) tow);
420      cout << "GLONASS " << ee->almanac_number <<  " Time Changed at " 
421           << currentTime.datestr()         << " " << currentTime.timestr() 
422           << endl
423           << "old: " << hTime.datestr()    << " " << hTime.timestr()       
424           << endl
425           << "new: " << newHTime.datestr() << " " << newHTime.timestr()   
426           << endl
427           << "eph: " << ee->GPSWeek << " " << ee->GPSTOW << " " << ee->tb
428           << endl
429           << "ww, tow (old): " << ww_old << " " << tow_old
430           << endl
431           << "ww, tow (new): " << ww     << " " << tow
432           << endl << endl;
433    }
434  }
435
436  bncTime hlpTime(ww, (double) tow);
437  unsigned year, month, day;
438  hlpTime.civil_date(year, month, day);
439  _gps_utc = gnumleap(year, month, day);
440
441  _TOC.set(ww, tow);
442  _E                 = ee->E;
443  _tau               = ee->tau;
444  _gamma             = ee->gamma;
445  _x_pos             = ee->x_pos;
446  _x_velocity        = ee->x_velocity;     
447  _x_acceleration    = ee->x_acceleration;
448  _y_pos             = ee->y_pos;         
449  _y_velocity        = ee->y_velocity;   
450  _y_acceleration    = ee->y_acceleration;
451  _z_pos             = ee->z_pos;         
452  _z_velocity        = ee->z_velocity;   
453  _z_acceleration    = ee->z_acceleration;
454  _health            = 0;
455  _frequency_number  = ee->frequency_number;
456  _tki               = ee->tk-3*60*60; if (_tki < 0) _tki += 86400;
457
458  // Initialize status vector
459  // ------------------------
460  _tt = _TOC;
461
462  _xv(1) = _x_pos * 1.e3; 
463  _xv(2) = _y_pos * 1.e3; 
464  _xv(3) = _z_pos * 1.e3; 
465  _xv(4) = _x_velocity * 1.e3; 
466  _xv(5) = _y_velocity * 1.e3; 
467  _xv(6) = _z_velocity * 1.e3; 
468
469  _ok = true;
470}
471
472// build up RTCM3 for GLONASS
473////////////////////////////////////////////////////////////////////////////
474#define GLONASSTOINT(type, value) static_cast<type>(round(value))
475
476#define GLONASSADDBITS(a, b) {bitbuffer = (bitbuffer<<(a)) \
477                       |(GLONASSTOINT(long long,b)&((1ULL<<(a))-1)); \
478                       numbits += (a); \
479                       while(numbits >= 8) { \
480                       buffer[size++] = bitbuffer>>(numbits-8);numbits -= 8;}}
481#define GLONASSADDBITSFLOATM(a,b,c) {int s; long long i; \
482                       if(b < 0.0) \
483                       { \
484                         s = 1; \
485                         i = GLONASSTOINT(long long,(-b)/(c)); \
486                         if(!i) s = 0; \
487                       } \
488                       else \
489                       { \
490                         s = 0; \
491                         i = GLONASSTOINT(long long,(b)/(c)); \
492                       } \
493                       GLONASSADDBITS(1,s) \
494                       GLONASSADDBITS(a-1,i)}
495
496int t_ephGlo::RTCM3(unsigned char *buffer)
497{
498
499  int size = 0;
500  int numbits = 0;
501  long long bitbuffer = 0;
502  unsigned char *startbuffer = buffer;
503  buffer= buffer+3;
504
505  GLONASSADDBITS(12, 1020)
506  GLONASSADDBITS(6, _prn.right((_prn.length()-1)).toInt())
507  GLONASSADDBITS(5, 7+_frequency_number)
508  GLONASSADDBITS(1, 0)
509  GLONASSADDBITS(1, 0)
510  GLONASSADDBITS(2, 0)
511  _tki=_tki+3*60*60;
512  GLONASSADDBITS(5, static_cast<int>(_tki)/(60*60))
513  GLONASSADDBITS(6, (static_cast<int>(_tki)/60)%60)
514  GLONASSADDBITS(1, (static_cast<int>(_tki)/30)%30)
515  GLONASSADDBITS(1, _health) 
516  GLONASSADDBITS(1, 0)
517  unsigned long long timeofday = (static_cast<int>(_tt.gpssec()+3*60*60-_gps_utc)%86400);
518  GLONASSADDBITS(7, timeofday/60/15)
519  GLONASSADDBITSFLOATM(24, _x_velocity*1000, 1000.0/static_cast<double>(1<<20))
520  GLONASSADDBITSFLOATM(27, _x_pos*1000, 1000.0/static_cast<double>(1<<11))
521  GLONASSADDBITSFLOATM(5, _x_acceleration*1000, 1000.0/static_cast<double>(1<<30))
522  GLONASSADDBITSFLOATM(24, _y_velocity*1000, 1000.0/static_cast<double>(1<<20))
523  GLONASSADDBITSFLOATM(27, _y_pos*1000, 1000.0/static_cast<double>(1<<11))
524  GLONASSADDBITSFLOATM(5, _y_acceleration*1000, 1000.0/static_cast<double>(1<<30))
525  GLONASSADDBITSFLOATM(24, _z_velocity*1000, 1000.0/static_cast<double>(1<<20))
526  GLONASSADDBITSFLOATM(27,_z_pos*1000, 1000.0/static_cast<double>(1<<11))
527  GLONASSADDBITSFLOATM(5, _z_acceleration*1000, 1000.0/static_cast<double>(1<<30))
528  GLONASSADDBITS(1, 0)
529  GLONASSADDBITSFLOATM(11, _gamma, 1.0/static_cast<double>(1<<30)
530  /static_cast<double>(1<<10))
531  GLONASSADDBITS(2, 0) /* GLONASS-M P */
532  GLONASSADDBITS(1, 0) /* GLONASS-M ln(3) */
533  GLONASSADDBITSFLOATM(22, _tau, 1.0/static_cast<double>(1<<30))
534  GLONASSADDBITS(5, 0) /* GLONASS-M delta tau */
535  GLONASSADDBITS(5, _E)
536  GLONASSADDBITS(1, 0) /* GLONASS-M P4 */
537  GLONASSADDBITS(4, 0) /* GLONASS-M FT */
538  GLONASSADDBITS(11, 0) /* GLONASS-M NT */
539  GLONASSADDBITS(2, 0) /* GLONASS-M active? */
540  GLONASSADDBITS(1, 0) /* GLONASS additional data */
541  GLONASSADDBITS(11, 0) /* GLONASS NA */
542  GLONASSADDBITS(32, 0) /* GLONASS tau C */
543  GLONASSADDBITS(5, 0) /* GLONASS-M N4 */
544  GLONASSADDBITS(22, 0) /* GLONASS-M tau GPS */
545  GLONASSADDBITS(1, 0) /* GLONASS-M ln(5) */
546  GLONASSADDBITS(7, 0) /* Reserved */
547
548  startbuffer[0]=0xD3;
549  startbuffer[1]=(size >> 8);
550  startbuffer[2]=size;
551  unsigned long i = CRC24(size+3, startbuffer);
552  buffer[size++] = i >> 16;
553  buffer[size++] = i >> 8;
554  buffer[size++] = i;
555  size += 3;
556  return size;
557}
558
559// Set Galileo Satellite Position
560////////////////////////////////////////////////////////////////////////////
561void t_ephGal::set(const galileoephemeris* ee) {
562
563  _receptDateTime = currentDateAndTimeGPS();
564
565  _prn = QString("E%1").arg(ee->satellite, 2, 10, QChar('0'));
566
567  _TOC.set(ee->Week, ee->TOC);
568  _clock_bias      = ee->clock_bias;
569  _clock_drift     = ee->clock_drift;
570  _clock_driftrate = ee->clock_driftrate;
571
572  _IODnav   = ee->IODnav;
573  _Crs      = ee->Crs;
574  _Delta_n  = ee->Delta_n;
575  _M0       = ee->M0;
576
577  _Cuc      = ee->Cuc;
578  _e        = ee->e;
579  _Cus      = ee->Cus;
580  _sqrt_A   = ee->sqrt_A;
581
582  _TOEsec   = ee->TOE;
583  _Cic      = ee->Cic;
584  _OMEGA0   = ee->OMEGA0;
585  _Cis      = ee->Cis;
586
587  _i0       = ee->i0;
588  _Crc      = ee->Crc;
589  _omega    = ee->omega;
590  _OMEGADOT = ee->OMEGADOT;
591
592  _IDOT     = ee->IDOT;
593  _TOEweek  = ee->Week;
594
595  _SISA     = ee->SISA;
596  _E5aHS    = ee->E5aHS;
597  _BGD_1_5A = ee->BGD_1_5A;
598  _BGD_1_5B = ee->BGD_1_5B;
599
600  _TOT      = 0.9999e9;
601
602  _ok = true;
603}
604
605// Compute Galileo Satellite Position (virtual)
606////////////////////////////////////////////////////////////////////////////
607void t_ephGal::position(int GPSweek, double GPSweeks, 
608                        double* xc,
609                        double* vv) const {
610
611  static const double omegaEarth = 7292115.1467e-11;
612  static const double gmWGS      = 398.6005e12;
613
614  memset(xc, 0, 4*sizeof(double));
615  memset(vv, 0, 3*sizeof(double));
616
617  double a0 = _sqrt_A * _sqrt_A;
618  if (a0 == 0) {
619    return;
620  }
621
622  double n0 = sqrt(gmWGS/(a0*a0*a0));
623
624  bncTime tt(GPSweek, GPSweeks);
625  double tk = tt - bncTime(_TOC.gpsw(), _TOEsec);
626
627  double n  = n0 + _Delta_n;
628  double M  = _M0 + n*tk;
629  double E  = M;
630  double E_last;
631  do {
632    E_last = E;
633    E = M + _e*sin(E);
634  } while ( fabs(E-E_last)*a0 > 0.001 );
635  double v      = 2.0*atan( sqrt( (1.0 + _e)/(1.0 - _e) )*tan( E/2 ) );
636  double u0     = v + _omega;
637  double sin2u0 = sin(2*u0);
638  double cos2u0 = cos(2*u0);
639  double r      = a0*(1 - _e*cos(E)) + _Crc*cos2u0 + _Crs*sin2u0;
640  double i      = _i0 + _IDOT*tk + _Cic*cos2u0 + _Cis*sin2u0;
641  double u      = u0 + _Cuc*cos2u0 + _Cus*sin2u0;
642  double xp     = r*cos(u);
643  double yp     = r*sin(u);
644  double OM     = _OMEGA0 + (_OMEGADOT - omegaEarth)*tk - 
645                  omegaEarth*_TOEsec;
646 
647  double sinom = sin(OM);
648  double cosom = cos(OM);
649  double sini  = sin(i);
650  double cosi  = cos(i);
651  xc[0] = xp*cosom - yp*cosi*sinom;
652  xc[1] = xp*sinom + yp*cosi*cosom;
653  xc[2] = yp*sini;                 
654 
655  double tc = tt - _TOC;
656  xc[3] = _clock_bias + _clock_drift*tc + _clock_driftrate*tc*tc;
657
658  // Velocity
659  // --------
660  double tanv2 = tan(v/2);
661  double dEdM  = 1 / (1 - _e*cos(E));
662  double dotv  = sqrt((1.0 + _e)/(1.0 - _e)) / cos(E/2)/cos(E/2) / (1 + tanv2*tanv2) 
663               * dEdM * n;
664  double dotu  = dotv + (-_Cuc*sin2u0 + _Cus*cos2u0)*2*dotv;
665  double dotom = _OMEGADOT - omegaEarth;
666  double doti  = _IDOT + (-_Cic*sin2u0 + _Cis*cos2u0)*2*dotv;
667  double dotr  = a0 * _e*sin(E) * dEdM * n
668                + (-_Crc*sin2u0 + _Crs*cos2u0)*2*dotv;
669  double dotx  = dotr*cos(u) - r*sin(u)*dotu;
670  double doty  = dotr*sin(u) + r*cos(u)*dotu;
671
672  vv[0]  = cosom   *dotx  - cosi*sinom   *doty      // dX / dr
673           - xp*sinom*dotom - yp*cosi*cosom*dotom   // dX / dOMEGA
674                       + yp*sini*sinom*doti;        // dX / di
675
676  vv[1]  = sinom   *dotx  + cosi*cosom   *doty
677           + xp*cosom*dotom - yp*cosi*sinom*dotom
678                          - yp*sini*cosom*doti;
679
680  vv[2]  = sini    *doty  + yp*cosi      *doti;
681
682  // Relativistic Correction
683  // -----------------------
684  //  xc(4) -= 4.442807633e-10 * _e * sqrt(a0) *sin(E);
685  xc[3] -= 2.0 * (xc[0]*vv[0] + xc[1]*vv[1] + xc[2]*vv[2]) / t_CST::c / t_CST::c;
686}
687
688// build up RTCM3 for Galileo
689////////////////////////////////////////////////////////////////////////////
690#define GALILEOTOINT(type, value) static_cast<type>(round(value))
691
692#define GALILEOADDBITS(a, b) {bitbuffer = (bitbuffer<<(a)) \
693                       |(GALILEOTOINT(long long,b)&((1LL<<a)-1)); \
694                       numbits += (a); \
695                       while(numbits >= 8) { \
696                       buffer[size++] = bitbuffer>>(numbits-8);numbits -= 8;}}
697#define GALILEOADDBITSFLOAT(a,b,c) {long long i = GALILEOTOINT(long long,(b)/(c)); \
698                             GALILEOADDBITS(a,i)};
699
700int t_ephGal::RTCM3(unsigned char *buffer) {
701  int size = 0;
702  int numbits = 0;
703  long long bitbuffer = 0;
704  unsigned char *startbuffer = buffer;
705  buffer= buffer+3;
706
707  GALILEOADDBITS(12, /*inav ? 1046 :*/ 1045)
708  GALILEOADDBITS(6, _prn.right((_prn.length()-1)).toInt())
709  GALILEOADDBITS(12, _TOC.gpsw())
710  GALILEOADDBITS(10, _IODnav)
711  GALILEOADDBITS(8, _SISA)
712  GALILEOADDBITSFLOAT(14, _IDOT, M_PI/static_cast<double>(1<<30)
713  /static_cast<double>(1<<13))
714  GALILEOADDBITS(14, _TOC.gpssec()/60)
715  GALILEOADDBITSFLOAT(6, _clock_driftrate, 1.0/static_cast<double>(1<<30)
716  /static_cast<double>(1<<29))
717  GALILEOADDBITSFLOAT(21, _clock_drift, 1.0/static_cast<double>(1<<30)
718  /static_cast<double>(1<<16))
719  GALILEOADDBITSFLOAT(31, _clock_bias, 1.0/static_cast<double>(1<<30)
720  /static_cast<double>(1<<4))
721  GALILEOADDBITSFLOAT(16, _Crs, 1.0/static_cast<double>(1<<5))
722  GALILEOADDBITSFLOAT(16, _Delta_n, M_PI/static_cast<double>(1<<30)
723  /static_cast<double>(1<<13))
724  GALILEOADDBITSFLOAT(32, _M0, M_PI/static_cast<double>(1<<30)/static_cast<double>(1<<1))
725  GALILEOADDBITSFLOAT(16, _Cuc, 1.0/static_cast<double>(1<<29))
726  GALILEOADDBITSFLOAT(32, _e, 1.0/static_cast<double>(1<<30)/static_cast<double>(1<<3))
727  GALILEOADDBITSFLOAT(16, _Cus, 1.0/static_cast<double>(1<<29))
728  GALILEOADDBITSFLOAT(32, _sqrt_A, 1.0/static_cast<double>(1<<19))
729  GALILEOADDBITS(14, _TOEsec/60)
730  GALILEOADDBITSFLOAT(16, _Cic, 1.0/static_cast<double>(1<<29))
731  GALILEOADDBITSFLOAT(32, _OMEGA0, M_PI/static_cast<double>(1<<30)
732  /static_cast<double>(1<<1))
733  GALILEOADDBITSFLOAT(16, _Cis, 1.0/static_cast<double>(1<<29))
734  GALILEOADDBITSFLOAT(32, _i0, M_PI/static_cast<double>(1<<30)/static_cast<double>(1<<1))
735  GALILEOADDBITSFLOAT(16, _Crc, 1.0/static_cast<double>(1<<5))
736  GALILEOADDBITSFLOAT(32, _omega, M_PI/static_cast<double>(1<<30)
737  /static_cast<double>(1<<1))
738  GALILEOADDBITSFLOAT(24, _OMEGADOT, M_PI/static_cast<double>(1<<30)
739  /static_cast<double>(1<<13))
740  GALILEOADDBITSFLOAT(10, _BGD_1_5A, 1.0/static_cast<double>(1<<30)
741  /static_cast<double>(1<<2))
742  /*if(inav)
743  {
744    GALILEOADDBITSFLOAT(10, _BGD_1_5B, 1.0/static_cast<double>(1<<30)
745    /static_cast<double>(1<<2))
746    GALILEOADDBITS(2, _E5bHS)
747    GALILEOADDBITS(1, flags & MNFGALEPHF_E5BDINVALID)
748  }
749  else*/
750  {
751    GALILEOADDBITS(2, _E5aHS)
752    GALILEOADDBITS(1, /*flags & MNFGALEPHF_E5ADINVALID*/0)
753  }
754  _TOEsec = 0.9999E9;
755  GALILEOADDBITS(20, _TOEsec)
756
757  GALILEOADDBITS(/*inav ? 1 :*/ 3, 0) /* fill up */
758
759  startbuffer[0]=0xD3;
760  startbuffer[1]=(size >> 8);
761  startbuffer[2]=size;
762  unsigned long i = CRC24(size+3, startbuffer);
763  buffer[size++] = i >> 16;
764  buffer[size++] = i >> 8;
765  buffer[size++] = i;
766  size += 3;
767  return size;
768}
769
770// Constructor
771//////////////////////////////////////////////////////////////////////////////
772t_ephGPS::t_ephGPS(float rnxVersion, const QStringList& lines) {
773
774  const int nLines = 8;
775
776  _ok = false;
777
778  if (lines.size() != nLines) {
779    return;
780  }
781
782  // RINEX Format
783  // ------------
784  int fieldLen = 19;
785
786  int pos[4];
787  pos[0] = (rnxVersion <= 2.12) ?  3 :  4;
788  pos[1] = pos[0] + fieldLen;
789  pos[2] = pos[1] + fieldLen;
790  pos[3] = pos[2] + fieldLen;
791
792  // Read eight lines
793  // ----------------
794  for (int iLine = 0; iLine < nLines; iLine++) {
795    QString line = lines[iLine];
796
797    if      ( iLine == 0 ) {
798      QTextStream in(line.left(pos[1]).toAscii());
799
800      int    year, month, day, hour, min;
801      double sec;
802     
803      in >> _prn >> year >> month >> day >> hour >> min >> sec;
804
805      if (_prn.at(0) != 'G') {
806        _prn = QString("G%1").arg(_prn.toInt(), 2, 10, QLatin1Char('0'));
807      }
808   
809      if      (year <  80) {
810        year += 2000;
811      }
812      else if (year < 100) {
813        year += 1900;
814      }
815
816      _TOC.set(year, month, day, hour, min, sec);
817
818      if ( readDbl(line, pos[1], fieldLen, _clock_bias     ) ||
819           readDbl(line, pos[2], fieldLen, _clock_drift    ) ||
820           readDbl(line, pos[3], fieldLen, _clock_driftrate) ) {
821        return;
822      }
823    }
824
825    else if      ( iLine == 1 ) {
826      if ( readDbl(line, pos[0], fieldLen, _IODE   ) ||
827           readDbl(line, pos[1], fieldLen, _Crs    ) ||
828           readDbl(line, pos[2], fieldLen, _Delta_n) ||
829           readDbl(line, pos[3], fieldLen, _M0     ) ) {
830        return;
831      }
832    }
833
834    else if ( iLine == 2 ) {
835      if ( readDbl(line, pos[0], fieldLen, _Cuc   ) ||
836           readDbl(line, pos[1], fieldLen, _e     ) ||
837           readDbl(line, pos[2], fieldLen, _Cus   ) ||
838           readDbl(line, pos[3], fieldLen, _sqrt_A) ) {
839        return;
840      }
841    }
842
843    else if ( iLine == 3 ) {
844      if ( readDbl(line, pos[0], fieldLen, _TOEsec)  ||
845           readDbl(line, pos[1], fieldLen, _Cic   )  ||
846           readDbl(line, pos[2], fieldLen, _OMEGA0)  ||
847           readDbl(line, pos[3], fieldLen, _Cis   ) ) {
848        return;
849      }
850    }
851
852    else if ( iLine == 4 ) {
853      if ( readDbl(line, pos[0], fieldLen, _i0      ) ||
854           readDbl(line, pos[1], fieldLen, _Crc     ) ||
855           readDbl(line, pos[2], fieldLen, _omega   ) ||
856           readDbl(line, pos[3], fieldLen, _OMEGADOT) ) {
857        return;
858      }
859    }
860
861    else if ( iLine == 5 ) {
862      if ( readDbl(line, pos[0], fieldLen, _IDOT   ) ||
863           readDbl(line, pos[1], fieldLen, _L2Codes) ||
864           readDbl(line, pos[2], fieldLen, _TOEweek  ) ||
865           readDbl(line, pos[3], fieldLen, _L2PFlag) ) {
866        return;
867      }
868    }
869
870    else if ( iLine == 6 ) {
871      if ( readDbl(line, pos[0], fieldLen, _ura   ) ||
872           readDbl(line, pos[1], fieldLen, _health) ||
873           readDbl(line, pos[2], fieldLen, _TGD   ) ||
874           readDbl(line, pos[3], fieldLen, _IODC  ) ) {
875        return;
876      }
877    }
878
879    else if ( iLine == 7 ) {
880      if ( readDbl(line, pos[0], fieldLen, _TOT) ) {
881        return;
882      }
883      readDbl(line, pos[1], fieldLen, _fitInterval); // _fitInterval optional
884    }
885  }
886
887  _ok = true;
888}
889
890// Constructor
891//////////////////////////////////////////////////////////////////////////////
892t_ephGlo::t_ephGlo(float rnxVersion, const QStringList& lines) {
893
894  const int nLines = 4;
895
896  _ok = false;
897
898  if (lines.size() != nLines) {
899    return;
900  }
901
902  // RINEX Format
903  // ------------
904  int fieldLen = 19;
905
906  int pos[4];
907  pos[0] = (rnxVersion <= 2.12) ?  3 :  4;
908  pos[1] = pos[0] + fieldLen;
909  pos[2] = pos[1] + fieldLen;
910  pos[3] = pos[2] + fieldLen;
911
912  // Read four lines
913  // ---------------
914  for (int iLine = 0; iLine < nLines; iLine++) {
915    QString line = lines[iLine];
916
917    if      ( iLine == 0 ) {
918      QTextStream in(line.left(pos[1]).toAscii());
919
920      int    year, month, day, hour, min;
921      double sec;
922     
923      in >> _prn >> year >> month >> day >> hour >> min >> sec;
924
925      if (_prn.at(0) != 'R') {
926        _prn = QString("R%1").arg(_prn.toInt(), 2, 10, QLatin1Char('0'));
927      }
928   
929      if      (year <  80) {
930        year += 2000;
931      }
932      else if (year < 100) {
933        year += 1900;
934      }
935
936      _gps_utc = gnumleap(year, month, day);
937
938      _TOC.set(year, month, day, hour, min, sec);
939      _TOC  = _TOC + _gps_utc;
940
941      if ( readDbl(line, pos[1], fieldLen, _tau  ) ||
942           readDbl(line, pos[2], fieldLen, _gamma) ||
943           readDbl(line, pos[3], fieldLen, _tki  ) ) {
944        return;
945      }
946
947      _tau = -_tau;
948    }
949
950    else if      ( iLine == 1 ) {
951      if ( readDbl(line, pos[0], fieldLen, _x_pos         ) ||
952           readDbl(line, pos[1], fieldLen, _x_velocity    ) ||
953           readDbl(line, pos[2], fieldLen, _x_acceleration) ||
954           readDbl(line, pos[3], fieldLen, _health        ) ) {
955        return;
956      }
957    }
958
959    else if ( iLine == 2 ) {
960      if ( readDbl(line, pos[0], fieldLen, _y_pos           ) ||
961           readDbl(line, pos[1], fieldLen, _y_velocity      ) ||
962           readDbl(line, pos[2], fieldLen, _y_acceleration  ) ||
963           readDbl(line, pos[3], fieldLen, _frequency_number) ) {
964        return;
965      }
966    }
967
968    else if ( iLine == 3 ) {
969      if ( readDbl(line, pos[0], fieldLen, _z_pos         )  ||
970           readDbl(line, pos[1], fieldLen, _z_velocity    )  ||
971           readDbl(line, pos[2], fieldLen, _z_acceleration)  ||
972           readDbl(line, pos[3], fieldLen, _E             ) ) {
973        return;
974      }
975    }
976  }
977
978  // Initialize status vector
979  // ------------------------
980  _tt = _TOC;
981  _xv.ReSize(6); 
982  _xv(1) = _x_pos * 1.e3; 
983  _xv(2) = _y_pos * 1.e3; 
984  _xv(3) = _z_pos * 1.e3; 
985  _xv(4) = _x_velocity * 1.e3; 
986  _xv(5) = _y_velocity * 1.e3; 
987  _xv(6) = _z_velocity * 1.e3; 
988
989  _ok = true;
990}
991
992// Constructor
993//////////////////////////////////////////////////////////////////////////////
994t_ephGal::t_ephGal(float rnxVersion, const QStringList& lines) {
995
996  const int nLines = 8;
997
998  _ok = false;
999
1000  if (lines.size() != nLines) {
1001    return;
1002  }
1003
1004  // RINEX Format
1005  // ------------
1006  int fieldLen = 19;
1007
1008  int pos[4];
1009  pos[0] = (rnxVersion <= 2.12) ?  3 :  4;
1010  pos[1] = pos[0] + fieldLen;
1011  pos[2] = pos[1] + fieldLen;
1012  pos[3] = pos[2] + fieldLen;
1013
1014  // Read eight lines
1015  // ----------------
1016  for (int iLine = 0; iLine < nLines; iLine++) {
1017    QString line = lines[iLine];
1018
1019    if      ( iLine == 0 ) {
1020      QTextStream in(line.left(pos[1]).toAscii());
1021
1022      int    year, month, day, hour, min;
1023      double sec;
1024     
1025      in >> _prn >> year >> month >> day >> hour >> min >> sec;
1026
1027      if (_prn.at(0) != 'E') {
1028        _prn = QString("E%1").arg(_prn.toInt(), 2, 10, QLatin1Char('0'));
1029      }
1030   
1031      if      (year <  80) {
1032        year += 2000;
1033      }
1034      else if (year < 100) {
1035        year += 1900;
1036      }
1037
1038      _TOC.set(year, month, day, hour, min, sec);
1039
1040      if ( readDbl(line, pos[1], fieldLen, _clock_bias     ) ||
1041           readDbl(line, pos[2], fieldLen, _clock_drift    ) ||
1042           readDbl(line, pos[3], fieldLen, _clock_driftrate) ) {
1043        return;
1044      }
1045    }
1046
1047    else if      ( iLine == 1 ) {
1048      if ( readDbl(line, pos[0], fieldLen, _IODnav ) ||
1049           readDbl(line, pos[1], fieldLen, _Crs    ) ||
1050           readDbl(line, pos[2], fieldLen, _Delta_n) ||
1051           readDbl(line, pos[3], fieldLen, _M0     ) ) {
1052        return;
1053      }
1054    }
1055
1056    else if ( iLine == 2 ) {
1057      if ( readDbl(line, pos[0], fieldLen, _Cuc   ) ||
1058           readDbl(line, pos[1], fieldLen, _e     ) ||
1059           readDbl(line, pos[2], fieldLen, _Cus   ) ||
1060           readDbl(line, pos[3], fieldLen, _sqrt_A) ) {
1061        return;
1062      }
1063    }
1064
1065    else if ( iLine == 3 ) {
1066      if ( readDbl(line, pos[0], fieldLen, _TOEsec)  ||
1067           readDbl(line, pos[1], fieldLen, _Cic   )  ||
1068           readDbl(line, pos[2], fieldLen, _OMEGA0)  ||
1069           readDbl(line, pos[3], fieldLen, _Cis   ) ) {
1070        return;
1071      }
1072    }
1073
1074    else if ( iLine == 4 ) {
1075      if ( readDbl(line, pos[0], fieldLen, _i0      ) ||
1076           readDbl(line, pos[1], fieldLen, _Crc     ) ||
1077           readDbl(line, pos[2], fieldLen, _omega   ) ||
1078           readDbl(line, pos[3], fieldLen, _OMEGADOT) ) {
1079        return;
1080      }
1081    }
1082
1083    else if ( iLine == 5 ) {
1084      if ( readDbl(line, pos[0], fieldLen, _IDOT    ) ||
1085           readDbl(line, pos[2], fieldLen, _TOEweek) ) {
1086        return;
1087      }
1088    }
1089
1090    else if ( iLine == 6 ) {
1091      if ( readDbl(line, pos[0], fieldLen, _SISA    ) ||
1092           readDbl(line, pos[1], fieldLen, _E5aHS   ) ||
1093           readDbl(line, pos[2], fieldLen, _BGD_1_5A) ||
1094           readDbl(line, pos[3], fieldLen, _BGD_1_5B) ) {
1095        return;
1096      }
1097    }
1098
1099    else if ( iLine == 7 ) {
1100      if ( readDbl(line, pos[0], fieldLen, _TOT) ) {
1101        return;
1102      }
1103    }
1104  }
1105
1106  _ok = true;
1107}
1108
1109//
1110//////////////////////////////////////////////////////////////////////////////
1111QString t_eph::rinexDateStr(const bncTime& tt, const QString& prn,
1112                            double version) {
1113
1114  QString datStr;
1115 
1116  unsigned year, month, day, hour, min;
1117  double   sec;
1118  tt.civil_date(year, month, day);
1119  tt.civil_time(hour, min, sec);
1120 
1121  QTextStream out(&datStr);
1122
1123  if (version < 3.0) {
1124    QString prnHlp = prn.mid(1,2); if (prnHlp[0] == '0') prnHlp[0] = ' ';
1125    out << prnHlp << QString(" %1 %2 %3 %4 %5%6")
1126      .arg(year % 100, 2, 10, QChar('0'))
1127      .arg(month,      2)
1128      .arg(day,        2)
1129      .arg(hour,       2)
1130      .arg(min,        2)
1131      .arg(sec, 5, 'f',1);
1132  }
1133  else {
1134    out << prn << QString(" %1 %2 %3 %4 %5 %6")
1135      .arg(year,     4)
1136      .arg(month,    2, 10, QChar('0'))
1137      .arg(day,      2, 10, QChar('0'))
1138      .arg(hour,     2, 10, QChar('0'))
1139      .arg(min,      2, 10, QChar('0'))
1140      .arg(int(sec), 2, 10, QChar('0'));
1141  }
1142
1143  return datStr;
1144}
1145
1146// RINEX Format String
1147//////////////////////////////////////////////////////////////////////////////
1148QString t_ephGPS::toString(double version) const {
1149
1150  QString rnxStr = rinexDateStr(_TOC, _prn, version);
1151 
1152  QTextStream out(&rnxStr);
1153 
1154  out << QString("%1%2%3\n")
1155    .arg(_clock_bias,      19, 'e', 12)
1156    .arg(_clock_drift,     19, 'e', 12)
1157    .arg(_clock_driftrate, 19, 'e', 12);
1158
1159  QString fmt = version < 3.0 ? "   %1%2%3%4\n" : "    %1%2%3%4\n";
1160
1161  out << QString(fmt)
1162    .arg(_IODE,    19, 'e', 12)
1163    .arg(_Crs,     19, 'e', 12)
1164    .arg(_Delta_n, 19, 'e', 12)
1165    .arg(_M0,      19, 'e', 12);
1166
1167  out << QString(fmt)
1168    .arg(_Cuc,    19, 'e', 12)
1169    .arg(_e,      19, 'e', 12)
1170    .arg(_Cus,    19, 'e', 12)
1171    .arg(_sqrt_A, 19, 'e', 12);
1172
1173  out << QString(fmt)
1174    .arg(_TOEsec, 19, 'e', 12)
1175    .arg(_Cic,    19, 'e', 12)
1176    .arg(_OMEGA0, 19, 'e', 12)
1177    .arg(_Cis,    19, 'e', 12);
1178
1179  out << QString(fmt)
1180    .arg(_i0,       19, 'e', 12)
1181    .arg(_Crc,      19, 'e', 12)
1182    .arg(_omega,    19, 'e', 12)
1183    .arg(_OMEGADOT, 19, 'e', 12);
1184
1185  out << QString(fmt)
1186    .arg(_IDOT,    19, 'e', 12)
1187    .arg(_L2Codes, 19, 'e', 12)
1188    .arg(_TOEweek, 19, 'e', 12)
1189    .arg(_L2PFlag, 19, 'e', 12);
1190
1191  out << QString(fmt)
1192    .arg(_ura,    19, 'e', 12)
1193    .arg(_health, 19, 'e', 12)
1194    .arg(_TGD,    19, 'e', 12)
1195    .arg(_IODC,   19, 'e', 12);
1196
1197  out << QString(fmt)
1198    .arg(_TOT,         19, 'e', 12)
1199    .arg(_fitInterval, 19, 'e', 12)
1200    .arg("",           19, QChar(' '))
1201    .arg("",           19, QChar(' '));
1202
1203  return rnxStr;
1204}
1205
1206// RINEX Format String
1207//////////////////////////////////////////////////////////////////////////////
1208QString t_ephGlo::toString(double version) const {
1209
1210  QString rnxStr = rinexDateStr(_TOC-_gps_utc, _prn, version);
1211
1212  QTextStream out(&rnxStr);
1213
1214  out << QString("%1%2%3\n")
1215    .arg(-_tau,  19, 'e', 12)
1216    .arg(_gamma, 19, 'e', 12)
1217    .arg(_tki,   19, 'e', 12);
1218
1219  QString fmt = version < 3.0 ? "   %1%2%3%4\n" : "    %1%2%3%4\n";
1220
1221  out << QString(fmt)
1222    .arg(_x_pos,          19, 'e', 12)
1223    .arg(_x_velocity,     19, 'e', 12)
1224    .arg(_x_acceleration, 19, 'e', 12)
1225    .arg(_health,         19, 'e', 12);
1226
1227  out << QString(fmt)
1228    .arg(_y_pos,            19, 'e', 12)
1229    .arg(_y_velocity,       19, 'e', 12)
1230    .arg(_y_acceleration,   19, 'e', 12)
1231    .arg(_frequency_number, 19, 'e', 12);
1232
1233  out << QString(fmt)
1234    .arg(_z_pos,          19, 'e', 12)
1235    .arg(_z_velocity,     19, 'e', 12)
1236    .arg(_z_acceleration, 19, 'e', 12)
1237    .arg(_E,              19, 'e', 12);
1238
1239  return rnxStr;
1240}
1241
1242// RINEX Format String
1243//////////////////////////////////////////////////////////////////////////////
1244QString t_ephGal::toString(double version) const {
1245
1246  QString rnxStr = rinexDateStr(_TOC, _prn, version);
1247
1248  QTextStream out(&rnxStr);
1249
1250  out << QString("%1%2%3\n")
1251    .arg(_clock_bias,      19, 'e', 12)
1252    .arg(_clock_drift,     19, 'e', 12)
1253    .arg(_clock_driftrate, 19, 'e', 12);
1254
1255  QString fmt = version < 3.0 ? "   %1%2%3%4\n" : "    %1%2%3%4\n";
1256
1257  out << QString(fmt)
1258    .arg(_IODnav,  19, 'e', 12)
1259    .arg(_Crs,     19, 'e', 12)
1260    .arg(_Delta_n, 19, 'e', 12)
1261    .arg(_M0,      19, 'e', 12);
1262
1263  out << QString(fmt)
1264    .arg(_Cuc,    19, 'e', 12)
1265    .arg(_e,      19, 'e', 12)
1266    .arg(_Cus,    19, 'e', 12)
1267    .arg(_sqrt_A, 19, 'e', 12);
1268
1269  out << QString(fmt)
1270    .arg(_TOEsec, 19, 'e', 12)
1271    .arg(_Cic,    19, 'e', 12)
1272    .arg(_OMEGA0, 19, 'e', 12)
1273    .arg(_Cis,    19, 'e', 12);
1274
1275  out << QString(fmt)
1276    .arg(_i0,       19, 'e', 12)
1277    .arg(_Crc,      19, 'e', 12)
1278    .arg(_omega,    19, 'e', 12)
1279    .arg(_OMEGADOT, 19, 'e', 12);
1280
1281  out << QString(fmt)
1282    .arg(_IDOT,    19, 'e', 12)
1283    .arg(0.0,      19, 'e', 12)
1284    .arg(_TOEweek, 19, 'e', 12)
1285    .arg(0.0,      19, 'e', 12);
1286
1287  out << QString(fmt)
1288    .arg(_SISA,     19, 'e', 12)
1289    .arg(_E5aHS,    19, 'e', 12)
1290    .arg(_BGD_1_5A, 19, 'e', 12)
1291    .arg(_BGD_1_5B, 19, 'e', 12);
1292
1293  out << QString(fmt)
1294    .arg(_TOT,    19, 'e', 12)
1295    .arg("",      19, QChar(' '))
1296    .arg("",      19, QChar(' '))
1297    .arg("",      19, QChar(' '));
1298
1299  return rnxStr;
1300}
1301
Note: See TracBrowser for help on using the repository browser.