Changeset 3182 in ntrip


Ignore:
Timestamp:
Mar 29, 2011, 7:15:07 PM (14 years ago)
Author:
mervart
Message:
 
Location:
trunk/BNC
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • trunk/BNC/bncephuser.h

    r3174 r3182  
    7878  ~bncEphUser();
    7979
    80  public slots:
    81   void slotNewEphGPS(gpsephemeris gpseph);
    82   void slotNewEphGlonass(glonassephemeris gloeph);
    83   void slotNewEphGalileo(galileoephemeris galeph);
    84 
    85  protected:
    86 
    8780  class t_ephPair {
    8881   public:
    … …  
    9588      delete prev;
    9689    }
    97     ColumnVector xx;
    9890    t_eph* last;
    9991    t_eph* prev;
    10092  };
    10193
     94
     95 public slots:
     96  void slotNewEphGPS(gpsephemeris gpseph);
     97  void slotNewEphGlonass(glonassephemeris gloeph);
     98  void slotNewEphGalileo(galileoephemeris galeph);
     99
     100 protected:
    102101  QMutex                    _mutex;
    103102  QMap<QString, t_ephPair*> _eph;

  • trunk/BNC/upload/bncrtnetdecoder.cpp

    r3179 r3182  
    4343#include "bncutils.h"
    4444#include "bncsettings.h"
    45 #include "bncclockrinex.h"
    46 #include "bncsp3.h"
    4745
    4846using namespace std;
    … …  
    5149////////////////////////////////////////////////////////////////////////
    5250bncRtnetDecoder::bncRtnetDecoder() {
     51
    5352  bncSettings settings;
     53
    5454  _year = 0;
    55   _append = Qt::CheckState(settings.value("rnxAppend").toInt()) == Qt::Checked;
    5655
    57   // RINEX writer
    58   // ------------
    59   if ( settings.value("rnxPath").toString().isEmpty() ) {
    60     _rnx = 0;
    61   }
    62   else {
    63     QString prep  = "BNC";
    64     QString ext   = ".clk";
    65     QString path  = settings.value("rnxPath").toString();
    66     QString intr  = settings.value("rnxIntr").toString();
    67     int     sampl = settings.value("rnxSampl").toInt();
    68     _rnx = new bncClockRinex(prep, ext, path, intr, sampl);
    69   }
     56  // List of upload casters
     57  // ----------------------
    7058
    71   // SP3 writer
    72   // ----------
    73   if ( settings.value("sp3Path").toString().isEmpty() ) {
    74     _sp3 = 0;
    75   }
    76   else {
    77     QString prep  = "BNC";
    78     QString ext   = ".sp3";
    79     QString path  = settings.value("sp3Path").toString();
    80     QString intr  = settings.value("sp3Intr").toString();
    81     int     sampl = settings.value("sp3Sampl").toInt();
    82     _sp3 = new bncSP3(prep, ext, path, intr, sampl);
    83   }
    8459}
    8560
    … …  
    8762////////////////////////////////////////////////////////////////////////
    8863bncRtnetDecoder::~bncRtnetDecoder() {
    89   delete _rnx;
    90   delete _sp3;
    9164}
    9265
    … …  
    127100  // -----------------------------------
    128101  if (lines.size() > 0) {
    129 
    130     QStringList prns;
    131 
    132102    for (int ic = 0; ic < _caster.size(); ic++) {
    133       _caster.at(ic)->open();
    134 
    135       struct ClockOrbit co;
    136       memset(&co, 0, sizeof(co));
    137       co.GPSEpochTime      = (int)_GPSweeks;
    138       co.GLONASSEpochTime  = (int)fmod(_GPSweeks, 86400.0)
    139                            + 3 * 3600 - gnumleap(_year, _month, _day);
    140       co.ClockDataSupplied = 1;
    141       co.OrbitDataSupplied = 1;
    142       co.SatRefDatum       = DATUM_ITRF;
    143      
    144       struct Bias bias;
    145       memset(&bias, 0, sizeof(bias));
    146       bias.GPSEpochTime      = (int)_GPSweeks;
    147       bias.GLONASSEpochTime  = (int)fmod(_GPSweeks, 86400.0)
    148                              + 3 * 3600 - gnumleap(_year, _month, _day);
    149 
    150       for (int ii = 0; ii < lines.size(); ii++) {
    151 
    152         QString      prn;
    153         ColumnVector xx(14); xx = 0.0;
    154         t_ephPair*   pair = 0;
    155      
    156         if (ic == 0) {
    157           QTextStream in(lines[ii].toAscii());
    158           in >> prn;
    159           prns << prn;
    160           if ( _eph.contains(prn) ) {
    161             in >> xx(1) >> xx(2) >> xx(3) >> xx(4) >> xx(5)
    162                >> xx(6) >> xx(7) >> xx(8) >> xx(9) >> xx(10)
    163                >> xx(11) >> xx(12) >> xx(13) >> xx(14);
    164             xx(1) *= 1e3;          // x-crd
    165             xx(2) *= 1e3;          // y-crd
    166             xx(3) *= 1e3;          // z-crd
    167             xx(4) *= 1e-6;         // clk
    168             xx(5) *= 1e-6;         // rel. corr.
    169                                    // xx(6), xx(7), xx(8) ... PhaseCent - CoM
    170                                    // xx(9)  ... P1-C1 DCB in meters
    171                                    // xx(10) ... P1-P2 DCB in meters
    172                                    // xx(11) ... dT
    173             xx(12) *= 1e3;         // x-crd at time + dT
    174             xx(13) *= 1e3;         // y-crd at time + dT
    175             xx(14) *= 1e3;         // z-crd at time + dT
    176 
    177             pair     = _eph[prn];
    178             pair->xx = xx;
    179           }
    180         }
    181         else {
    182           prn = prns[ii];
    183           if ( _eph.contains(prn) ) {
    184             pair = _eph[prn];
    185             xx   = pair->xx;
    186           }
    187         }
    188 
    189         // Use old ephemeris if the new one is too recent
    190         // ----------------------------------------------
    191         t_eph* ep = 0;
    192         if (pair) {
    193           ep = pair->last;
    194           if (pair->prev && ep &&
    195               ep->receptDateTime().secsTo(QDateTime::currentDateTime()) < 60) {
    196             ep = pair->prev;
    197           }
    198         }
    199 
    200         if (ep != 0) {
    201           struct ClockOrbit::SatData* sd = 0;
    202           if      (prn[0] == 'G') {
    203             sd = co.Sat + co.NumberOfGPSSat;
    204             ++co.NumberOfGPSSat;
    205           }
    206           else if (prn[0] == 'R') {
    207             sd = co.Sat + CLOCKORBIT_NUMGPS + co.NumberOfGLONASSSat;
    208             ++co.NumberOfGLONASSSat;
    209           }
    210           if (sd) {
    211             QString outLine;
    212             processSatellite(ic, _caster.at(ic)->crdTrafo(),
    213                              _caster.at(ic)->CoM(), ep,
    214                              _GPSweek, _GPSweeks, prn, xx, sd, outLine);
    215             _caster.at(ic)->printAscii(outLine);
    216           }
    217 
    218           struct Bias::BiasSat* biasSat = 0;
    219           if      (prn[0] == 'G') {
    220             biasSat = bias.Sat + bias.NumberOfGPSSat;
    221             ++bias.NumberOfGPSSat;
    222           }
    223           else if (prn[0] == 'R') {
    224             biasSat = bias.Sat + CLOCKORBIT_NUMGPS + bias.NumberOfGLONASSSat;
    225             ++bias.NumberOfGLONASSSat;
    226           }
    227 
    228           // Coefficient of Ionosphere-Free LC
    229           // ---------------------------------
    230           const static double a_L1_GPS =  2.54572778;
    231           const static double a_L2_GPS = -1.54572778;
    232           const static double a_L1_Glo =  2.53125000;
    233           const static double a_L2_Glo = -1.53125000;
    234 
    235           if (biasSat) {
    236             biasSat->ID = prn.mid(1).toInt();
    237             biasSat->NumberOfCodeBiases = 3;
    238             if      (prn[0] == 'G') {
    239               biasSat->Biases[0].Type = CODETYPEGPS_L1_Z;
    240               biasSat->Biases[0].Bias = - a_L2_GPS * xx(10);
    241               biasSat->Biases[1].Type = CODETYPEGPS_L1_CA;
    242               biasSat->Biases[1].Bias = - a_L2_GPS * xx(10) + xx(9);
    243               biasSat->Biases[2].Type = CODETYPEGPS_L2_Z;
    244               biasSat->Biases[2].Bias = a_L1_GPS * xx(10);
    245             }
    246             else if (prn[0] == 'R') {
    247               biasSat->Biases[0].Type = CODETYPEGLONASS_L1_P;
    248               biasSat->Biases[0].Bias = - a_L2_Glo * xx(10);
    249               biasSat->Biases[1].Type = CODETYPEGLONASS_L1_CA;
    250               biasSat->Biases[1].Bias = - a_L2_Glo * xx(10) + xx(9);
    251               biasSat->Biases[2].Type = CODETYPEGLONASS_L2_P;
    252               biasSat->Biases[2].Bias = a_L1_Glo * xx(10);
    253             }
    254           }
    255         }
    256       }
    257      
    258       if ( _caster.at(ic)->usedSocket() &&
    259            (co.NumberOfGPSSat > 0 || co.NumberOfGLONASSSat > 0) ) {
    260         char obuffer[CLOCKORBIT_BUFFERSIZE];
    261 
    262         int len = MakeClockOrbit(&co, COTYPE_AUTO, 0, obuffer, sizeof(obuffer));
    263         if (len > 0) {
    264           _caster.at(ic)->write(obuffer, len);
    265         }
    266       }
    267 
    268       if ( _caster.at(ic)->usedSocket() &&
    269            (bias.NumberOfGPSSat > 0 || bias.NumberOfGLONASSSat > 0) ) {
    270         char obuffer[CLOCKORBIT_BUFFERSIZE];
    271         int len = MakeBias(&bias, BTYPE_AUTO, 0, obuffer, sizeof(obuffer));
    272         if (len > 0) {
    273           _caster.at(ic)->write(obuffer, len);
    274         }
    275       }
     103      _caster.at(ic)->uploadClockOrbitBias(lines, _eph, _year, _month, _day,
     104                                           _GPSweek, _GPSweeks);
    276105    }
    277106  }
    … …  
    280109}
    281110
    282 //
    283 ////////////////////////////////////////////////////////////////////////////
    284 void bncRtnetDecoder::processSatellite(int iCaster, const QString trafo,
    285                                        bool CoM, t_eph* ep, int GPSweek,
    286                                        double GPSweeks, const QString& prn,
    287                                        const ColumnVector& xx,
    288                                        struct ClockOrbit::SatData* sd,
    289                                        QString& outLine) {
    290 
    291   const double secPerWeek = 7.0 * 86400.0;
    292 
    293   ColumnVector rsw(3);
    294   ColumnVector rsw2(3);
    295   double dClk;
    296 
    297   for (int ii = 1; ii <= 2; ++ii) {
    298 
    299     int    GPSweek12  = GPSweek;
    300     double GPSweeks12 = GPSweeks;
    301     if (ii == 2) {
    302       GPSweeks12 += xx(11);
    303       if (GPSweeks12 > secPerWeek) {
    304         GPSweek12  += 1;
    305         GPSweeks12 -= secPerWeek;
    306       }
    307     }
    308 
    309     ColumnVector xB(4);
    310     ColumnVector vv(3);
    311 
    312     ep->position(GPSweek12, GPSweeks12, xB.data(), vv.data());
    313    
    314     ColumnVector xyz;
    315     if (ii == 1) {
    316       xyz = xx.Rows(1,3);
    317     }
    318     else {
    319       xyz = xx.Rows(12,14);
    320     }
    321    
    322     // Correction Center of Mass -> Antenna Phase Center
    323     // -------------------------------------------------
    324     if (! CoM) {
    325       xyz(1) += xx(6);
    326       xyz(2) += xx(7);
    327       xyz(3) += xx(8);
    328     }
    329    
    330     if (trafo != "IGS05") {
    331       crdTrafo(GPSweek12, xyz, trafo);
    332     }
    333    
    334     ColumnVector dx = xB.Rows(1,3) - xyz ;
    335    
    336     if (ii == 1) {
    337       XYZ_to_RSW(xB.Rows(1,3), vv, dx, rsw);
    338       dClk = (xx(4) + xx(5) - xB(4)) * 299792458.0;
    339     }
    340     else {
    341       XYZ_to_RSW(xB.Rows(1,3), vv, dx, rsw2);
    342     }
    343   }
    344 
    345   if (sd) {
    346     sd->ID                    = prn.mid(1).toInt();
    347     sd->IOD                   = ep->IOD();
    348     sd->Clock.DeltaA0         = dClk;
    349     sd->Orbit.DeltaRadial     = rsw(1);
    350     sd->Orbit.DeltaAlongTrack = rsw(2);
    351     sd->Orbit.DeltaCrossTrack = rsw(3);
    352     sd->Orbit.DotDeltaRadial     = (rsw2(1) - rsw(1)) / xx(11);
    353     sd->Orbit.DotDeltaAlongTrack = (rsw2(2) - rsw(2)) / xx(11);
    354     sd->Orbit.DotDeltaCrossTrack = (rsw2(3) - rsw(3)) / xx(11);
    355   }
    356 
    357   outLine.sprintf("%d %.1f %s  %3d  %10.3f  %8.3f %8.3f %8.3f\n",
    358                   GPSweek, GPSweeks, ep->prn().c_str(),
    359                   ep->IOD(), dClk, rsw(1), rsw(2), rsw(3));
    360 
    361   if (iCaster == 0) {
    362     if (_rnx) {
    363       _rnx->write(GPSweek, GPSweeks, prn, xx);
    364     }
    365     if (_sp3) {
    366       _sp3->write(GPSweek, GPSweeks, prn, xx, _append);
    367     }
    368   }
    369 }
    370 
    371 // Transform Coordinates
    372 ////////////////////////////////////////////////////////////////////////////
    373 void bncRtnetDecoder::crdTrafo(int GPSWeek, ColumnVector& xyz,
    374                                const QString& trafo) {
    375 
    376   bncSettings settings;
    377 
    378   if (trafo == "ETRF2000") {
    379     _dx  =    0.0541;
    380     _dy  =    0.0502;
    381     _dz  =   -0.0538;
    382     _dxr =   -0.0002;
    383     _dyr =    0.0001;
    384     _dzr =   -0.0018;
    385     _ox  =  0.000891;
    386     _oy  =  0.005390;
    387     _oz  = -0.008712;
    388     _oxr =  0.000081;
    389     _oyr =  0.000490;
    390     _ozr = -0.000792;
    391     _sc  =      0.40;
    392     _scr =      0.08;
    393     _t0  =    2000.0;
    394   }
    395   else if (trafo == "NAD83") {
    396     _dx  =    0.9963;
    397     _dy  =   -1.9024;
    398     _dz  =   -0.5210;
    399     _dxr =    0.0005;
    400     _dyr =   -0.0006;
    401     _dzr =   -0.0013;
    402     _ox  =  0.025915;
    403     _oy  =  0.009426;
    404     _oz  =  0.011599;
    405     _oxr =  0.000067;
    406     _oyr = -0.000757;
    407     _ozr = -0.000051;
    408     _sc  =      0.78;
    409     _scr =     -0.10;
    410     _t0  =    1997.0;
    411   }
    412   else if (trafo == "GDA94") {
    413     _dx  =   -0.07973;
    414     _dy  =   -0.00686;
    415     _dz  =    0.03803;
    416     _dxr =    0.00225;
    417     _dyr =   -0.00062;
    418     _dzr =   -0.00056;
    419     _ox  =  0.0000351;
    420     _oy  = -0.0021211;
    421     _oz  = -0.0021411;
    422     _oxr = -0.0014707;
    423     _oyr = -0.0011443;
    424     _ozr = -0.0011701;
    425     _sc  =      6.636;
    426     _scr =      0.294;
    427     _t0  =     1994.0;
    428   }
    429   else if (trafo == "SIRGAS2000") {
    430     _dx  =   -0.0051;
    431     _dy  =   -0.0065;
    432     _dz  =   -0.0099;
    433     _dxr =    0.0000;
    434     _dyr =    0.0000;
    435     _dzr =    0.0000;
    436     _ox  =  0.000150;
    437     _oy  =  0.000020;
    438     _oz  =  0.000021;
    439     _oxr =  0.000000;
    440     _oyr =  0.000000;
    441     _ozr =  0.000000;
    442     _sc  =     0.000;
    443     _scr =     0.000;
    444     _t0  =    0000.0;
    445   }
    446   else if (trafo == "SIRGAS95") {
    447     _dx  =    0.0077;
    448     _dy  =    0.0058;
    449     _dz  =   -0.0138;
    450     _dxr =    0.0000;
    451     _dyr =    0.0000;
    452     _dzr =    0.0000;
    453     _ox  =  0.000000;
    454     _oy  =  0.000000;
    455     _oz  = -0.000030;
    456     _oxr =  0.000000;
    457     _oyr =  0.000000;
    458     _ozr =  0.000000;
    459     _sc  =     1.570;
    460     _scr =     0.000;
    461     _t0  =    0000.0;
    462   }
    463   else if (trafo == "Custom") {
    464     _dx  = settings.value("trafo_dx").toDouble();
    465     _dy  = settings.value("trafo_dy").toDouble();
    466     _dz  = settings.value("trafo_dz").toDouble();
    467     _dxr = settings.value("trafo_dxr").toDouble();
    468     _dyr = settings.value("trafo_dyr").toDouble();
    469     _dzr = settings.value("trafo_dzr").toDouble();
    470     _ox  = settings.value("trafo_ox").toDouble();
    471     _oy  = settings.value("trafo_oy").toDouble();
    472     _oz  = settings.value("trafo_oz").toDouble();
    473     _oxr = settings.value("trafo_oxr").toDouble();
    474     _oyr = settings.value("trafo_oyr").toDouble();
    475     _ozr = settings.value("trafo_ozr").toDouble();
    476     _sc  = settings.value("trafo_sc").toDouble();
    477     _scr = settings.value("trafo_scr").toDouble();
    478     _t0  = settings.value("trafo_t0").toDouble();
    479   }
    480 
    481   // Current epoch minus 2000.0 in years
    482   // ------------------------------------
    483   double dt = (GPSWeek - (1042.0+6.0/7.0)) / 365.2422 * 7.0 + 2000.0 - _t0;
    484 
    485   ColumnVector dx(3);
    486 
    487   dx(1) = _dx + dt * _dxr;
    488   dx(2) = _dy + dt * _dyr;
    489   dx(3) = _dz + dt * _dzr;
    490 
    491   static const double arcSec = 180.0 * 3600.0 / M_PI;
    492 
    493   double ox = (_ox + dt * _oxr) / arcSec;
    494   double oy = (_oy + dt * _oyr) / arcSec;
    495   double oz = (_oz + dt * _ozr) / arcSec;
    496 
    497   double sc = 1.0 + _sc * 1e-9 + dt * _scr * 1e-9;
    498 
    499   Matrix rMat(3,3);
    500   rMat(1,1) = 1.0;
    501   rMat(1,2) = -oz;
    502   rMat(1,3) =  oy;
    503   rMat(2,1) =  oz;
    504   rMat(2,2) = 1.0;
    505   rMat(2,3) = -ox;
    506   rMat(3,1) = -oy;
    507   rMat(3,2) =  ox;
    508   rMat(3,3) = 1.0;
    509 
    510   xyz = sc * rMat * xyz + dx;
    511 }

  • trunk/BNC/upload/bncrtnetdecoder.h

    r3179 r3182  
    3232#include "RTCM/GPSDecoder.h"
    3333
    34 class bncClockRinex;
    35 class bncSP3;
    36 
    3734class bncRtnetDecoder: public GPSDecoder, public bncEphUser {
    3835 public:
    … …  
    4239 private:
    4340  void readEpochTime(const QString& line);
    44   void processSatellite(int iCaster, const QString trafo,
    45                         bool CoM, t_eph* ep, int GPSweek,
    46                         double GPSweeks, const QString& prn,
    47                         const ColumnVector& xx,
    48                         struct ClockOrbit::SatData* sd,
    49                         QString& outLine);
    50   void crdTrafo(int GPSWeek, ColumnVector& xyz,
    51                 const QString& trafo);
    52 
    5341  QList<bncUploadCaster*> _caster;
    54   bool                    _append;
    5542  QString                 _buffer;
    5643  int                     _GPSweek;
    … …  
    6249  int                     _min;
    6350  double                  _sec;
    64   double                  _dx;
    65   double                  _dy;
    66   double                  _dz;
    67   double                  _dxr;
    68   double                  _dyr;
    69   double                  _dzr;
    70   double                  _ox;
    71   double                  _oy;
    72   double                  _oz;
    73   double                  _oxr;
    74   double                  _oyr;
    75   double                  _ozr;
    76   double                  _sc;
    77   double                  _scr;
    78   double                  _t0;
    79   bncClockRinex*          _rnx;
    80   bncSP3*                 _sp3;
    8151};
    8252

  • trunk/BNC/upload/bncuploadcaster.cpp

    r3174 r3182  
    2020#include "bncversion.h"
    2121#include "bncapp.h"
     22#include "bncclockrinex.h"
     23#include "bncsp3.h"
    2224
    2325using namespace std;
    … …  
    2931                                 const QString& password,
    3032                                 const QString& crdTrafo, bool  CoM,
     33                                 const QString& rnxFileName,
     34                                 const QString& sp3FileName,
    3135                                 const QString& outFileName) {
    3236
    3337  bncSettings settings;
     38
     39  connect(this, SIGNAL(newMessage(QByteArray,bool)),
     40          ((bncApp*)qApp), SLOT(slotMessage(const QByteArray,bool)));
    3441
    3542  _mountpoint = mountpoint;
    … …  
    4350  _sOpenTrial = 0;
    4451
     52  _append = Qt::CheckState(settings.value("rnxAppend").toInt()) == Qt::Checked;
     53
    4554  if (outFileName.isEmpty()) {
    4655    _outFile   = 0;
    … …  
    5059    _outFile = new QFile(outFileName);
    5160    QIODevice::OpenMode oMode;
    52     if (Qt::CheckState(settings.value("rnxAppend").toInt()) == Qt::Checked) {
     61    if (_append) {
    5362      oMode = QIODevice::WriteOnly | QIODevice::Unbuffered | QIODevice::Append;
    5463    }
    … …  
    6170    }
    6271  }
    63   connect(this, SIGNAL(newMessage(QByteArray,bool)),
    64           ((bncApp*)qApp), SLOT(slotMessage(const QByteArray,bool)));
     72
     73  // RINEX writer
     74  // ------------
     75  if ( settings.value("rnxPath").toString().isEmpty() ) {
     76    _rnx = 0;
     77  }
     78  else {
     79    QString prep  = "BNC";
     80    QString ext   = ".clk";
     81    QString path  = settings.value("rnxPath").toString();
     82    QString intr  = settings.value("rnxIntr").toString();
     83    int     sampl = settings.value("rnxSampl").toInt();
     84    _rnx = new bncClockRinex(prep, ext, path, intr, sampl);
     85  }
     86
     87  // SP3 writer
     88  // ----------
     89  if ( settings.value("sp3Path").toString().isEmpty() ) {
     90    _sp3 = 0;
     91  }
     92  else {
     93    QString prep  = "BNC";
     94    QString ext   = ".sp3";
     95    QString path  = settings.value("sp3Path").toString();
     96    QString intr  = settings.value("sp3Intr").toString();
     97    int     sampl = settings.value("sp3Sampl").toInt();
     98    _sp3 = new bncSP3(prep, ext, path, intr, sampl);
     99  }
     100
     101  // Set Transformation Parameters
     102  // -----------------------------
     103  if      (_crdTrafo == "ETRF2000") {
     104    _dx  =    0.0541;
     105    _dy  =    0.0502;
     106    _dz  =   -0.0538;
     107    _dxr =   -0.0002;
     108    _dyr =    0.0001;
     109    _dzr =   -0.0018;
     110    _ox  =  0.000891;
     111    _oy  =  0.005390;
     112    _oz  = -0.008712;
     113    _oxr =  0.000081;
     114    _oyr =  0.000490;
     115    _ozr = -0.000792;
     116    _sc  =      0.40;
     117    _scr =      0.08;
     118    _t0  =    2000.0;
     119  }
     120  else if (_crdTrafo == "NAD83") {
     121    _dx  =    0.9963;
     122    _dy  =   -1.9024;
     123    _dz  =   -0.5210;
     124    _dxr =    0.0005;
     125    _dyr =   -0.0006;
     126    _dzr =   -0.0013;
     127    _ox  =  0.025915;
     128    _oy  =  0.009426;
     129    _oz  =  0.011599;
     130    _oxr =  0.000067;
     131    _oyr = -0.000757;
     132    _ozr = -0.000051;
     133    _sc  =      0.78;
     134    _scr =     -0.10;
     135    _t0  =    1997.0;
     136  }
     137  else if (_crdTrafo == "GDA94") {
     138    _dx  =   -0.07973;
     139    _dy  =   -0.00686;
     140    _dz  =    0.03803;
     141    _dxr =    0.00225;
     142    _dyr =   -0.00062;
     143    _dzr =   -0.00056;
     144    _ox  =  0.0000351;
     145    _oy  = -0.0021211;
     146    _oz  = -0.0021411;
     147    _oxr = -0.0014707;
     148    _oyr = -0.0011443;
     149    _ozr = -0.0011701;
     150    _sc  =      6.636;
     151    _scr =      0.294;
     152    _t0  =     1994.0;
     153  }
     154  else if (_crdTrafo == "SIRGAS2000") {
     155    _dx  =   -0.0051;
     156    _dy  =   -0.0065;
     157    _dz  =   -0.0099;
     158    _dxr =    0.0000;
     159    _dyr =    0.0000;
     160    _dzr =    0.0000;
     161    _ox  =  0.000150;
     162    _oy  =  0.000020;
     163    _oz  =  0.000021;
     164    _oxr =  0.000000;
     165    _oyr =  0.000000;
     166    _ozr =  0.000000;
     167    _sc  =     0.000;
     168    _scr =     0.000;
     169    _t0  =    0000.0;
     170  }
     171  else if (_crdTrafo == "SIRGAS95") {
     172    _dx  =    0.0077;
     173    _dy  =    0.0058;
     174    _dz  =   -0.0138;
     175    _dxr =    0.0000;
     176    _dyr =    0.0000;
     177    _dzr =    0.0000;
     178    _ox  =  0.000000;
     179    _oy  =  0.000000;
     180    _oz  = -0.000030;
     181    _oxr =  0.000000;
     182    _oyr =  0.000000;
     183    _ozr =  0.000000;
     184    _sc  =     1.570;
     185    _scr =     0.000;
     186    _t0  =    0000.0;
     187  }
     188  else if (_crdTrafo == "Custom") {
     189    _dx  = settings.value("trafo_dx").toDouble();
     190    _dy  = settings.value("trafo_dy").toDouble();
     191    _dz  = settings.value("trafo_dz").toDouble();
     192    _dxr = settings.value("trafo_dxr").toDouble();
     193    _dyr = settings.value("trafo_dyr").toDouble();
     194    _dzr = settings.value("trafo_dzr").toDouble();
     195    _ox  = settings.value("trafo_ox").toDouble();
     196    _oy  = settings.value("trafo_oy").toDouble();
     197    _oz  = settings.value("trafo_oz").toDouble();
     198    _oxr = settings.value("trafo_oxr").toDouble();
     199    _oyr = settings.value("trafo_oyr").toDouble();
     200    _ozr = settings.value("trafo_ozr").toDouble();
     201    _sc  = settings.value("trafo_sc").toDouble();
     202    _scr = settings.value("trafo_scr").toDouble();
     203    _t0  = settings.value("trafo_t0").toDouble();
     204  }
    65205}
    66206
    … …  
    150290  }
    151291}
     292
     293// Encode and Upload Clocks, Orbits, and Biases
     294////////////////////////////////////////////////////////////////////////////
     295void bncUploadCaster::uploadClockOrbitBias(const QStringList& lines,
     296                           const QMap<QString, bncEphUser::t_ephPair*>& ephMap,
     297                           int year, int month, int day,
     298                           int GPSweek, double GPSweeks) {
     299
     300  this->open();
     301 
     302  struct ClockOrbit co;
     303  memset(&co, 0, sizeof(co));
     304  co.GPSEpochTime      = (int)GPSweeks;
     305  co.GLONASSEpochTime  = (int)fmod(GPSweeks, 86400.0)
     306                       + 3 * 3600 - gnumleap(year, month, day);
     307  co.ClockDataSupplied = 1;
     308  co.OrbitDataSupplied = 1;
     309  co.SatRefDatum       = DATUM_ITRF;
     310 
     311  struct Bias bias;
     312  memset(&bias, 0, sizeof(bias));
     313  bias.GPSEpochTime      = (int)GPSweeks;
     314  bias.GLONASSEpochTime  = (int)fmod(GPSweeks, 86400.0)
     315                         + 3 * 3600 - gnumleap(year, month, day);
     316 
     317  for (int ii = 0; ii < lines.size(); ii++) {
     318 
     319    QString      prn;
     320    ColumnVector xx(14); xx = 0.0;
     321    bncEphUser::t_ephPair*   pair = 0;
     322 
     323    QTextStream in(lines[ii].toAscii());
     324    in >> prn;
     325    if ( ephMap.contains(prn) ) {
     326      in >> xx(1) >> xx(2) >> xx(3) >> xx(4) >> xx(5)
     327         >> xx(6) >> xx(7) >> xx(8) >> xx(9) >> xx(10)
     328         >> xx(11) >> xx(12) >> xx(13) >> xx(14);
     329      xx(1) *= 1e3;          // x-crd
     330      xx(2) *= 1e3;          // y-crd
     331      xx(3) *= 1e3;          // z-crd
     332      xx(4) *= 1e-6;         // clk
     333      xx(5) *= 1e-6;         // rel. corr.
     334                             // xx(6), xx(7), xx(8) ... PhaseCent - CoM
     335                             // xx(9)  ... P1-C1 DCB in meters
     336                             // xx(10) ... P1-P2 DCB in meters
     337                             // xx(11) ... dT
     338      xx(12) *= 1e3;         // x-crd at time + dT
     339      xx(13) *= 1e3;         // y-crd at time + dT
     340      xx(14) *= 1e3;         // z-crd at time + dT
     341 
     342      pair = ephMap[prn];
     343    }
     344 
     345    // Use old ephemeris if the new one is too recent
     346    // ----------------------------------------------
     347    t_eph* ep = 0;
     348    if (pair) {
     349      ep = pair->last;
     350      if (pair->prev && ep &&
     351          ep->receptDateTime().secsTo(QDateTime::currentDateTime()) < 60) {
     352        ep = pair->prev;
     353      }
     354    }
     355 
     356    if (ep != 0) {
     357      struct ClockOrbit::SatData* sd = 0;
     358      if      (prn[0] == 'G') {
     359        sd = co.Sat + co.NumberOfGPSSat;
     360        ++co.NumberOfGPSSat;
     361      }
     362      else if (prn[0] == 'R') {
     363        sd = co.Sat + CLOCKORBIT_NUMGPS + co.NumberOfGLONASSSat;
     364        ++co.NumberOfGLONASSSat;
     365      }
     366      if (sd) {
     367        QString outLine;
     368        processSatellite(ep, GPSweek, GPSweeks, prn, xx, sd, outLine);
     369        this->printAscii(outLine);
     370      }
     371 
     372      struct Bias::BiasSat* biasSat = 0;
     373      if      (prn[0] == 'G') {
     374        biasSat = bias.Sat + bias.NumberOfGPSSat;
     375        ++bias.NumberOfGPSSat;
     376      }
     377      else if (prn[0] == 'R') {
     378        biasSat = bias.Sat + CLOCKORBIT_NUMGPS + bias.NumberOfGLONASSSat;
     379        ++bias.NumberOfGLONASSSat;
     380      }
     381 
     382      // Coefficient of Ionosphere-Free LC
     383      // ---------------------------------
     384      const static double a_L1_GPS =  2.54572778;
     385      const static double a_L2_GPS = -1.54572778;
     386      const static double a_L1_Glo =  2.53125000;
     387      const static double a_L2_Glo = -1.53125000;
     388 
     389      if (biasSat) {
     390        biasSat->ID = prn.mid(1).toInt();
     391        biasSat->NumberOfCodeBiases = 3;
     392        if      (prn[0] == 'G') {
     393          biasSat->Biases[0].Type = CODETYPEGPS_L1_Z;
     394          biasSat->Biases[0].Bias = - a_L2_GPS * xx(10);
     395          biasSat->Biases[1].Type = CODETYPEGPS_L1_CA;
     396          biasSat->Biases[1].Bias = - a_L2_GPS * xx(10) + xx(9);
     397          biasSat->Biases[2].Type = CODETYPEGPS_L2_Z;
     398          biasSat->Biases[2].Bias = a_L1_GPS * xx(10);
     399        }
     400        else if (prn[0] == 'R') {
     401          biasSat->Biases[0].Type = CODETYPEGLONASS_L1_P;
     402          biasSat->Biases[0].Bias = - a_L2_Glo * xx(10);
     403          biasSat->Biases[1].Type = CODETYPEGLONASS_L1_CA;
     404          biasSat->Biases[1].Bias = - a_L2_Glo * xx(10) + xx(9);
     405          biasSat->Biases[2].Type = CODETYPEGLONASS_L2_P;
     406          biasSat->Biases[2].Bias = a_L1_Glo * xx(10);
     407        }
     408      }
     409    }
     410  }
     411 
     412  if ( this->usedSocket() &&
     413       (co.NumberOfGPSSat > 0 || co.NumberOfGLONASSSat > 0) ) {
     414    char obuffer[CLOCKORBIT_BUFFERSIZE];
     415 
     416    int len = MakeClockOrbit(&co, COTYPE_AUTO, 0, obuffer, sizeof(obuffer));
     417    if (len > 0) {
     418      this->write(obuffer, len);
     419    }
     420  }
     421 
     422  if ( this->usedSocket() &&
     423       (bias.NumberOfGPSSat > 0 || bias.NumberOfGLONASSSat > 0) ) {
     424    char obuffer[CLOCKORBIT_BUFFERSIZE];
     425    int len = MakeBias(&bias, BTYPE_AUTO, 0, obuffer, sizeof(obuffer));
     426    if (len > 0) {
     427      this->write(obuffer, len);
     428    }
     429  }
     430}
     431
     432//
     433////////////////////////////////////////////////////////////////////////////
     434void bncUploadCaster::processSatellite(t_eph* eph, int GPSweek,
     435                                       double GPSweeks, const QString& prn,
     436                                       const ColumnVector& xx,
     437                                       struct ClockOrbit::SatData* sd,
     438                                       QString& outLine) {
     439
     440  const double secPerWeek = 7.0 * 86400.0;
     441
     442  ColumnVector rsw(3);
     443  ColumnVector rsw2(3);
     444  double dClk;
     445
     446  for (int ii = 1; ii <= 2; ++ii) {
     447
     448    int    GPSweek12  = GPSweek;
     449    double GPSweeks12 = GPSweeks;
     450    if (ii == 2) {
     451      GPSweeks12 += xx(11);
     452      if (GPSweeks12 > secPerWeek) {
     453        GPSweek12  += 1;
     454        GPSweeks12 -= secPerWeek;
     455      }
     456    }
     457
     458    ColumnVector xB(4);
     459    ColumnVector vv(3);
     460
     461    eph->position(GPSweek12, GPSweeks12, xB.data(), vv.data());
     462   
     463    ColumnVector xyz;
     464    if (ii == 1) {
     465      xyz = xx.Rows(1,3);
     466    }
     467    else {
     468      xyz = xx.Rows(12,14);
     469    }
     470   
     471    // Correction Center of Mass -> Antenna Phase Center
     472    // -------------------------------------------------
     473    if (! _CoM) {
     474      xyz(1) += xx(6);
     475      xyz(2) += xx(7);
     476      xyz(3) += xx(8);
     477    }
     478   
     479    if (_crdTrafo != "IGS05") {
     480      crdTrafo(GPSweek12, xyz);
     481    }
     482   
     483    ColumnVector dx = xB.Rows(1,3) - xyz ;
     484   
     485    if (ii == 1) {
     486      XYZ_to_RSW(xB.Rows(1,3), vv, dx, rsw);
     487      dClk = (xx(4) + xx(5) - xB(4)) * 299792458.0;
     488    }
     489    else {
     490      XYZ_to_RSW(xB.Rows(1,3), vv, dx, rsw2);
     491    }
     492  }
     493
     494  if (sd) {
     495    sd->ID                    = prn.mid(1).toInt();
     496    sd->IOD                   = eph->IOD();
     497    sd->Clock.DeltaA0         = dClk;
     498    sd->Orbit.DeltaRadial     = rsw(1);
     499    sd->Orbit.DeltaAlongTrack = rsw(2);
     500    sd->Orbit.DeltaCrossTrack = rsw(3);
     501    sd->Orbit.DotDeltaRadial     = (rsw2(1) - rsw(1)) / xx(11);
     502    sd->Orbit.DotDeltaAlongTrack = (rsw2(2) - rsw(2)) / xx(11);
     503    sd->Orbit.DotDeltaCrossTrack = (rsw2(3) - rsw(3)) / xx(11);
     504  }
     505
     506  outLine.sprintf("%d %.1f %s  %3d  %10.3f  %8.3f %8.3f %8.3f\n",
     507                  GPSweek, GPSweeks, eph->prn().c_str(),
     508                  eph->IOD(), dClk, rsw(1), rsw(2), rsw(3));
     509
     510  if (_rnx) {
     511    _rnx->write(GPSweek, GPSweeks, prn, xx);
     512  }
     513  if (_sp3) {
     514    _sp3->write(GPSweek, GPSweeks, prn, xx, _append);
     515  }
     516}
     517
     518// Transform Coordinates
     519////////////////////////////////////////////////////////////////////////////
     520void bncUploadCaster::crdTrafo(int GPSWeek, ColumnVector& xyz) {
     521
     522  // Current epoch minus 2000.0 in years
     523  // ------------------------------------
     524  double dt = (GPSWeek - (1042.0+6.0/7.0)) / 365.2422 * 7.0 + 2000.0 - _t0;
     525
     526  ColumnVector dx(3);
     527
     528  dx(1) = _dx + dt * _dxr;
     529  dx(2) = _dy + dt * _dyr;
     530  dx(3) = _dz + dt * _dzr;
     531
     532  static const double arcSec = 180.0 * 3600.0 / M_PI;
     533
     534  double ox = (_ox + dt * _oxr) / arcSec;
     535  double oy = (_oy + dt * _oyr) / arcSec;
     536  double oz = (_oz + dt * _ozr) / arcSec;
     537
     538  double sc = 1.0 + _sc * 1e-9 + dt * _scr * 1e-9;
     539
     540  Matrix rMat(3,3);
     541  rMat(1,1) = 1.0;
     542  rMat(1,2) = -oz;
     543  rMat(1,3) =  oy;
     544  rMat(2,1) =  oz;
     545  rMat(2,2) = 1.0;
     546  rMat(2,3) = -ox;
     547  rMat(3,1) = -oy;
     548  rMat(3,2) =  ox;
     549  rMat(3,3) = 1.0;
     550
     551  xyz = sc * rMat * xyz + dx;
     552}

  • trunk/BNC/upload/bncuploadcaster.h

    r3179 r3182  
    33
    44#include <QtNetwork>
     5#include "bncephuser.h"
     6
     7class bncClockRinex;
     8class bncSP3;
    59
    610class bncUploadCaster : public QObject {
    … …  
    1115                  const QString& password,
    1216                  const QString& crdTrafo, bool  CoM,
     17                  const QString& rnxFileName,
     18                  const QString& sp3FileName,
    1319                  const QString& outFileName);
    1420  virtual ~bncUploadCaster();
    … …  
    1723  void printAscii(const QString& line);
    1824  bool usedSocket() const {return _outSocket;}
    19   QString crdTrafo() const {return _crdTrafo;}
    20   bool CoM() const {return _CoM;}
     25  void uploadClockOrbitBias(const QStringList& lines,
     26                            const QMap<QString, bncEphUser::t_ephPair*>& ephMap,
     27                            int year, int month, int day,
     28                            int GPSweek, double GPSweeks);
    2129
    2230 signals:
    … …  
    2533
    2634 private:
     35  void processSatellite(t_eph* eph, int GPSweek,
     36                        double GPSweeks, const QString& prn,
     37                        const ColumnVector& xx,
     38                        struct ClockOrbit::SatData* sd,
     39                        QString& outLine);
     40  void crdTrafo(int GPSWeek, ColumnVector& xyz);
    2741  QString        _mountpoint;
    2842  QString        _outHost;
    … …  
    3650  QFile*         _outFile;
    3751  QTextStream*   _outStream;
     52  bool           _append;
     53  double         _dx;
     54  double         _dy;
     55  double         _dz;
     56  double         _dxr;
     57  double         _dyr;
     58  double         _dzr;
     59  double         _ox;
     60  double         _oy;
     61  double         _oz;
     62  double         _oxr;
     63  double         _oyr;
     64  double         _ozr;
     65  double         _sc;
     66  double         _scr;
     67  double         _t0;
     68  bncClockRinex* _rnx;
     69  bncSP3*        _sp3;
    3870};
    3971
Note: See TracChangeset for help on using the changeset viewer.