Changeset 10791 in ntrip for trunk/BNC/src/PPP


Ignore:
Timestamp:
Dec 3, 2025, 5:37:16 PM (5 months ago)
Author:
mervart
Message:

BNC Multifrequency and PPPAR Client (initial version)

Location:
trunk/BNC/src/PPP
Files:
4 added
11 edited

Legend:

Unmodified
Added
Removed

  • trunk/BNC/src/PPP/pppClient.cpp

    r10765 r10791  
    117117//////////////////////////////////////////////////////////////////////////////
    118118void t_pppClient::putOrbCorrections(const vector<t_orbCorr*>& corr) {
     119  if (OPT->_logMode == t_pppOptions::normal && corr.size() > 0) {
     120    LOG << "orbCorrections " << string(corr[0]->_time) << ' ' << corr.size() << endl;
     121  }
    119122  for (unsigned ii = 0; ii < corr.size(); ii++) {
    120123    _ephPool->putOrbCorrection(new t_orbCorr(*corr[ii]));
    … …  
    125128//////////////////////////////////////////////////////////////////////////////
    126129void t_pppClient::putClkCorrections(const vector<t_clkCorr*>& corr) {
     130  if (OPT->_logMode == t_pppOptions::normal && corr.size() > 0) {
     131    LOG << "clkCorrections " << string(corr[0]->_time) << ' ' << corr.size() << endl;
     132  }
    127133  for (unsigned ii = 0; ii < corr.size(); ii++) {
    128134    _ephPool->putClkCorrection(new t_clkCorr(*corr[ii]));
    … …  
    133139//////////////////////////////////////////////////////////////////////////////
    134140void t_pppClient::putCodeBiases(const vector<t_satCodeBias*>& biases) {
     141  if (OPT->_logMode == t_pppOptions::normal && biases.size() > 0) {
     142    LOG << "codeBiases     " << string(biases[0]->_time) << ' ' << biases.size() << endl;
     143  }
     144  set<char> systems;
    135145  for (unsigned ii = 0; ii < biases.size(); ii++) {
    136     _obsPool->putCodeBias(new t_satCodeBias(*biases[ii]));
     146    t_satCodeBias* newBias = new t_satCodeBias(*biases[ii]);
     147    char sys = newBias->_prn.system();
     148    if (systems.find(sys) == systems.end()) {
     149      _obsPool->clearCodeBiases(sys);
     150      systems.insert(sys);
     151    }
     152    _obsPool->putCodeBias(newBias);
    137153  }
    138154}
    … …  
    141157//////////////////////////////////////////////////////////////////////////////
    142158void t_pppClient::putPhaseBiases(const vector<t_satPhaseBias*>& biases) {
     159  if (OPT->_logMode == t_pppOptions::normal && biases.size() > 0) {
     160    LOG << "phaseBiases    " << string(biases[0]->_time) << ' ' << biases.size() << endl;
     161  }
     162  set<char> systems;
    143163  for (unsigned ii = 0; ii < biases.size(); ii++) {
    144     _obsPool->putPhaseBias(new t_satPhaseBias(*biases[ii]));
     164    t_satPhaseBias* newBias = new t_satPhaseBias(*biases[ii]);
     165    char sys = newBias->_prn.system();
     166    if (systems.find(sys) == systems.end()) {
     167      _obsPool->clearPhaseBiases(sys);
     168      systems.insert(sys);
     169    }
     170    _obsPool->putPhaseBias(newBias);
    145171  }
    146172}
    … …  
    211237    }
    212238  }
    213 /*
    214   vector<t_pppSatObs*>::iterator it = obsVector.begin();
    215   while (it != obsVector.end()) {
    216     t_pppSatObs* satObs = *it;
    217     satObs->printObsMinusComputed();
    218     it++;
    219   }
    220 */
     239
     240  if (OPT->_logMode == t_pppOptions::all) {
     241    vector<t_pppSatObs*>::iterator it = obsVector.begin();
     242    while (it != obsVector.end()) {
     243      t_pppSatObs* satObs = *it;
     244      satObs->printObsMinusComputed();
     245      it++;
     246    }
     247  }
    221248  return pseudoObsIono;
    222249}
    223 
    224 //
    225 //////////////////////////////////////////////////////////////////////////////
    226 void t_pppClient::useObsWithCodeBiasesOnly(std::vector<t_pppSatObs*>& obsVector) {
    227 
    228   vector<t_pppSatObs*>::iterator it = obsVector.begin();
    229   while (it != obsVector.end()) {
    230     t_pppSatObs* pppSatObs = *it;
    231     bool codeBiasesAvailable = false;
    232     if (pppSatObs->getCodeBias(pppSatObs->fType1()) &&
    233         pppSatObs->getCodeBias(pppSatObs->fType2())) {
    234         codeBiasesAvailable = true;
    235     }
    236     if (codeBiasesAvailable) {
    237       ++it;
    238     }
    239     else {
    240       it = obsVector.erase(it);
    241       delete pppSatObs;
    242     }
    243   }
    244   return;
    245 }
    246 
    247250
    248251// Compute the Bancroft position, check for blunders
    … …  
    251254                                  vector<t_pppSatObs*>& obsVector,
    252255                                  ColumnVector& xyzc, bool print) {
    253   t_lc::type tLC = t_lc::dummy;
     256
    254257  int  numBancroft = obsVector.size();
    255258
    … …  
    259262    for (unsigned ii = 0; ii < obsVector.size(); ii++) {
    260263      const t_pppSatObs* satObs = obsVector.at(ii);
    261       if (tLC == t_lc::dummy) {
    262         if (satObs->isValid(t_lc::cIF)) {
    263           tLC = t_lc::cIF;
    264         }
    265         else if (satObs->isValid(t_lc::c1)) {
    266           tLC = t_lc::c1;
    267         }
    268         else if (satObs->isValid(t_lc::c2)) {
    269           tLC = t_lc::c2;
    270         }
    271       }
    272       if ( satObs->isValid(tLC) &&
     264      t_lc LC = satObs->rangeLC();
     265      if ( satObs->isValid(LC) &&
    273266          (!satObs->modelSet() || satObs->eleSat() >= _opt->_minEle) ) {
    274267        ++iObs;
    … …  
    276269        BB[iObs][1] = satObs->xc()[1];
    277270        BB[iObs][2] = satObs->xc()[2];
    278         BB[iObs][3] = satObs->obsValue(tLC) - satObs->cmpValueForBanc(tLC);
     271        BB[iObs][3] = satObs->obsValue(LC) - satObs->cmpValueForBanc(LC);
    279272      }
    280273    }
    … …  
    297290    for (unsigned ii = 0; ii < obsVector.size(); ii++) {
    298291      const t_pppSatObs* satObs = obsVector.at(ii);
    299       if (satObs->isValid() &&
     292      t_lc LC = satObs->rangeLC();
     293      if (satObs->isValid(LC) &&
    300294          (!satObs->modelSet() || satObs->eleSat() >= _opt->_minEle) ) {
    301295        ColumnVector rr = satObs->xc().Rows(1,3) - xyzc.Rows(1,3);
    302         double res = rr.NormFrobenius() - satObs->obsValue(tLC)
     296        double res = rr.NormFrobenius() - satObs->obsValue(LC)
    303297                   - (satObs->xc()[3] - xyzc[3]) * t_CST::c;
    304298        if (fabs(res) > maxRes) {
    … …  
    354348//
    355349//////////////////////////////////////////////////////////////////////////////
    356 void t_pppClient::finish(t_irc irc, int ind) {
     350void t_pppClient::finish(t_irc irc, int /*ind*/) {
    357351#ifdef BNC_DEBUG_PPP
    358352  LOG << "t_pppClient::finish(" << ind << "): " << irc << endl;
    … …  
    463457      }
    464458    }
    465     // use observations only if satellite code biases are available
    466     /* ------------------------------------------------------------
    467     if (!_opt->_corrMount.empty() {
    468         useObsWithCodeBiasesOnly(_obsRover);
    469     }*/
     459
     460    // Use observations only if satellite code biases are available
     461    // ------------------------------------------------------------
     462    if (_opt->ambRes()) {
     463      useObsWithBiasesOnly(_obsRover);
     464    }
    470465
    471466    if (int(_obsRover.size()) < _opt->_minObs) {
    … …  
    613608
    614609}
     610
     611//
     612//////////////////////////////////////////////////////////////////////////////
     613void t_pppClient::useObsWithBiasesOnly(vector<t_pppSatObs*>& obsVector) const {
     614  vector<t_pppSatObs*>::iterator it = obsVector.begin();
     615  while (it != obsVector.end()) {
     616    t_pppSatObs* pppSatObs = *it;
     617    if (pppSatObs->hasBiases()) {
     618      ++it;
     619    }
     620    else {
     621      it = obsVector.erase(it);
     622      delete pppSatObs;
     623    }
     624  }
     625}

  • trunk/BNC/src/PPP/pppClient.h

    r10599 r10791  
    5252                   std::vector<t_pppSatObs*>& obsVector, bncTime& epoTime);
    5353  bool  preparePseudoObs(std::vector<t_pppSatObs*>& obsVector);
    54   void  useObsWithCodeBiasesOnly(std::vector<t_pppSatObs*>& obsVector);
    5554  t_irc cmpModel(t_pppStation* station, const ColumnVector& xyzc,
    5655                 std::vector<t_pppSatObs*>& obsVector);
    … …  
    6160  double cmpOffGal(std::vector<t_pppSatObs*>& obsVector);
    6261  double cmpOffBds(std::vector<t_pppSatObs*>& obsVector);
    63 
     62  void   useObsWithBiasesOnly(std::vector<t_pppSatObs*>& obsVector) const;
    6463
    6564  t_output*                 _output;
    … …  
    7574  t_tides*                  _tides;
    7675  bool                      _pseudoObsIono;
    77   QMap<char, int>           _usedSystems;
    7876  bool                      _running;
    7977};

  • trunk/BNC/src/PPP/pppEphPool.cpp

    r10330 r10791  
    1616
    1717#include <iostream>
     18#include <iomanip>
    1819#include "pppEphPool.h"
    1920#include "pppInclude.h"
    … …  
    4142  if (corr) {
    4243    _satEphPool[corr->_prn.toInt()].putOrbCorrection(corr);
     44    if (OPT->_logMode > t_pppOptions::normal) {
     45      LOG << "orbCorr " << string(corr->_time) << ' ' << corr->_prn.toString() << endl;
     46    }
    4347  }
    4448}
    … …  
    4953  if (corr) {
    5054    _satEphPool[corr->_prn.toInt()].putClkCorrection(corr);
     55    if (OPT->_logMode > t_pppOptions::normal) {
     56      LOG.setf(ios::fixed);
     57      LOG << "clkCorr " << string(corr->_time) << ' ' << corr->_prn.toString() << ' '
     58          << setw(7) << setprecision(3) << corr->_dClk * t_CST::c << endl;
     59    }
    5160  }
    5261}

  • trunk/BNC/src/PPP/pppFilter.cpp

    r10631 r10791  
    1616
    1717#include <iostream>
     18#include <sstream>
    1819#include <iomanip>
    1920#include <cmath>
    … …  
    2829#include "pppStation.h"
    2930#include "pppClient.h"
     31#include "ambres.h"
    3032
    3133using namespace BNC_PPP;
    … …  
    107109  }
    108110
    109   // close epoch processing
     111  // Dillution of Precision
    110112  // ----------------------
    111113  cmpDOP(allObs);
    112   _parlist->printResult(_epoTime, _QFlt, _xFlt);
     114 
     115  // Ambiguity Resolution
     116  // --------------------
     117  if (OPT->ambRes()) {
     118    AmbRes ambRes;
     119    ostringstream msg;
     120    ColumnVector    xFix     = _xFlt;
     121    SymmetricMatrix QFix     = _QFlt;
     122    double          fixRatio = 0.0;
     123    ambRes.run(_epoTime, _parlist->params(), QFix, xFix, fixRatio, msg);
     124    LOG << msg.str();
     125    _parlist->printResult(_epoTime, QFix, xFix, fixRatio);
     126  }
     127
     128  // Float Solution
     129  // --------------
     130  else {
     131    _parlist->printResult(_epoTime, _QFlt, _xFlt);
     132  }
     133
    113134  _lastEpoTimeOK = _epoTime; // remember time of last successful epoch processing
     135
    114136  return success;
    115137}
    … …  
    117139// Process Selected LCs
    118140////////////////////////////////////////////////////////////////////////////
    119 t_irc t_pppFilter::processSystem(const vector<t_lc::type> &LCs,
     141t_irc t_pppFilter::processSystem(const vector<t_lc>& LCs,
    120142                                 const vector<t_pppSatObs*> &obsVector,
    121143                                 bool pseudoObsIonoAvailable) {
    … …  
    158180    int iObs = -1;
    159181    vector<t_pppSatObs*> usedObs;
    160     vector<t_lc::type> usedTypes;
     182    vector<t_lc>        usedTypes;
    161183
    162184    // Real Observations
    … …  
    171193        nSat++;
    172194        for (unsigned jj = 0; jj < usedLCs; jj++) {
    173           const t_lc::type tLC = LCs[jj];
    174           if (tLC == t_lc::GIM) {
     195          const t_lc LC = LCs[jj];
     196          if (LC._type == t_lc::GIM) {
     197            continue;
     198          }
     199          if (LC._frq1 == t_frequency::G5 && !obs->isValid(LC)) {
    175200            continue;
    176201          }
    177202          ++iObs;
    178203          usedObs.push_back(obs);
    179           usedTypes.push_back(tLC);
     204          usedTypes.push_back(LC);
    180205          for (unsigned iPar = 0; iPar < nPar; iPar++) {
    181206            const t_pppParam *par = params[iPar];
    182             AA[iObs][iPar] = par->partial(_epoTime, obs, tLC);
    183           }
    184 
    185           ll[iObs] = obs->obsValue(tLC) - obs->cmpValue(tLC) - DotProduct(_x0, AA.Row(iObs + 1));
    186           PP[iObs] = 1.0 / (obs->sigma(tLC) * obs->sigma(tLC));
     207            AA[iObs][iPar] = par->partial(_epoTime, obs, LC);
     208          }
     209
     210          ll[iObs] = obs->obsValue(LC) - obs->cmpValue(LC) - DotProduct(_x0, AA.Row(iObs + 1));
     211          PP[iObs] = 1.0 / (obs->sigma(LC) * obs->sigma(LC));
    187212        }
    188213      }
    … …  
    202227        if (!obs->outlier()) {
    203228          for (unsigned jj = 0; jj < usedLCs; jj++) {
    204             const t_lc::type tLC = LCs[jj];
    205             if (tLC == t_lc::GIM) {
     229            const t_lc LC = LCs[jj];
     230            if (LC._type == t_lc::GIM) {
    206231              ++iObs;
    207232            } else {
    … …  
    209234            }
    210235            usedObs.push_back(obs);
    211             usedTypes.push_back(tLC);
     236            usedTypes.push_back(LC);
    212237            for (unsigned iPar = 0; iPar < nPar; iPar++) {
    213238              const t_pppParam *par = params[iPar];
    214               AA[iObs][iPar] = par->partial(_epoTime, obs, tLC);
     239              AA[iObs][iPar] = par->partial(_epoTime, obs, LC);
    215240            }
    216             ll[iObs] = obs->obsValue(tLC) - obs->cmpValue(tLC) - DotProduct(_x0, AA.Row(iObs + 1));
    217             PP[iObs] = 1.0 / (obs->sigma(tLC) * obs->sigma(tLC));
     241            ll[iObs] = obs->obsValue(LC) - obs->cmpValue(LC) - DotProduct(_x0, AA.Row(iObs + 1));
     242            PP[iObs] = 1.0 / (obs->sigma(LC) * obs->sigma(LC));
    218243          }
    219244        }
    … …  
    236261    double maxOutlier = 0.0;
    237262    int maxOutlierIndex = -1;
    238     t_lc::type maxOutlierLC = t_lc::dummy;
     263    t_lc maxOutlierLC;
    239264    for (unsigned ii = 0; ii < usedObs.size(); ii++) {
    240       const t_lc::type tLC = usedTypes[ii];
     265      const t_lc LC = usedTypes[ii];
    241266      double res = fabs(vv[ii]);
    242       if (res > usedObs[ii]->maxRes(tLC)) {
     267      if (res > usedObs[ii]->maxRes(LC)) {
    243268        if (res > fabs(maxOutlier)) {
    244269          maxOutlier = vv[ii];
    245270          maxOutlierIndex = ii;
    246           maxOutlierLC = tLC;
     271          maxOutlierLC = LC;
    247272        }
    248273      }
    … …  
    254279      t_pppSatObs *obs = usedObs[maxOutlierIndex];
    255280      t_pppParam *par = 0;
    256       LOG << epoTimeStr << " Outlier " << t_lc::toString(maxOutlierLC) << ' '
     281      LOG << epoTimeStr << " Outlier " << maxOutlierLC.toString() << ' '
    257282          << obs->prn().toString() << ' ' << setw(8) << setprecision(4)
    258283          << maxOutlier << endl;
    259284      for (unsigned iPar = 0; iPar < nPar; iPar++) {
    260         t_pppParam *hlp = params[iPar];
     285        t_pppParam* hlp = params[iPar];
    261286        if (hlp->type() == t_pppParam::amb &&
    262287            hlp->prn()  == obs->prn() &&
    263             hlp->tLC()  == usedTypes[maxOutlierIndex]) {
     288            hlp->LC()   == usedTypes[maxOutlierIndex]) {
    264289          par = hlp;
    265290        }
    … …  
    277302        for (unsigned jj = 0; jj < LCs.size(); jj++) {
    278303          for (unsigned ii = 0; ii < usedObs.size(); ii++) {
    279             const t_lc::type tLC = usedTypes[ii];
    280             t_pppSatObs *obs = usedObs[ii];
    281             if (tLC == LCs[jj]) {
    282               obs->setRes(tLC, vv[ii]);
     304            const t_lc LC = usedTypes[ii];
     305            t_pppSatObs* obs = usedObs[ii];
     306            if (LC == LCs[jj]) {
     307              obs->setRes(LC, vv[ii]);
    283308              LOG << epoTimeStr << " RES " << left << setw(3)
    284                   << t_lc::toString(tLC) << right << ' '
     309                  << LC.toString() << right << ' '
    285310                  << obs->prn().toString() << ' '
    286311                  << setw(8) << setprecision(4) << vv[ii] << endl;
    … …  
    296321// Cycle-Slip Detection
    297322////////////////////////////////////////////////////////////////////////////
    298 t_irc t_pppFilter::detectCycleSlips(const vector<t_lc::type> &LCs,
     323t_irc t_pppFilter::detectCycleSlips(const vector<t_lc>& LCs,
    299324                                    const vector<t_pppSatObs*> &obsVector) {
    300325
    … …  
    309334  SLIP *= fac;
    310335  string epoTimeStr = string(_epoTime);
    311   const vector<t_pppParam*> &params = _parlist->params();
    312336
    313337  for (unsigned ii = 0; ii < LCs.size(); ii++) {
    314     const t_lc::type &tLC = LCs[ii];
    315     if (t_lc::includesPhase(tLC)) {
     338    const t_lc& LC = LCs[ii];
     339    if (LC.includesPhase()) {
    316340      for (unsigned iObs = 0; iObs < obsVector.size(); iObs++) {
    317341        const t_pppSatObs *obs = obsVector[iObs];
    … …  
    343367        // --------
    344368        if (slip) {
    345           resetAmb(obs->prn(), obsVector, tLC);
    346         }
    347 
    348         // Check Pre-Fit Residuals
    349         /* -----------------------
    350         else {
    351           ColumnVector AA(params.size());
    352           for (unsigned iPar = 0; iPar < params.size(); iPar++) {
    353             const t_pppParam* par = params[iPar];
    354             AA[iPar] = par->partial(_epoTime, obs, tLC);
    355           }
    356           double ll = obs->obsValue(tLC) - obs->cmpValue(tLC) - DotProduct(_x0, AA);
    357           double vv = DotProduct(AA, _xFlt) - ll;
    358 
    359           if (fabs(vv) > SLIP) {
    360             LOG << epoTimeStr << " cycle slip detected " << t_lc::toString(tLC) << ' '
    361                 << obs->prn().toString() << ' ' << setw(8) << setprecision(4) << vv << endl;
    362             resetAmb(obs->prn(), obsVector, tLC);
    363           }
    364         }*/
     369          resetAmb(obs->prn(), obsVector, t_lc());
     370        }
    365371      }
    366372    }
    … …  
    371377// Reset Ambiguity Parameter (cycle slip)
    372378////////////////////////////////////////////////////////////////////////////
    373 t_irc t_pppFilter::resetAmb(const t_prn prn, const vector<t_pppSatObs*> &obsVector, t_lc::type lc,
     379t_irc t_pppFilter::resetAmb(const t_prn prn, const vector<t_pppSatObs*> &obsVector, t_lc lc,
    374380                            SymmetricMatrix *QSav, ColumnVector *xSav) {
    375381
    … …  
    386392      (par->firstObsTime().undef()) ?
    387393        firstObsTime = lastObsTime : firstObsTime = par->firstObsTime();
    388       t_lc::type tLC = par->tLC();
    389       if (tLC != lc) {continue;}
     394      t_lc LC = par->LC();
     395      if (lc.valid() && LC != lc) {continue;}
    390396      LOG << string(_epoTime) << " RESET " << par->toString() << endl;
    391       delete par; par = new t_pppParam(t_pppParam::amb, prn, tLC, &obsVector);
     397      delete par; par = new t_pppParam(t_pppParam::amb, prn, LC, &obsVector);
    392398      par->setIndex(ind);
    393399      par->setFirstObsTime(firstObsTime);
    … …  
    436442        for (unsigned iPar = 0; iPar < numPar; iPar++) {
    437443          t_pppParam* par = _parlist->params()[iPar];
    438           AA[_numSat - 1][iPar] = par->partial(_epoTime, obs, t_lc::c1);
     444          AA[_numSat - 1][iPar] = par->partial(_epoTime, obs, obs->rangeLC());
    439445          if      (par->type() == t_pppParam::crdX) {
    440446            parX = par;

  • trunk/BNC/src/PPP/pppFilter.h

    r10256 r10791  
    7373  };
    7474
    75   t_irc processSystem(const std::vector<t_lc::type>& LCs,
     75  t_irc processSystem(const std::vector<t_lc>& LCs,
    7676                      const std::vector<t_pppSatObs*>& obsVector,
    7777                      bool pseudoObsIonoAvailable);
    7878
    79   t_irc detectCycleSlips(const std::vector<t_lc::type>& LCs,
     79  t_irc detectCycleSlips(const std::vector<t_lc>& LCs,
    8080                         const std::vector<t_pppSatObs*>& obsVector);
    8181
    82   t_irc resetAmb(t_prn prn, const std::vector<t_pppSatObs*>& obsVector, t_lc::type lc,
     82  t_irc resetAmb(t_prn prn, const std::vector<t_pppSatObs*>& obsVector, t_lc lc,
    8383                 SymmetricMatrix* QSav = 0, ColumnVector* xSav = 0);
    8484

  • trunk/BNC/src/PPP/pppObsPool.cpp

    r10034 r10791  
    1515 * -----------------------------------------------------------------------*/
    1616
     17#include <iomanip>
    1718#include "pppObsPool.h"
     19#include "pppClient.h"
    1820
    1921using namespace BNC_PPP;
    … …  
    7476  delete _satCodeBiases[iPrn];
    7577  _satCodeBiases[iPrn] = satCodeBias;
     78  if (OPT->_logMode > t_pppOptions::normal) {
     79    LOG << "codeBias " << string(satCodeBias->_time) << ' ' << satCodeBias->_prn.toString();
     80    for (const auto& bias : satCodeBias->_bias) {
     81      LOG << "  " << bias._rnxType2ch << ' ' << setw(7) << setprecision(3) << bias._value;
     82    }
     83    LOG << endl;
     84  }
    7685}
    7786
    … …  
    8291  delete _satPhaseBiases[iPrn];
    8392  _satPhaseBiases[iPrn] = satPhaseBias;
     93  if (OPT->_logMode > t_pppOptions::normal) {
     94    LOG.setf(ios::fixed);
     95    LOG << "phaseBias " << string(satPhaseBias->_time) << ' ' << satPhaseBias->_prn.toString()
     96        << ' ' << setw(7) << setprecision(2) << satPhaseBias->_yaw     * 180.0 / M_PI
     97        << ' ' << setw(7) << setprecision(3) << satPhaseBias->_yawRate * 180.0 / M_PI;
     98    for (const auto& bias : satPhaseBias->_bias) {
     99      LOG << "  " << bias._rnxType2ch << ' ' << setw(2) << bias._jumpCounter;
     100      if (bias._fixIndicator) {
     101        LOG << " x ";
     102      }
     103      else {
     104        LOG << " . ";
     105      }
     106      LOG << setw(6) << setprecision(3) << bias._value;
     107    }
     108    LOG << endl;
     109  }
     110}
     111
     112//
     113/////////////////////////////////////////////////////////////////////////////
     114void t_pppObsPool::clearCodeBiases(char sys) {
     115  for (unsigned iPrn = 1; iPrn <= t_prn::MAXPRN; ++iPrn) {
     116    t_prn prn; prn.set(iPrn);
     117    if (prn.system() == sys) {
     118      delete _satCodeBiases[iPrn];
     119      _satCodeBiases[iPrn] = 0;
     120    }
     121  }
     122}
     123
     124//
     125/////////////////////////////////////////////////////////////////////////////
     126void t_pppObsPool::clearPhaseBiases(char sys) {
     127  for (unsigned iPrn = 1; iPrn <= t_prn::MAXPRN; ++iPrn) {
     128    t_prn prn; prn.set(iPrn);
     129    if (prn.system() == sys) {
     130      delete _satPhaseBiases[iPrn];
     131      _satPhaseBiases[iPrn] = 0;
     132    }
     133  }
    84134}
    85135

  • trunk/BNC/src/PPP/pppObsPool.h

    r10034 r10791  
    3333  void putPhaseBias(t_satPhaseBias* satPhaseBias);
    3434  void putTec(t_vTec* _vTec);
     35  void clearCodeBiases(char sys);
     36  void clearPhaseBiases(char sys);
    3537
    3638  void putEpoch(const bncTime& epoTime, std::vector<t_pppSatObs*>& obsVector,

  • trunk/BNC/src/PPP/pppParlist.cpp

    r10623 r10791  
    3434// Constructor
    3535////////////////////////////////////////////////////////////////////////////
    36 t_pppParam::t_pppParam(e_type type, const t_prn& prn, t_lc::type tLC,
     36t_pppParam::t_pppParam(e_type type, const t_prn& prn, t_lc LC,
    3737                       const vector<t_pppSatObs*>* obsVector) {
    3838
    3939  _type     = type;
    4040  _prn      = prn;
    41   _tLC      = tLC;
     41  _sys      = '?';
     42  _LC       = LC;
    4243  _x0       = 0.0;
    4344  _indexOld = -1;
    … …  
    4748
    4849  switch (_type) {
    49    case crdX:
    50      _epoSpec = false;
    51      _sigma0  = OPT->_aprSigCrd[0];
    52      _noise   = OPT->_noiseCrd[0];
    53      break;
    54    case crdY:
    55      _epoSpec = false;
    56      _sigma0  = OPT->_aprSigCrd[1];
    57      _noise   = OPT->_noiseCrd[1];
    58      break;
    59    case crdZ:
    60      _epoSpec = false;
    61      _sigma0  = OPT->_aprSigCrd[2];
    62      _noise   = OPT->_noiseCrd[2];
    63      break;
    64    case rClkG:
    65      _epoSpec = true;
    66      _sigma0  = OPT->_aprSigClk;
    67      break;
    68    case rClkR:
    69      _epoSpec = true;
    70      _sigma0  = OPT->_aprSigClk;
    71      break;
    72    case rClkE:
    73      _epoSpec = true;
    74      _sigma0  = OPT->_aprSigClk;
    75      break;
    76    case rClkC:
    77      _epoSpec = true;
    78      _sigma0  = OPT->_aprSigClk;
    79      break;
    80    case amb:
    81      _epoSpec = false;
    82      _sigma0  = OPT->_aprSigAmb;
    83      if (obsVector) {
    84        for (unsigned ii = 0; ii < obsVector->size(); ii++) {
    85          const t_pppSatObs* obs = obsVector->at(ii);
    86          if (obs->prn() == _prn) {
    87            _x0 = floor((obs->obsValue(tLC) - obs->cmpValue(tLC)) / obs->lambda(tLC) + 0.5);
    88            break;
    89          }
    90        }
    91      }
    92      break;
    93    case trp:
    94      _epoSpec = false;
    95      _sigma0  = OPT->_aprSigTrp;
    96      _noise   = OPT->_noiseTrp;
    97      break;
    98     case ion:
    99       _epoSpec = false;
    100       _sigma0  = OPT->_aprSigIon;
    101       _noise   = OPT->_noiseIon;
    102      break;
    103     case cBiasG1:   case cBiasR1:   case cBiasE1:   case cBiasC1:
    104     case cBiasG2:   case cBiasR2:   case cBiasE2:   case cBiasC2:
    105       _epoSpec = false;
    106       _sigma0  = OPT->_aprSigCodeBias;
    107       _noise   = OPT->_noiseCodeBias;
    108       break;
    109     case pBiasG1:   case pBiasR1:   case pBiasE1:   case pBiasC1:
    110     case pBiasG2:   case pBiasR2:   case pBiasE2:   case pBiasC2:
    111       _epoSpec = false;
    112       _sigma0  = OPT->_aprSigPhaseBias;
    113       _noise   = OPT->_noisePhaseBias;
    114       break;
     50  case crdX:
     51    _epoSpec = false;
     52    _sigma0  = OPT->_aprSigCrd[0];
     53    _noise   = OPT->_noiseCrd[0];
     54    break;
     55  case crdY:
     56    _epoSpec = false;
     57    _sigma0  = OPT->_aprSigCrd[1];
     58    _noise   = OPT->_noiseCrd[1];
     59    break;
     60  case crdZ:
     61    _epoSpec = false;
     62    _sigma0  = OPT->_aprSigCrd[2];
     63    _noise   = OPT->_noiseCrd[2];
     64    break;
     65  case rClk:
     66    _epoSpec = true;
     67    _sigma0  = OPT->_aprSigClk;
     68    break;
     69  case amb:
     70    _epoSpec = false;
     71    _sigma0  = OPT->_aprSigAmb;
     72    if (obsVector) {
     73      for (unsigned ii = 0; ii < obsVector->size(); ii++) {
     74        const t_pppSatObs* obs = obsVector->at(ii);
     75        if (obs->prn() == _prn) {
     76          _x0 = floor((obs->obsValue(LC) - obs->cmpValue(LC)) / obs->lambda(LC) + 0.5);
     77          break;
     78        }
     79      }
     80    }
     81    break;
     82  case trp:
     83    _epoSpec = false;
     84    _sigma0  = OPT->_aprSigTrp;
     85    _noise   = OPT->_noiseTrp;
     86    break;
     87  case ion:
     88    _epoSpec = true;
     89    _sigma0  = OPT->_aprSigIon;
     90    break;
     91  case bias:
     92    _epoSpec = true;
     93    _sigma0  = OPT->_aprSigCodeBias;
     94    break;
    11595  }
    11696}
    … …  
    124104////////////////////////////////////////////////////////////////////////////
    125105double t_pppParam::partial(const bncTime& /* epoTime */, const t_pppSatObs* obs,
    126                            const t_lc::type& tLC) const {
     106                           const t_lc& obsLC) const {
    127107
    128108  // Special Case - Melbourne-Wuebbena
    129109  // ---------------------------------
    130   if (tLC == t_lc::MW && _type != amb) {
     110  if (obsLC._type == t_lc::MW && _type != amb) {
    131111    return 0.0;
     112  }
     113
     114  // Special Case - GIM
     115  // ------------------
     116  if (obsLC._type == t_lc::GIM) {
     117    if (_type == ion) {
     118      return 1.0;
     119    }
     120    else {
     121      return 0.0;
     122    }
    132123  }
    133124
    … …  
    138129  map<t_frequency::type, double> phaseCoeff;
    139130  map<t_frequency::type, double> ionoCoeff;
    140   obs->lcCoeff(tLC, codeCoeff, phaseCoeff, ionoCoeff);
    141 
     131  obs->lcCoeff(obsLC, codeCoeff, phaseCoeff, ionoCoeff);
     132 
    142133  switch (_type) {
    143134  case crdX:
    144     if (tLC == t_lc::GIM) {return 0.0;}
    145135    return (sta->xyzApr()[0] - obs->xc()[0]) / rhoV.NormFrobenius();
    146136  case crdY:
    147     if (tLC == t_lc::GIM) {return 0.0;}
    148137    return (sta->xyzApr()[1] - obs->xc()[1]) / rhoV.NormFrobenius();
    149138  case crdZ:
    150     if (tLC == t_lc::GIM) {return 0.0;}
    151139    return (sta->xyzApr()[2] - obs->xc()[2]) / rhoV.NormFrobenius();
    152   case rClkG:
    153     if (tLC == t_lc::GIM) {return 0.0;}
    154     return (obs->prn().system() == 'G') ? 1.0 : 0.0;
    155   case rClkR:
    156     if (tLC == t_lc::GIM) {return 0.0;}
    157     return (obs->prn().system() == 'R') ? 1.0 : 0.0;
    158   case rClkE:
    159     if (tLC == t_lc::GIM) {return 0.0;}
    160     return (obs->prn().system() == 'E') ? 1.0 : 0.0;
    161   case rClkC:
    162     if (tLC == t_lc::GIM) {return 0.0;}
    163     return (obs->prn().system() == 'C') ? 1.0 : 0.0;
     140  case rClk:
     141    return (_sys != obsLC.system() || obsLC.isGeometryFree()) ? 0.0 : 1.0;
    164142  case amb:
    165     if (tLC == t_lc::GIM) {
    166       return 0.0;
    167     }
    168     else {
    169       if (obs->prn() == _prn) {
    170         if      (tLC == _tLC) {
    171           return (obs->lambda(tLC));
    172         }
    173         else if (tLC == t_lc::lIF && _tLC == t_lc::MW) {
    174           return obs->lambda(t_lc::lIF) * obs->lambda(t_lc::MW) / obs->lambda(t_lc::l2);
    175         }
    176         else {
    177           if      (_tLC == t_lc::l1) {
    178             return obs->lambda(t_lc::l1) * phaseCoeff[obs->fType1()];
    179           }
    180           else if (_tLC == t_lc::l2) {
    181             return obs->lambda(t_lc::l2) * phaseCoeff[obs->fType2()];
    182           }
    183         }
     143    if (obs->prn() == _prn) {
     144      if      (obsLC == _LC) {
     145        return (obs->lambda(obsLC));
     146      }
     147      else if (_LC._type == t_lc::phase) {
     148        return obs->lambda(_LC) * phaseCoeff[_LC._frq1];
     149      }
     150      else if (obsLC._type == t_lc::phaseIF && _LC._type == t_lc::MW) {
     151        return obs->lambda(obsLC) * obs->lambda(_LC) / obs->lambda(t_lc(t_lc::phase, _LC._frq2));
    184152      }
    185153    }
    186154    break;
    187155  case trp:
    188     if      (tLC == t_lc::GIM) {
    189       return 0.0;
    190     }
    191     else {
    192       return 1.0 / sin(obs->eleSat());
    193     }
     156    return 1.0 / sin(obs->eleSat());
    194157  case ion:
    195158    if (obs->prn() == _prn) {
    196       if      (tLC == t_lc::c1) {
    197         return ionoCoeff[obs->fType1()];
    198       }
    199       else if (tLC == t_lc::c2) {
    200         return ionoCoeff[obs->fType2()];
    201       }
    202       else if (tLC == t_lc::l1) {
    203         return ionoCoeff[obs->fType1()];
    204       }
    205       else if (tLC == t_lc::l2) {
    206         return ionoCoeff[obs->fType2()];
    207       }
    208       else if (tLC == t_lc::GIM) {
    209         return 1.0;
    210       }
    211     }
    212     break;
    213   case cBiasG1:
    214     if ((obs->prn().system() == 'G') && (tLC == t_lc::c1)) {return 1.0;} else {return 0.0;}
    215     break;
    216   case cBiasR1:
    217     if ((obs->prn().system() == 'R') && (tLC == t_lc::c1)) {return 1.0;} else {return 0.0;}
    218     break;
    219   case cBiasE1:
    220     if ((obs->prn().system() == 'E') && (tLC == t_lc::c1)) {return 1.0;} else {return 0.0;}
    221     break;
    222   case cBiasC1:
    223     if ((obs->prn().system() == 'C') && (tLC == t_lc::c1)) {return 1.0;} else {return 0.0;}
    224     break;
    225   case cBiasG2:
    226     if ((obs->prn().system() == 'G') && (tLC == t_lc::c2)) {return 1.0;} else {return 0.0;}
    227     break;
    228   case cBiasR2:
    229     if ((obs->prn().system() == 'R') && (tLC == t_lc::c2)) {return 1.0;} else {return 0.0;}
    230     break;
    231   case cBiasE2:
    232     if ((obs->prn().system() == 'E') && (tLC == t_lc::c2)) {return 1.0;} else {return 0.0;}
    233         break;
    234   case cBiasC2:
    235     if ((obs->prn().system() == 'C') && (tLC == t_lc::c2)) {return 1.0;} else {return 0.0;}
    236     break;
    237   case pBiasG1:
    238     if ((obs->prn().system() == 'G') && (tLC == t_lc::l1)) {return 1.0;} else {return 0.0;}
    239     break;
    240   case pBiasR1:
    241     if ((obs->prn().system() == 'R') && (tLC == t_lc::l1)) {return 1.0;} else {return 0.0;}
    242     break;
    243   case pBiasE1:
    244     if ((obs->prn().system() == 'E') && (tLC == t_lc::l1)) {return 1.0;} else {return 0.0;}
    245     break;
    246   case pBiasC1:
    247     if ((obs->prn().system() == 'C') && (tLC == t_lc::l1)) {return 1.0;} else {return 0.0;}
    248     break;
    249   case pBiasG2:
    250     if ((obs->prn().system() == 'G') && (tLC == t_lc::l2)) {return 1.0;} else {return 0.0;}
    251     break;
    252   case pBiasR2:
    253     if ((obs->prn().system() == 'R') && (tLC == t_lc::l2)) {return 1.0;} else {return 0.0;}
    254     break;
    255   case pBiasE2:
    256     if ((obs->prn().system() == 'E') && (tLC == t_lc::l2)) {return 1.0;} else {return 0.0;}
    257     break;
    258   case pBiasC2:
    259     if ((obs->prn().system() == 'C') && (tLC == t_lc::l2)) {return 1.0;} else {return 0.0;}
    260     break;
    261 
     159      if (obsLC._type == t_lc::code || obsLC._type == t_lc::phase) {
     160        return ionoCoeff[obsLC._frq1];
     161      }
     162    }
     163    break;
     164  case bias:
     165    return (_LC == obsLC ? 1.0 : 0.0);
    262166  }
    263167  return 0.0;
    … …  
    278182    ss << "CRD_Z";
    279183    break;
    280   case rClkG:
    281     ss << "REC_CLK  G  ";
    282     break;
    283   case rClkR:
    284     ss << "REC_CLK  R  ";
    285     break;
    286   case rClkE:
    287     ss << "REC_CLK  E  ";
    288     break;
    289   case rClkC:
    290     ss << "REC_CLK  C  ";
     184  case rClk:
     185    ss << "REC_CLK  " << _sys << "  ";
    291186    break;
    292187  case trp:
    … …  
    294189    break;
    295190  case amb:
    296     ss << "AMB  " << left << setw(3) << t_lc::toString(_tLC) << right << ' ' << _prn.toString();
     191    ss << "AMB  " << left << setw(3) << _LC.toString() << right << ' ' << _prn.toString();
    297192    break;
    298193  case ion:
    299     ss << "ION  " << left << setw(3) << t_lc::toString(_tLC) << right << ' ' << _prn.toString();
    300     break;
    301   case cBiasG1:  case pBiasG1:
    302   case cBiasG2:  case pBiasG2:
    303     ss << "BIA  " << left << setw(3) << t_lc::toString(_tLC) << right << " G  ";
    304     break;
    305   case cBiasR1:  case pBiasR1:
    306   case cBiasR2:  case pBiasR2:
    307     ss << "BIA  " << left << setw(3) << t_lc::toString(_tLC) << right << " R  ";
    308     break;
    309   case cBiasE1:  case pBiasE1:
    310   case cBiasE2:  case pBiasE2:
    311     ss << "BIA  " << left << setw(3) << t_lc::toString(_tLC) << right << " E  ";
    312     break;
    313   case cBiasC1:  case pBiasC1:
    314   case cBiasC2:  case pBiasC2:
    315     ss << "BIA  " << left << setw(3) << t_lc::toString(_tLC) << right << " C  ";
     194    ss << "ION  " << left << setw(3) << _LC.toString() << right << ' ' << _prn.toString();
     195    break;
     196  case bias:
     197    char sys = t_frequency::toSystem(_LC._frq1);
     198    ss << "BIA  " << left << setw(3) << _LC.toString() << right << ' ' << sys << "  ";
    316199    break;
    317200  }
    … …  
    353236             par->type() == t_pppParam::crdY ||
    354237             par->type() == t_pppParam::crdZ ||
    355              par->type() == t_pppParam::ion  ||
    356              par->type() == t_pppParam::cBiasC1 ||
    357              par->type() == t_pppParam::cBiasC2 ||
    358              par->type() == t_pppParam::cBiasE1 ||
    359              par->type() == t_pppParam::cBiasE2 ||
    360              par->type() == t_pppParam::cBiasR1 ||
    361              par->type() == t_pppParam::cBiasR2 ||
    362              par->type() == t_pppParam::pBiasC1 ||
    363              par->type() == t_pppParam::pBiasC2 ||
    364              par->type() == t_pppParam::pBiasE1 ||
    365              par->type() == t_pppParam::pBiasE2 ||
    366              par->type() == t_pppParam::pBiasR1 ||
    367              par->type() == t_pppParam::pBiasR2) {
     238             par->type() == t_pppParam::ion  ) {
    368239      if (par->lastObsTime().valid() && (epoTime - par->lastObsTime() > 60.0)) {
    369240        remove = true;
    … …  
    383254  }
    384255
    385   // check which systems have observations
     256  // Check which systems have observations
    386257  // -------------------------------------
    387   _usedSystems['G'] = _usedSystems['R'] = _usedSystems['E'] = _usedSystems['C'] = 0;
    388   for (unsigned jj = 0; jj < obsVector.size(); jj++) {
    389     const t_pppSatObs* satObs = obsVector[jj];
    390     char sys = satObs->prn().system();
    391     if (OPT->LCs(sys).size()) {
    392       _usedSystems[sys]++;
     258  vector<char> systems = OPT->systems();
     259  for (char sys : systems) {
     260    _usedSystems[sys] = 0;
     261    for (unsigned jj = 0; jj < obsVector.size(); jj++) {
     262      const t_pppSatObs* satObs = obsVector[jj];
     263      if (satObs->prn().system() == sys) {
     264        _usedSystems[sys]++;
     265      }
    393266    }
    394267  };
    395 
    396268
    397269  // Check whether parameters have observations
    … …  
    425297  // Coordinates
    426298  // -----------
    427   required.push_back(new t_pppParam(t_pppParam::crdX, t_prn(), t_lc::dummy));
    428   required.push_back(new t_pppParam(t_pppParam::crdY, t_prn(), t_lc::dummy));
    429   required.push_back(new t_pppParam(t_pppParam::crdZ, t_prn(), t_lc::dummy));
     299  required.push_back(new t_pppParam(t_pppParam::crdX, t_prn()));
     300  required.push_back(new t_pppParam(t_pppParam::crdY, t_prn()));
     301  required.push_back(new t_pppParam(t_pppParam::crdZ, t_prn()));
    430302
    431303  // Receiver Clocks
    432304  // ---------------
    433    if (_usedSystems['G']) {
    434    //if (OPT->useSystem('G')) {
    435      required.push_back(new t_pppParam(t_pppParam::rClkG, t_prn(), t_lc::dummy));
    436    }
    437    if (_usedSystems['R']) {
    438    //if (OPT->useSystem('R')) {
    439      required.push_back(new t_pppParam(t_pppParam::rClkR, t_prn(), t_lc::dummy));
    440    }
    441    if (_usedSystems['E']) {
    442    //if (OPT->useSystem('E')) {
    443      required.push_back(new t_pppParam(t_pppParam::rClkE, t_prn(), t_lc::dummy));
    444    }
    445    if (_usedSystems['C']) {
    446    //if (OPT->useSystem('C')) {
    447      required.push_back(new t_pppParam(t_pppParam::rClkC, t_prn(), t_lc::dummy));
    448    }
     305  for (const auto& [sys, numObs] : _usedSystems) {
     306    if (numObs > 0) {
     307      t_pppParam* clk = new t_pppParam(t_pppParam::rClk, t_prn());
     308      clk->setSystem(sys);
     309      required.push_back(clk);
     310    }
     311  }
    449312
    450313  // Troposphere
    451314  // -----------
    452315  if (OPT->estTrp()) {
    453     required.push_back(new t_pppParam(t_pppParam::trp, t_prn(), t_lc::dummy));
     316    required.push_back(new t_pppParam(t_pppParam::trp, t_prn()));
    454317  }
    455318
    … …  
    459322    for (unsigned jj = 0; jj < obsVector.size(); jj++) {
    460323      const t_pppSatObs* satObs = obsVector[jj];
    461       char sys = satObs->prn().system();
    462       std::vector<t_lc::type> LCs = OPT->LCs(sys);
    463       if (std::find(LCs.begin(), LCs.end(), t_lc::cIF) == LCs.end() &&
    464           std::find(LCs.begin(), LCs.end(), t_lc::lIF) == LCs.end()) {
    465         required.push_back(new t_pppParam(t_pppParam::ion, satObs->prn(), t_lc::dummy));
    466       }
    467     }
    468   }
     324      std::vector<t_lc> LCs = OPT->LCs(satObs->prn().system());
     325      for (auto it = LCs.begin(); it != LCs.end(); ++it) {
     326        const t_lc& lc = *it;
     327        if (!lc.isIonoFree()) {
     328          required.push_back(new t_pppParam(t_pppParam::ion, satObs->prn()));
     329          break;
     330        }
     331      }
     332    }
     333  }
     334 
    469335  // Ambiguities
    470336  // -----------
    … …  
    472338    const t_pppSatObs*  satObs = obsVector[jj];
    473339    char sys = satObs->prn().system();
    474     const vector<t_lc::type>& ambLCs = OPT->ambLCs(sys);
     340    const vector<t_lc>& ambLCs = OPT->ambLCs(sys);
    475341    for (unsigned ii = 0; ii < ambLCs.size(); ii++) {
     342      if (ambLCs[ii]._frq1 == t_frequency::G5 && !satObs->isValid(ambLCs[ii])) {
     343        continue;
     344      }
    476345      required.push_back(new t_pppParam(t_pppParam::amb, satObs->prn(), ambLCs[ii], &obsVector));
    477346    }
    478347  }
    479348
    480   // Receiver Code Biases
    481   // --------------------
    482     if (OPT->_ionoModelType == OPT->PPP_RTK) {
    483     std::vector<t_lc::type> lc;
    484     if (OPT->useSystem('G')) {
    485       lc = OPT->LCs('G');
    486       if (std::find(lc.begin(), lc.end(), t_lc::c1) != lc.end()) {
    487         required.push_back(new t_pppParam(t_pppParam::cBiasG1, t_prn(), t_lc::c1));
    488       }
    489       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    490         required.push_back(new t_pppParam(t_pppParam::cBiasG2, t_prn(), t_lc::c2));
    491       }
    492     }
    493     if (OPT->useSystem('R')) {
    494       lc = OPT->LCs('R');
    495       if (std::find(lc.begin(), lc.end(), t_lc::c1) != lc.end()) {
    496         required.push_back(new t_pppParam(t_pppParam::cBiasR1, t_prn(), t_lc::c1));
    497       }
    498       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    499         required.push_back(new t_pppParam(t_pppParam::cBiasR2, t_prn(), t_lc::c2));
    500       }
    501     }
    502     if (OPT->useSystem('E')) {
    503       lc = OPT->LCs('E');
    504       if (std::find(lc.begin(), lc.end(), t_lc::c1) != lc.end()) {
    505         required.push_back(new t_pppParam(t_pppParam::cBiasE1, t_prn(), t_lc::c1));
    506       }
    507       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    508         required.push_back(new t_pppParam(t_pppParam::cBiasE2, t_prn(), t_lc::c2));
    509       }
    510     }
    511     if (OPT->useSystem('C')) {
    512       lc = OPT->LCs('C');
    513       if (std::find(lc.begin(), lc.end(), t_lc::c1) != lc.end()) {
    514         required.push_back(new t_pppParam(t_pppParam::cBiasC1, t_prn(), t_lc::c1));
    515       }
    516       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    517         required.push_back(new t_pppParam(t_pppParam::cBiasC2, t_prn(), t_lc::c2));
    518       }
    519     }
    520   }
    521 
    522   if (OPT->_pseudoObsIono) {
    523     std::vector<t_lc::type> lc;
    524     if (OPT->useSystem('G')) {
    525       lc = OPT->LCs('G');
    526       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    527         required.push_back(new t_pppParam(t_pppParam::cBiasG2, t_prn(), t_lc::c2));
    528       }
    529     }
    530     if (OPT->useSystem('R')) {
    531       lc = OPT->LCs('R');
    532       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    533         required.push_back(new t_pppParam(t_pppParam::cBiasR2, t_prn(), t_lc::c2));
    534       }
    535     }
    536     if (OPT->useSystem('E')) {
    537       lc = OPT->LCs('E');
    538       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    539         required.push_back(new t_pppParam(t_pppParam::cBiasE2, t_prn(), t_lc::c2));
    540       }
    541     }
    542     if (OPT->useSystem('C')) {
    543       lc = OPT->LCs('C');
    544       if (std::find(lc.begin(), lc.end(), t_lc::c2) != lc.end()) {
    545         required.push_back(new t_pppParam(t_pppParam::cBiasC2, t_prn(), t_lc::c2));
    546       }
    547     }
    548   }
    549 
    550   if (OPT->_ionoModelType == OPT->PPP_RTK) {
    551     std::vector<t_lc::type> lc;
    552     if (OPT->useSystem('G')) {
    553       lc = OPT->LCs('G');
    554       if (std::find(lc.begin(), lc.end(), t_lc::l1) != lc.end()) {
    555         required.push_back(new t_pppParam(t_pppParam::pBiasG1, t_prn(), t_lc::l1));
    556       }
    557       if (std::find(lc.begin(), lc.end(), t_lc::l2) != lc.end()) {
    558         required.push_back(new t_pppParam(t_pppParam::pBiasG2, t_prn(), t_lc::l2));
    559       }
    560     }
    561     if (OPT->useSystem('R')) {
    562       lc = OPT->LCs('R');
    563       if (std::find(lc.begin(), lc.end(), t_lc::l1) != lc.end()) {
    564         required.push_back(new t_pppParam(t_pppParam::pBiasR1, t_prn(), t_lc::l1));
    565       }
    566       if (std::find(lc.begin(), lc.end(), t_lc::l2) != lc.end()) {
    567         required.push_back(new t_pppParam(t_pppParam::pBiasR2, t_prn(), t_lc::l2));
    568       }
    569     }
    570     if (OPT->useSystem('E')) {
    571       lc = OPT->LCs('E');
    572       if (std::find(lc.begin(), lc.end(), t_lc::l1) != lc.end()) {
    573         required.push_back(new t_pppParam(t_pppParam::pBiasE1, t_prn(), t_lc::l1));
    574       }
    575       if (std::find(lc.begin(), lc.end(), t_lc::l2) != lc.end()) {
    576         required.push_back(new t_pppParam(t_pppParam::pBiasE2, t_prn(), t_lc::l2));
    577       }
    578     }
    579     if (OPT->useSystem('C')) {
    580       lc = OPT->LCs('C');
    581       if (std::find(lc.begin(), lc.end(), t_lc::l1) != lc.end()) {
    582         required.push_back(new t_pppParam(t_pppParam::pBiasC1, t_prn(), t_lc::l1));
    583       }
    584       if (std::find(lc.begin(), lc.end(), t_lc::l2) != lc.end()) {
    585         required.push_back(new t_pppParam(t_pppParam::pBiasC2, t_prn(), t_lc::l2));
     349  // Biases
     350  // ------
     351  int maxSkip = (OPT->_pseudoObsIono ? 1 : 2);
     352  for (const auto& [sys, numObs] : _usedSystems) {
     353    if (numObs > 0) {
     354      bool ar   = OPT->arSystem(sys);
     355      int  skip = 0;
     356      vector<t_lc> LCs = OPT->LCs(sys);
     357      for (const t_lc& lc : LCs) {
     358        if (ar) {
     359          if (skip < maxSkip && lc.includesPhase()) {
     360            skip += 1;
     361          }
     362          else {
     363            required.push_back(new t_pppParam(t_pppParam::bias, t_prn(), lc));
     364          }
     365        }
     366        else {
     367          if (lc.includesCode()) {
     368            if (skip < maxSkip) {
     369              skip += 1;
     370            }
     371            else {
     372              required.push_back(new t_pppParam(t_pppParam::bias, t_prn(), lc));
     373            }
     374          }
     375        }
    586376      }
    587377    }
    … …  
    605395    }
    606396    else {
    607 #ifdef BNC_DEBUG_PPP
    608       //LOG << "push_back  parReq " << parReq->toString() << std::endl;
    609 #endif
    610397      _params.push_back(parReq);
    611398    }
    … …  
    634421////////////////////////////////////////////////////////////////////////////
    635422void t_pppParlist::printResult(const bncTime& epoTime, const SymmetricMatrix& QQ,
    636                                const ColumnVector& xx) const {
     423                               const ColumnVector& xx, double fixRatio) const {
    637424
    638425  string epoTimeStr = string(epoTime);
    … …  
    703490
    704491        << " dU = " << setprecision(4) << neu[2] << " +- "
    705         << setprecision(4) << sqrt(QQneu[2][2])
    706 
    707         << endl;
     492        << setprecision(4) << sqrt(QQneu[2][2]);
     493
     494    if (fixRatio > 0.0) {
     495      LOG << " fix " <<  int(100*fixRatio) << " %";
     496    }
     497    else {
     498      LOG << " flt ";
     499    }
     500
     501    LOG << endl;
    708502  }
    709503  return;

  • trunk/BNC/src/PPP/pppParlist.h

    r10583 r10791  
    1414class t_pppParam {
    1515 public:
    16   enum e_type {crdX, crdY, crdZ, rClkG, rClkR, rClkE, rClkC, trp, ion, amb,
    17                cBiasG1, cBiasR1, cBiasE1, cBiasC1, pBiasG1, pBiasR1, pBiasE1, pBiasC1,
    18                cBiasG2, cBiasR2, cBiasE2, cBiasC2, pBiasG2, pBiasR2, pBiasE2, pBiasC2};
     16  enum e_type {crdX, crdY, crdZ, rClk, trp, ion, amb, bias};
    1917
    20   t_pppParam(e_type type, const t_prn& prn, t_lc::type tLC, const std::vector<t_pppSatObs*>* obsVector = 0);
     18  t_pppParam(e_type type, const t_prn& prn, t_lc LC = t_lc(), const std::vector<t_pppSatObs*>* obsVector = 0);
    2119  ~t_pppParam();
    2220
    2321  e_type type() const {return _type;}
    2422  double x0()  const {return _x0;}
    25   double partial(const bncTime& epoTime, const t_pppSatObs* obs,
    26                  const t_lc::type& tLC) const;
     23  double partial(const bncTime& epoTime, const t_pppSatObs* obs, const t_lc& LC) const;
    2724  bool   epoSpec() const {return _epoSpec;}
    2825  bool   isEqual(const t_pppParam* par2) const {
    29     return (_type == par2->_type && _prn == par2->_prn && _tLC == par2->_tLC);
     26    return (_type == par2->_type && _prn == par2->_prn && _LC == par2->_LC && _sys == par2->_sys);
    3027  }
    3128  void   setIndex(int indexNew) {
    … …  
    3633    _indexOld = -1;
    3734  }
    38   int indexOld() const {return _indexOld;}
    39   int indexNew() const {return _indexNew;}
    40   double sigma0() const {return _sigma0;}
    41   double noise() const {return _noise;}
    42   t_lc::type tLC() const {return _tLC;}
    43   t_prn prn() const {return _prn;}
     35  int         indexOld() const {return _indexOld;}
     36  int         indexNew() const {return _indexNew;}
     37  double      sigma0() const {return _sigma0;}
     38  double      noise() const {return _noise;}
     39  t_lc        LC() const {return _LC;}
     40  t_prn       prn() const {return _prn;}
    4441  std::string toString() const;
    4542
    … …  
    5552  unsigned ambNumEpo() const           {return _ambInfo ? _ambInfo->_numEpo : 0;}
    5653  void     stepAmbNumEpo()             {if (_ambInfo) _ambInfo->_numEpo += 1;}
     54  char     system() const {
     55    if (_prn.valid()) {
     56      return _prn.system();
     57    }
     58    else if (_LC.valid()) {
     59      return _LC.system();
     60    }
     61    else {
     62      return _sys;
     63    }
     64  }
     65  void     setSystem(char sys) {_sys = sys;}
    5766
    5867  static bool sortFunction(const t_pppParam* p1, const t_pppParam* p2) {
    … …  
    6069      return p1->_type < p2->_type;
    6170    }
    62     else if (p1->_tLC != p2->_tLC) {
    63       return p1->_tLC < p2->_tLC;
     71    else if (p1->_LC != p2->_LC) {
     72      return p1->_LC < p2->_LC;
    6473    }
    6574    else if (p1->_prn != p2->_prn) {
    … …  
    8493  e_type       _type;
    8594  t_prn        _prn;
    86   t_lc::type   _tLC;
     95  char         _sys;
     96  t_lc         _LC;
    8797  double       _x0;
    8898  bool         _epoSpec;
    … …  
    104114  const std::vector<t_pppParam*>& params() const {return _params;}
    105115  std::vector<t_pppParam*>& params() {return _params;}
    106   const QMap<char, int>& usedSystems() const {return _usedSystems;}
    107116  void printResult(const bncTime& epoTime, const SymmetricMatrix& QQ,
    108                    const ColumnVector& xx) const;
     117                   const ColumnVector& xx, double fixRatio = 0.0) const;
    109118  void printParams(const bncTime& epoTime);
    110119
    111120 private:
    112121  std::vector<t_pppParam*> _params;
    113   QMap<char, int>          _usedSystems;
     122  std::map<char, int>      _usedSystems;
    114123};
    115124

  • trunk/BNC/src/PPP/pppSatObs.cpp

    r10626 r10791  
    1919#include <iomanip>
    2020#include <cmath>
     21#include <algorithm>
     22#include <set>
    2123#include <newmatio.h>
    2224
    … …  
    4345  _reference  = false;
    4446  _stecSat    = 0.0;
    45   _signalPriorities = QString::fromStdString(OPT->_signalPriorities);
    4647  for (unsigned ii = 0; ii < t_frequency::max; ii++) {
    4748    _obs[ii] = 0;
    … …  
    6061//
    6162////////////////////////////////////////////////////////////////////////////
    62 void t_pppSatObs::prepareObs(const t_satObs& pppSatObs) {
     63bool t_pppSatObs::isBetter(const t_frqObs* aa, t_frqObs* bb, const string& trkModes) const {
     64
     65  if (!trkModes.empty()) {
     66    size_t posA = trkModes.find(aa->trkChar());
     67    size_t posB = trkModes.find(bb->trkChar());
     68    if (posA != posB) {
     69      if      (posA == string::npos) {
     70        return false;
     71      }
     72      else if (posB == string::npos) {
     73        return true;
     74      }
     75      else {
     76        return posA < posB;
     77      }
     78    }
     79  }
     80     
     81  unsigned numValA = 0;
     82  if (aa->_codeValid)  numValA += 1;
     83  if (aa->_phaseValid) numValA += 1;
     84
     85  unsigned numValB = 0;
     86  if (bb->_codeValid)  numValB += 1;
     87  if (bb->_phaseValid) numValB += 1;
     88
     89  if (numValA != numValB) {
     90    return numValA > numValB;
     91  }
     92
     93  return false;
     94}
     95
     96//
     97////////////////////////////////////////////////////////////////////////////
     98void t_pppSatObs::prepareObs(const t_satObs& satObs) {
    6399
    64100  _model.reset();
    65101
    66   std::vector<char> bb = OPT->frqBands(_prn.system());
    67   char frqNum1 = '0';
    68   if (bb.size() >= 1) {
    69     frqNum1 = bb[0];
    70   }
    71   char frqNum2 = '0';
    72   if (bb.size() == 2) {
    73     frqNum2 = bb[1];
    74   }
    75 
    76   // Select pseudo-ranges and phase observations
    77   // -------------------------------------------
    78   QStringList priorList = _signalPriorities.split(" ", Qt::SkipEmptyParts);
    79   string preferredAttrib;
    80   for (unsigned iFreq = 1; iFreq < t_frequency::max; iFreq++) {
    81     t_frequency::type frqType = static_cast<t_frequency::type>(iFreq);
    82     char frqSys = t_frequency::toString(frqType)[0]; //cout << "frqSys: " << frqSys << endl;
    83     char frqNum = t_frequency::toString(frqType)[1]; //cout << "frqNum: " << frqNum << endl;
    84     if (frqSys != _prn.system()) {
    85       continue;
    86     }
    87     if (frqNum != frqNum1 &&
    88         frqNum != frqNum2 ) {
    89       continue;
    90     }
    91     QStringList hlp;
    92     for (int ii = 0; ii < priorList.size(); ii++) {
    93       if (priorList[ii].indexOf(":") != -1) {
    94         hlp = priorList[ii].split(":", Qt::SkipEmptyParts);
    95         if (hlp.size() == 2 && hlp[0].length() == 1 && hlp[0][0] == frqSys) {
    96           hlp = hlp[1].split("&", Qt::SkipEmptyParts);
    97         }
    98         if (hlp.size() == 2 && hlp[0].indexOf(frqNum) != -1) {
    99           preferredAttrib = hlp[1].toStdString(); //cout << "preferredAttrib: " << preferredAttrib << endl;
    100         }
    101       }
    102       for (unsigned iPref = 0; iPref < preferredAttrib.length(); iPref++) {
    103         QString obsType = QString("%1").arg(frqNum) + preferredAttrib[iPref];  //cout << "obstype: " << obsType.toStdString().c_str() << endl;
    104         if (_obs[iFreq] == 0) {
    105           for (unsigned ii = 0; ii < pppSatObs._obs.size(); ii++) {
    106             const t_frqObs* obs = pppSatObs._obs[ii];
    107             //cout << "observation2char: " << obs->_rnxType2ch << " vs. " << obsType.toStdString().c_str()<< endl;
    108             if (obs->_rnxType2ch == obsType.toStdString() &&
    109                 obs->_codeValid  && obs->_code &&
    110                 obs->_phaseValid && obs->_phase) {
    111               _obs[iFreq] = new t_frqObs(*obs); //cout << "================> newObs: " << obs->_rnxType2ch << " obs->_lockTime: " << obs->_lockTime << endl;
    112             }
    113           }
    114         }
    115       }
    116     }
    117   }
    118 
    119   // Used frequency types
    120   // --------------------
    121   _fType1 = t_frqBand::toFreq(_prn.system(), frqNum1);
    122   _fType2 = t_frqBand::toFreq(_prn.system(), frqNum2);
    123 
     102  const t_pppOptions::SysTrkModes* sysTrkModes = OPT->sysTrkModes(_prn.system());
     103  using FrqTrkModes = t_pppOptions::SysTrkModes::FrqTrkModes;
     104 
     105  for (const t_frqObs* obs : satObs._obs) {
     106    if ( (obs->_codeValid || obs->_phaseValid)) {
     107      t_frequency::type frq = t_frequency::toFreq(_prn.system(), obs->frqChar());
     108      if (frq == t_frequency::dummy) {
     109        continue;
     110      }
     111      string trkModes;
     112      if (sysTrkModes) {
     113        const vector<FrqTrkModes>& frqTrkModes = sysTrkModes->_frqTrkModes;
     114        auto it = find_if(frqTrkModes.begin(), frqTrkModes.end(),
     115                          [frq](const FrqTrkModes& mm){return mm._frq == frq;});
     116        if (it != frqTrkModes.end()) {
     117          trkModes = it->_trkModes;
     118        }
     119      }
     120      if (_obs[frq] == 0 || isBetter(obs, _obs[frq], trkModes)) {
     121        delete _obs[frq]; 
     122        _obs[frq] = new t_frqObs(*obs);
     123      }
     124    }
     125  }
     126 
    124127  // Check whether all required frequencies available
    125128  // ------------------------------------------------
    126   for (unsigned ii = 0; ii < OPT->LCs(_prn.system()).size(); ii++) {
    127     t_lc::type tLC = OPT->LCs(_prn.system())[ii];
    128     if (tLC == t_lc::GIM) {continue;}
    129     if (!isValid(tLC)) {
     129  const std::vector<t_lc>& LCs = OPT->LCs(_prn.system());
     130  for (unsigned ii = 0; ii < LCs.size(); ii++) {
     131    t_lc LC = LCs[ii];
     132    if (LC._type == t_lc::GIM) {
     133      continue;
     134    }
     135    if (LC._frq1 != t_frequency::G5 && !isValid(LC)) {
    130136      _valid = false;
    131137      return;
    … …  
    148154  bool totOK  = false;
    149155  ColumnVector satPosOld(6); satPosOld = 0.0;
    150   t_lc::type tLC = t_lc::dummy;
    151   if (isValid(t_lc::cIF)) {
    152     tLC = t_lc::cIF;
    153   }
    154   if (tLC == t_lc::dummy && isValid(t_lc::c1)) {
    155       tLC = t_lc::c1;
    156   }
    157   if (tLC == t_lc::dummy && isValid(t_lc::c2)) {
    158       tLC = t_lc::c2;
    159   }
    160   if (tLC == t_lc::dummy) {
    161     _valid = false;
     156
     157  _valid = false;
     158  t_frequency::type frq1 = t_frequency::dummy;
     159  t_frequency::type frq2 = t_frequency::dummy;
     160  OPT->defaultFrqs(_prn.system(), frq1, frq2);
     161  if (frq1 != t_frequency::dummy) {
     162    t_lc lc1(t_lc::code, frq1);
     163    if (isValid(lc1)) {
     164      _valid = true;
     165      _rangeLC = lc1;
     166    }
     167    if (frq2 != t_frequency::dummy) {
     168      t_lc lcIF(t_lc::codeIF, frq1, frq2);
     169      if (isValid(lcIF)) {
     170        _valid = true;
     171        _rangeLC = lcIF;
     172      }
     173    }
     174  }
     175  if (!_valid) {
    162176    return;
    163177  }
    164   double prange = obsValue(tLC);
     178 
     179  double prange = obsValue(_rangeLC);
    165180  for (int ii = 1; ii <= 10; ii++) {
    166181    bncTime ToT = _time - prange / t_CST::c - _xcSat[3];
    … …  
    188203//
    189204////////////////////////////////////////////////////////////////////////////
    190 void t_pppSatObs::lcCoeff(t_lc::type tLC,
     205void t_pppSatObs::lcCoeff(t_lc LC,
    191206                          map<t_frequency::type, double>& codeCoeff,
    192207                          map<t_frequency::type, double>& phaseCoeff,
    … …  
    197212  ionoCoeff.clear();
    198213
    199   double f1 = t_CST::freq(_fType1, _channel);
    200   double f2 = t_CST::freq(_fType2, _channel);
     214  double f1 = t_CST::freq(LC._frq1, _channel);
     215  double f2 = t_CST::freq(LC._frq2, _channel);
    201216  double f1GPS = t_CST::freq(t_frequency::G1, 0);
    202217
    203   switch (tLC) {
    204   case t_lc::l1:
    205     phaseCoeff[_fType1] =  1.0;
    206     ionoCoeff [_fType1] = -1.0 * pow(f1GPS, 2) / pow(f1, 2);
     218  switch (LC._type) {
     219  case t_lc::phase:
     220    phaseCoeff[LC._frq1] =  1.0;
     221    ionoCoeff [LC._frq1] = -1.0 * pow(f1GPS, 2) / pow(f1, 2);
    207222    return;
    208   case t_lc::l2:
    209     phaseCoeff[_fType2] = 1.0;
    210     ionoCoeff [_fType2] = -1.0 * pow(f1GPS, 2) / pow(f2, 2);
     223  case t_lc::code:
     224    codeCoeff[LC._frq1] = 1.0;
     225    ionoCoeff[LC._frq1] = pow(f1GPS, 2) / pow(f1, 2);
    211226    return;
    212   case t_lc::lIF:
    213     phaseCoeff[_fType1] =  f1 * f1 / (f1 * f1 - f2 * f2);
    214     phaseCoeff[_fType2] = -f2 * f2 / (f1 * f1 - f2 * f2);
     227  case t_lc::phaseIF:
     228    phaseCoeff[LC._frq1] =  f1 * f1 / (f1 * f1 - f2 * f2);
     229    phaseCoeff[LC._frq2] = -f2 * f2 / (f1 * f1 - f2 * f2);
     230    return;
     231  case t_lc::codeIF:
     232    codeCoeff[LC._frq1] =  f1 * f1 / (f1 * f1 - f2 * f2);
     233    codeCoeff[LC._frq2] = -f2 * f2 / (f1 * f1 - f2 * f2);
    215234    return;
    216235  case t_lc::MW:
    217     phaseCoeff[_fType1] =  f1 / (f1 - f2);
    218     phaseCoeff[_fType2] = -f2 / (f1 - f2);
    219     codeCoeff[_fType1] = -f1 / (f1 + f2);
    220     codeCoeff[_fType2] = -f2 / (f1 + f2);
     236    phaseCoeff[LC._frq1] =  f1 / (f1 - f2);
     237    phaseCoeff[LC._frq2] = -f2 / (f1 - f2);
     238    codeCoeff [LC._frq1] = -f1 / (f1 + f2);
     239    codeCoeff [LC._frq2] = -f2 / (f1 + f2);
    221240    return;
    222241  case t_lc::CL:
    223     phaseCoeff[_fType1] =  0.5;
    224     codeCoeff [_fType1] =  0.5;
    225     return;
    226   case t_lc::c1:
    227     codeCoeff[_fType1] = 1.0;
    228     ionoCoeff[_fType1] = pow(f1GPS, 2) / pow(f1, 2);
    229     return;
    230   case t_lc::c2:
    231     codeCoeff[_fType2] = 1.0;
    232     ionoCoeff[_fType2] = pow(f1GPS, 2) / pow(f2, 2);
    233     return;
    234   case t_lc::cIF:
    235     codeCoeff[_fType1] =  f1 * f1 / (f1 * f1 - f2 * f2);
    236     codeCoeff[_fType2] = -f2 * f2 / (f1 * f1 - f2 * f2);
     242    phaseCoeff[LC._frq1] =  0.5;
     243    codeCoeff [LC._frq1] =  0.5;
    237244    return;
    238245  case t_lc::GIM:
    … …  
    245252//
    246253////////////////////////////////////////////////////////////////////////////
    247 bool t_pppSatObs::isValid(t_lc::type tLC) const {
     254bool t_pppSatObs::isValid(t_lc LC) const {
    248255  bool valid = true;
    249   obsValue(tLC, &valid);
     256  obsValue(LC, &valid);
    250257
    251258  return valid;
    … …  
    253260//
    254261////////////////////////////////////////////////////////////////////////////
    255 double t_pppSatObs::obsValue(t_lc::type tLC, bool* valid) const {
     262double t_pppSatObs::obsValue(t_lc LC, bool* valid) const {
    256263
    257264  double retVal = 0.0;
    … …  
    259266
    260267  // Pseudo observations
    261   if (tLC == t_lc::GIM) {
     268  if (LC._type == t_lc::GIM) {
    262269    if (_stecSat == 0.0) {
    263270      if (valid) *valid = false;
    … …  
    272279  map<t_frequency::type, double> phaseCoeff;
    273280  map<t_frequency::type, double> ionoCoeff;
    274   lcCoeff(tLC, codeCoeff, phaseCoeff, ionoCoeff);
     281  lcCoeff(LC, codeCoeff, phaseCoeff, ionoCoeff);
    275282
    276283  map<t_frequency::type, double>::const_iterator it;
    … …  
    303310//
    304311////////////////////////////////////////////////////////////////////////////
    305 double t_pppSatObs::lambda(t_lc::type tLC) const {
    306 
    307   double f1 = t_CST::freq(_fType1, _channel);
    308   double f2 = t_CST::freq(_fType2, _channel);
    309 
    310   if      (tLC == t_lc::l1) {
     312double t_pppSatObs::lambda(t_lc LC) const {
     313
     314  double f1 = t_CST::freq(LC._frq1, _channel);
     315  double f2 = t_CST::freq(LC._frq2, _channel);
     316
     317  if      (LC._type == t_lc::phase) {
    311318    return t_CST::c / f1;
    312319  }
    313   else if (tLC == t_lc::l2) {
    314     return t_CST::c / f2;
    315   }
    316   else if (tLC == t_lc::lIF) {
     320  else if (LC._type == t_lc::phaseIF) {
    317321    return t_CST::c / (f1 + f2);
    318322  }
    319   else if (tLC == t_lc::MW) {
     323  else if (LC._type == t_lc::MW) {
    320324    return t_CST::c / (f1 - f2);
    321325  }
    322   else if (tLC == t_lc::CL) {
     326  else if (LC._type == t_lc::CL) {
    323327    return t_CST::c / f1 / 2.0;
    324328  }
    … …  
    329333//
    330334////////////////////////////////////////////////////////////////////////////
    331 double t_pppSatObs::sigma(t_lc::type tLC) const {
     335double t_pppSatObs::sigma(t_lc LC) const {
    332336
    333337  double retVal = 0.0;
    … …  
    335339  map<t_frequency::type, double> phaseCoeff;
    336340  map<t_frequency::type, double> ionoCoeff;
    337   lcCoeff(tLC, codeCoeff, phaseCoeff, ionoCoeff);
    338 
    339   if (tLC == t_lc::GIM) {
     341  lcCoeff(LC, codeCoeff, phaseCoeff, ionoCoeff);
     342
     343  if (LC._type == t_lc::GIM) {
    340344    retVal = OPT->_sigmaGIM * OPT->_sigmaGIM;
    341345  }
    … …  
    354358  // De-Weight R
    355359  // -----------
    356   if (_prn.system() == 'R'&& t_lc::includesCode(tLC)) {
     360  if (_prn.system() == 'R'&& LC.includesCode()) {
    357361    retVal *= 5.0;
    358362  }
    … …  
    361365  // -----------------------------
    362366  double cEle = 1.0;
    363   if ( (OPT->_eleWgtCode  && t_lc::includesCode(tLC)) ||
    364        (OPT->_eleWgtPhase && t_lc::includesPhase(tLC)) ) {
     367  if ( (OPT->_eleWgtCode  && LC.includesCode()) ||
     368       (OPT->_eleWgtPhase && LC.includesPhase()) ) {
    365369    double eleD = eleSat()*180.0/M_PI;
    366370    double hlp  = fabs(90.0 - eleD);
    … …  
    373377//
    374378////////////////////////////////////////////////////////////////////////////
    375 double t_pppSatObs::maxRes(t_lc::type tLC) const {
     379double t_pppSatObs::maxRes(t_lc LC) const {
    376380  double retVal = 0.0;
    377381
    … …  
    379383  map<t_frequency::type, double> phaseCoeff;
    380384  map<t_frequency::type, double> ionoCoeff;
    381   lcCoeff(tLC, codeCoeff, phaseCoeff, ionoCoeff);
     385  lcCoeff(LC, codeCoeff, phaseCoeff, ionoCoeff);
    382386
    383387  map<t_frequency::type, double>::const_iterator it;
    … …  
    388392    retVal += it->second * it->second * OPT->_maxResL1 * OPT->_maxResL1;
    389393  }
    390   if (tLC == t_lc::GIM) {
     394  if (LC._type == t_lc::GIM) {
    391395    retVal = OPT->_maxResGIM * OPT->_maxResGIM + OPT->_maxResGIM * OPT->_maxResGIM;
    392396  }
    … …  
    517521        }
    518522        const t_frqObs* obs = _obs[iFreq];
    519         if (obs &&
    520             obs->_rnxType2ch == bias._rnxType2ch) {
    521           _model._codeBias[iFreq]  = bias._value;
     523        if (obs && obs->_rnxType2ch == bias._rnxType2ch) {
     524          _model._codeBias[iFreq] = (bias._value != 0.0 ? bias._value : ZEROVALUE);
    522525        }
    523526      }
    … …  
    527530  // Phase Biases
    528531  // -----------
    529   const t_satPhaseBias* satPhaseBias = PPP_CLIENT->obsPool()->satPhaseBias(_prn);
    530   double yaw = 0.0;
    531   bool ssr = false;
    532   if (satPhaseBias) {
    533     double dt = station->epochTime() - satPhaseBias->_time;
    534     if (satPhaseBias->_updateInt) {
    535       dt -= (0.5 * ssrUpdateInt[satPhaseBias->_updateInt]);
    536     }
    537     yaw = satPhaseBias->_yaw + satPhaseBias->_yawRate * dt;
    538     ssr = true;
    539     for (unsigned ii = 0; ii < satPhaseBias->_bias.size(); ii++) {
    540       const t_frqPhaseBias& bias = satPhaseBias->_bias[ii];
    541       for (unsigned iFreq = 1; iFreq < t_frequency::max; iFreq++) {
    542         string frqStr = t_frequency::toString(t_frequency::type(iFreq));
    543         if (frqStr[0] != _prn.system()) {
    544           continue;
    545         }
    546         const t_frqObs* obs = _obs[iFreq];
    547         if (obs &&
    548             obs->_rnxType2ch == bias._rnxType2ch) {
    549           _model._phaseBias[iFreq]  = bias._value;
     532  double yaw    = 0.0;
     533  bool   useYaw = false;
     534  if (OPT->arSystem(_prn.system())) {
     535    const t_satPhaseBias* satPhaseBias = PPP_CLIENT->obsPool()->satPhaseBias(_prn);
     536    if (satPhaseBias) {
     537      if (OPT->_ar._useYaw) {
     538        double dt = station->epochTime() - satPhaseBias->_time;
     539        if (satPhaseBias->_updateInt) {
     540          dt -= (0.5 * ssrUpdateInt[satPhaseBias->_updateInt]);
     541        }
     542        yaw = satPhaseBias->_yaw + satPhaseBias->_yawRate * dt;
     543        useYaw = true;
     544      }
     545      for (unsigned ii = 0; ii < satPhaseBias->_bias.size(); ii++) {
     546        const t_frqPhaseBias& bias = satPhaseBias->_bias[ii];
     547        if (bias._fixIndicator) {  // if AR, biases without fixIndicator not used
     548          for (unsigned iFreq = 1; iFreq < t_frequency::max; iFreq++) {
     549            string frqStr = t_frequency::toString(t_frequency::type(iFreq));
     550            if (frqStr[0] != _prn.system()) {
     551              continue;
     552            }
     553            t_frqObs* obs = _obs[iFreq];
     554            if (obs && obs->_rnxType2ch[0] == bias._rnxType2ch[0]) { // allow different tracking mode
     555              _model._phaseBias[iFreq] = (bias._value != 0.0 ? bias._value : ZEROVALUE);
     556              obs->_biasJumpCounter = bias._jumpCounter;
     557            }
     558          }
    550559        }
    551560      }
    … …  
    555564  // Phase Wind-Up
    556565  // -------------
    557   _model._windUp = station->windUp(_time, _prn, rSat, ssr, yaw, vSat) ;
     566  _model._windUp = station->windUp(_time, _prn, rSat, useYaw, yaw, vSat) ;
    558567
    559568  // Relativistic effect due to earth gravity
    … …  
    574583  bool vTecUsage = true;
    575584  for (unsigned ii = 0; ii < OPT->LCs(_prn.system()).size(); ii++) {
    576     t_lc::type tLC = OPT->LCs(_prn.system())[ii];
    577     if (tLC == t_lc::cIF || tLC == t_lc::lIF) {
     585    t_lc LC = OPT->LCs(_prn.system())[ii];
     586    if (LC._type == t_lc::codeIF || LC._type == t_lc::phaseIF) {
    578587      vTecUsage = false;
    579588    }
    … …  
    601610  _model._set = true;
    602611
    603   //printModel();
     612  if (OPT->_logMode == t_pppOptions::all) {
     613    printModel();
     614  }
    604615
    605616  return success;
    … …  
    636647      if (_prn.system() == frqStr[0]) {
    637648      LOG << "PCO           : " << frqStr << setw(12) << setprecision(3) << _model._antPCO[iFreq]       << endl
    638           << "BIAS CODE     : " << frqStr << setw(12) << setprecision(3) << _model._codeBias[iFreq]     << "\t(" << _obs[iFreq]->_rnxType2ch[1] << ") " << endl
    639           << "BIAS PHASE    : " << frqStr << setw(12) << setprecision(3) << _model._phaseBias[iFreq]    << "\t(" << _obs[iFreq]->_rnxType2ch[1] << ") " << endl
     649          << "BIAS CODE     : " << frqStr << setw(12) << setprecision(3) << _model._codeBias[iFreq]     << "\t(" << _obs[iFreq]->trkChar() << ") " << endl
     650          << "BIAS PHASE    : " << frqStr << setw(12) << setprecision(3) << _model._phaseBias[iFreq]    << "\t(" << _obs[iFreq]->trkChar() << ") " << endl
    640651          << "IONO CODEDELAY: " << frqStr << setw(12) << setprecision(3) << _model._ionoCodeDelay[iFreq]<< endl;
    641652      }
    … …  
    652663  char sys = _prn.system();
    653664  for (unsigned ii = 0; ii < OPT->LCs(sys).size(); ii++) {
    654     t_lc::type tLC = OPT->LCs(sys)[ii];
    655     LOG << "OBS-CMP " << setw(4) << t_lc::toString(tLC) << ": " << _prn.toString() << " "
    656         << setw(12) << setprecision(3) << obsValue(tLC) << " "
    657         << setw(12) << setprecision(3) << cmpValue(tLC) << " "
    658         << setw(12) << setprecision(3) << obsValue(tLC)  - cmpValue(tLC) << endl;
    659   }
    660 }
    661 
    662 //
    663 ////////////////////////////////////////////////////////////////////////////
    664 double t_pppSatObs::cmpValueForBanc(t_lc::type tLC) const {
    665   return cmpValue(tLC) - _model._rho - _model._sagnac - _model._recClkM;
    666 }
    667 
    668 //
    669 ////////////////////////////////////////////////////////////////////////////
    670 double t_pppSatObs::cmpValue(t_lc::type tLC) const {
     665    t_lc LC = OPT->LCs(sys)[ii];
     666    LOG << "OBS-CMP " << setw(4) << LC.toString() << ": " << _prn.toString() << " "
     667        << setw(12) << setprecision(3) << obsValue(LC) << " "
     668        << setw(12) << setprecision(3) << cmpValue(LC) << " "
     669        << setw(12) << setprecision(3) << obsValue(LC)  - cmpValue(LC) << endl;
     670  }
     671}
     672
     673//
     674////////////////////////////////////////////////////////////////////////////
     675double t_pppSatObs::cmpValueForBanc(t_lc LC) const {
     676  return cmpValue(LC) - _model._rho - _model._sagnac - _model._recClkM;
     677}
     678
     679//
     680////////////////////////////////////////////////////////////////////////////
     681double t_pppSatObs::cmpValue(t_lc LC) const {
    671682  double cmpValue;
    672683
    673   if      (!isValid(tLC)) {
     684  if      (!isValid(LC)) {
    674685    cmpValue =  0.0;
    675686  }
    676   else if (tLC == t_lc::GIM) {
     687  else if (LC._type == t_lc::GIM) {
    677688    cmpValue =  _stecSat;
    678689  }
    … …  
    691702    map<t_frequency::type, double> phaseCoeff;
    692703    map<t_frequency::type, double> ionoCoeff;
    693     lcCoeff(tLC, codeCoeff, phaseCoeff, ionoCoeff);
     704    lcCoeff(LC, codeCoeff, phaseCoeff, ionoCoeff);
    694705    map<t_frequency::type, double>::const_iterator it;
    695706    for (it = codeCoeff.begin(); it != codeCoeff.end(); it++) {
    696707      t_frequency::type tFreq = it->first;
    697708      dispPart += it->second * (_model._antPCO[tFreq] - _model._codeBias[tFreq]);
    698       if (OPT->PPP_RTK) {
    699         dispPart += it->second * (_model._ionoCodeDelay[tFreq]);
    700       }
    701709    }
    702710    for (it = phaseCoeff.begin(); it != phaseCoeff.end(); it++) {
    … …  
    704712      dispPart += it->second * (_model._antPCO[tFreq] - _model._phaseBias[tFreq] +
    705713                                _model._windUp * t_CST::lambda(tFreq, _channel));
    706       if (OPT->PPP_RTK) {
    707         dispPart += it->second * (- _model._ionoCodeDelay[tFreq]);
    708       }
    709714    }
    710715    cmpValue = nonDisp + dispPart;
    … …  
    716721//
    717722////////////////////////////////////////////////////////////////////////////
    718 void t_pppSatObs::setRes(t_lc::type tLC, double res) {
    719   _res[tLC] = res;
    720 }
    721 
    722 //
    723 ////////////////////////////////////////////////////////////////////////////
    724 double t_pppSatObs::getRes(t_lc::type tLC) const {
    725   map<t_lc::type, double>::const_iterator it = _res.find(tLC);
     723void t_pppSatObs::setRes(t_lc LC, double res) {
     724  _res[LC] = res;
     725}
     726
     727//
     728////////////////////////////////////////////////////////////////////////////
     729double t_pppSatObs::getRes(t_lc LC) const {
     730  map<t_lc, double>::const_iterator it = _res.find(LC);
    726731  if (it != _res.end()) {
    727732    return it->second;
    … …  
    742747  return pseudoObsIono;
    743748}
     749
     750//
     751////////////////////////////////////////////////////////////////////////////
     752bool t_pppSatObs::hasBiases() const {
     753  bool ar = OPT->arSystem(_prn.system());
     754  set<t_frequency::type> frqs;
     755  for (const auto& lc : OPT->LCs(_prn.system())) {
     756    if (lc._frq1 != t_frequency::dummy) frqs.insert(lc._frq1);
     757    if (lc._frq2 != t_frequency::dummy) frqs.insert(lc._frq2);
     758  }
     759  for (int iFreq : frqs) {
     760    if (_obs[iFreq] != 0) {
     761      if (_model._codeBias[iFreq] == 0 || (ar && _model._phaseBias[iFreq] == 0)) {
     762        return false;
     763      }
     764    }
     765  }
     766  return true;
     767}

  • trunk/BNC/src/PPP/pppSatObs.h

    r10251 r10791  
    1111namespace BNC_PPP {
    1212
    13 class t_pppStation;
     13  class t_pppStation;
    1414
    15 class t_pppSatObs {
    16  public:
    17   t_pppSatObs(const t_satObs& satObs);
    18   ~t_pppSatObs();
    19   bool                isValid() const {return _valid;};
    20   bool                isValid(t_lc::type tLC) const;
    21   bool                isReference() const {return _reference;};
    22   void                setAsReference() {_reference = true;};
    23   void                resetReference() {_reference = false;};
    24   const t_prn&        prn() const {return _prn;}
    25   const ColumnVector& xc() const {return _xcSat;}
    26   const bncTime&      time() const {return _time;}
    27   t_irc               cmpModel(const t_pppStation* station);
    28   double              obsValue(t_lc::type tLC, bool* valid = 0) const;
    29   double              cmpValue(t_lc::type tLC) const;
    30   double              cmpValueForBanc(t_lc::type tLC) const;
    31   double              rho() const {return _model._rho;}
    32   double              sagnac() const {return _model._sagnac;}
    33   double              eleSat() const {return _model._eleSat;}
    34   bool                modelSet() const {return _model._set;}
    35   void                printModel() const;
    36   void                printObsMinusComputed() const;
    37   void                lcCoeff(t_lc::type tLC,
    38                               std::map<t_frequency::type, double>& codeCoeff,
    39                               std::map<t_frequency::type, double>& phaseCoeff,
    40                               std::map<t_frequency::type, double>& ionoCoeff) const;
    41   double              lambda(t_lc::type tLC) const;
    42   double              sigma(t_lc::type tLC) const;
    43   double              maxRes(t_lc::type tLC) const;
    44   bool                outlier() const {return _outlier;}
    45   void                setOutlier() {_outlier = true;}
    46   void                resetOutlier() {_outlier = false;}
    47   void                setRes(t_lc::type tLC, double res);
    48   double              getRes(t_lc::type tLC) const;
    49   bool                setPseudoObsIono(t_frequency::type freq);
    50   double              getIonoCodeDelay(t_frequency::type freq) {return _model._ionoCodeDelay[freq];}
    51   double              getCodeBias(t_frequency::type freq) {return _model._codeBias[freq];}
    52   t_frequency::type   fType1() const {return _fType1;}
    53   t_frequency::type   fType2() const {return _fType2;}
     15  class t_pppSatObs {
     16  public:
     17    t_pppSatObs(const t_satObs& satObs);
     18    ~t_pppSatObs();
     19    bool                isValid() const {return _valid;};
     20    bool                isValid(t_lc LC) const;
     21    t_lc                rangeLC() const {return _rangeLC;};
     22    bool                isReference() const {return _reference;};
     23    void                setAsReference() {_reference = true;};
     24    void                resetReference() {_reference = false;};
     25    const t_prn&        prn() const {return _prn;}
     26    const ColumnVector& xc() const {return _xcSat;}
     27    const bncTime&      time() const {return _time;}
     28    t_irc               cmpModel(const t_pppStation* station);
     29    double              obsValue(t_lc LC, bool* valid = 0) const;
     30    double              cmpValue(t_lc LC) const;
     31    double              cmpValueForBanc(t_lc LC) const;
     32    double              rho() const {return _model._rho;}
     33    double              sagnac() const {return _model._sagnac;}
     34    double              eleSat() const {return _model._eleSat;}
     35    bool                modelSet() const {return _model._set;}
     36    void                printModel() const;
     37    void                printObsMinusComputed() const;
     38    void                lcCoeff(t_lc LC,
     39                                std::map<t_frequency::type, double>& codeCoeff,
     40                                std::map<t_frequency::type, double>& phaseCoeff,
     41                                std::map<t_frequency::type, double>& ionoCoeff) const;
     42    double              lambda(t_lc LC) const;
     43    double              sigma(t_lc LC) const;
     44    double              maxRes(t_lc LC) const;
     45    bool                outlier() const {return _outlier;}
     46    void                setOutlier() {_outlier = true;}
     47    void                resetOutlier() {_outlier = false;}
     48    void                setRes(t_lc LC, double res);
     49    double              getRes(t_lc LC) const;
     50    bool                setPseudoObsIono(t_frequency::type freq);
     51    double              getIonoCodeDelay(t_frequency::type freq) {return _model._ionoCodeDelay[freq];}
     52    double              getCodeBias(t_frequency::type freq) {return _model._codeBias[freq];}
     53    bool                hasBiases() const;
    5454
    55   // RINEX
    56   bool slip() const {
    57     for (unsigned ii = 1; ii < t_frequency::max; ii++) {
    58       if (_obs[ii] && _obs[ii]->_slip) {
    59         return true;
     55    // RINEX
     56    bool slip() const {
     57      for (unsigned ii = 1; ii < t_frequency::max; ii++) {
     58        if (_obs[ii] && _obs[ii]->_slip) {
     59          return true;
     60        }
    6061      }
     62      return false;
    6163    }
    62     return false;
    63   }
    6464
    65   // RTCM
    66   int slipCounter() const {
    67     int cnt = -1;
    68     for (unsigned ii = 1; ii < t_frequency::max; ii++) {
    69       if (_obs[ii] && _obs[ii]->_slipCounter > cnt) {
    70         cnt = _obs[ii]->_slipCounter;
     65    // RTCM
     66    int slipCounter() const {
     67      int cnt = -1;
     68      for (unsigned ii = 1; ii < t_frequency::max; ii++) {
     69        if (_obs[ii] && _obs[ii]->_slipCounter > cnt) {
     70          cnt = _obs[ii]->_slipCounter;
     71        }
    7172      }
     73      return cnt;
    7274    }
    73     return cnt;
    74   }
    7575
    76   int biasJumpCounter() const {
    77     int jmp = -1;
    78     for (unsigned ii = 1; ii < t_frequency::max; ii++) {
    79       if (_obs[ii] && _obs[ii]->_biasJumpCounter > jmp) {
    80         jmp = _obs[ii]->_biasJumpCounter;
     76    int biasJumpCounter() const {
     77      int jmp = -1;
     78      for (unsigned ii = 1; ii < t_frequency::max; ii++) {
     79        if (_obs[ii] && _obs[ii]->_biasJumpCounter > jmp) {
     80          jmp = _obs[ii]->_biasJumpCounter;
     81        }
    8182      }
     83      return jmp;
    8284    }
    83     return jmp;
    84   }
    8585
    86  private:
    87   class t_model {
    88    public:
    89     t_model() {reset();}
    90     ~t_model() {}
    91     void reset() {
    92       _set       = false;
    93       _rho       = 0.0;
    94       _eleSat    = 0.0;
    95       _azSat     = 0.0;
    96       _elTx      = 0.0;
    97       _azTx      = 0.0;
    98       _recClkM   = 0.0;
    99       _satClkM   = 0.0;
    100       _sagnac    = 0.0;
    101       _antEcc    = 0.0;
    102       _tropo     = 0.0;
    103       _tropo0    = 0.0;
    104       _tideEarth = 0.0;
    105       _tideOcean = 0.0;
    106       _windUp    = 0.0;
    107       _rel       = 0.0;
    108       for (unsigned ii = 0; ii < t_frequency::max; ii++) {
    109         _antPCO[ii]        = 0.0;
    110         _codeBias[ii]      = 0.0;
    111         _phaseBias[ii]     = 0.0;
    112         _ionoCodeDelay[ii] = 0.0;
     86  private:
     87    class t_model {
     88    public:
     89      t_model() {reset();}
     90      ~t_model() {}
     91      void reset() {
     92        _set       = false;
     93        _rho       = 0.0;
     94        _eleSat    = 0.0;
     95        _azSat     = 0.0;
     96        _elTx      = 0.0;
     97        _azTx      = 0.0;
     98        _recClkM   = 0.0;
     99        _satClkM   = 0.0;
     100        _sagnac    = 0.0;
     101        _antEcc    = 0.0;
     102        _tropo     = 0.0;
     103        _tropo0    = 0.0;
     104        _tideEarth = 0.0;
     105        _tideOcean = 0.0;
     106        _windUp    = 0.0;
     107        _rel       = 0.0;
     108        for (unsigned ii = 0; ii < t_frequency::max; ii++) {
     109          _antPCO[ii]        = 0.0;
     110          _codeBias[ii]      = 0.0;
     111          _phaseBias[ii]     = 0.0;
     112          _ionoCodeDelay[ii] = 0.0;
     113        }
    113114      }
    114     }
    115     bool   _set;
    116     double _rho;
    117     double _eleSat;
    118     double _azSat;
    119     double _elTx;
    120     double _azTx;
    121     double _recClkM;
    122     double _satClkM;
    123     double _sagnac;
    124     double _antEcc;
    125     double _tropo;
    126     double _tropo0;
    127     double _tideEarth;
    128     double _tideOcean;
    129     double _windUp;
    130     double _rel;
    131     double _antPCO[t_frequency::max];
    132     double _codeBias[t_frequency::max];
    133     double _phaseBias[t_frequency::max];
    134     double _ionoCodeDelay[t_frequency::max];
     115      bool   _set;
     116      double _rho;
     117      double _eleSat;
     118      double _azSat;
     119      double _elTx;
     120      double _azTx;
     121      double _recClkM;
     122      double _satClkM;
     123      double _sagnac;
     124      double _antEcc;
     125      double _tropo;
     126      double _tropo0;
     127      double _tideEarth;
     128      double _tideOcean;
     129      double _windUp;
     130      double _rel;
     131      double _antPCO[t_frequency::max];
     132      double _codeBias[t_frequency::max];
     133      double _phaseBias[t_frequency::max];
     134      double _ionoCodeDelay[t_frequency::max];
     135    };
     136
     137    void prepareObs(const t_satObs& satObs);
     138
     139    bool isBetter(const t_frqObs* aa, t_frqObs* bb, const std::string& trkModes) const;
     140
     141    bool                   _valid;
     142    bool                   _reference;
     143    t_lc                   _rangeLC;
     144    t_prn                  _prn;
     145    bncTime                _time;
     146    int                    _channel;
     147    t_frqObs*              _obs[t_frequency::max];
     148    ColumnVector           _xcSat;
     149    ColumnVector           _vvSat;
     150    t_model                _model;
     151    bool                   _outlier;
     152    std::map<t_lc, double> _res;
     153    double                 _signalPropagationTime;
     154    double                 _stecSat;
     155    double                 _tropo0;
    135156  };
    136 
    137   void prepareObs(const t_satObs& satObs);
    138 
    139   bool                         _valid;
    140   bool                         _reference;
    141   t_frequency::type            _fType1;
    142   t_frequency::type            _fType2;
    143   t_prn                        _prn;
    144   bncTime                      _time;
    145   int                          _channel;
    146   t_frqObs*                    _obs[t_frequency::max];
    147   ColumnVector                 _xcSat;
    148   ColumnVector                 _vvSat;
    149   t_model                      _model;
    150   bool                         _outlier;
    151   std::map<t_lc::type, double> _res;
    152   double                       _signalPropagationTime;
    153   double                       _stecSat;
    154   double                       _tropo0;
    155   QString                      _signalPriorities;
    156         };
    157157
    158158}
Note: See TracChangeset for help on using the changeset viewer.