Changeset 3405 in ntrip for trunk/BNC


Ignore:
Timestamp:
Sep 2, 2011, 7:00:35 PM (13 years ago)
Author:
mervart
Message:
 
Location:
trunk/BNC
Files:
7 edited

Legend:

Unmodified
Added
Removed

  • trunk/BNC/bncmapview.cpp

    r3395 r3405  
    9191
    9292// -------------
    93 void BncMapView::mouseReleaseEvent(QMouseEvent* /* event */)
     93void BncMapView::mouseReleaseEvent(QMouseEvent* event)
    9494{   
    9595   setCursor(Qt::OpenHandCursor);

  • trunk/BNC/bncmodel.cpp

    r3404 r3405  
    5252#include "bnctides.h"
    5353#include "bncantex.h"
    54 #include "bnccomb.h"
    5554
    5655using namespace std;
    5756
    58 const unsigned MINOBS                = 5;
    59 const double   MINELE                = 10.0 * M_PI / 180.0;
    60 const double   MAXRES_CODE           = 10.0;
    61 const double   MAXRES_PHASE          = 0.02;
    62 const double   GLONASS_WEIGHT_FACTOR = 5.0;
     57const unsigned MINOBS           =    5;
     58const double   MINELE_GPS       = 10.0 * M_PI / 180.0;
     59const double   MINELE_GLO       = 10.0 * M_PI / 180.0;
     60const double   MINELE_GAL       = 10.0 * M_PI / 180.0;
     61const double   MAXRES_CODE_GPS  = 10.0;
     62const double   MAXRES_PHASE_GPS = 0.05;
     63const double   MAXRES_PHASE_GLO = 0.05;
     64const double   MAXRES_CODE_GAL  = 10.0;
     65const double   MAXRES_PHASE_GAL = 0.05;
    6366
    6467// Constructor
    … …  
    254257  _xcBanc.ReSize(4);  _xcBanc  = 0.0;
    255258  _ellBanc.ReSize(3); _ellBanc = 0.0;
    256 
    257   // Save copy of data (used in outlier detection)
    258   // ---------------------------------------------
    259   _epoData_sav = new t_epoData();
    260259}
    261260
    … …  
    272271    delete _params[iPar-1];
    273272  }
    274   delete _epoData_sav;
    275273}
    276274
    … …  
    324322  Tracer tracer("bncModel::cmpBancroft");
    325323
    326   if (epoData->sizeSys('G') < MINOBS) {
     324  if (epoData->sizeGPS() < MINOBS) {
    327325    _log += "bncModel::cmpBancroft: not enough data\n";
    328326    return failure;
    329327  }
    330328
    331   Matrix BB(epoData->sizeSys('G'), 4);
    332 
    333   QMapIterator<QString, t_satData*> it(epoData->satData);
     329  Matrix BB(epoData->sizeGPS(), 4);
     330
     331  QMapIterator<QString, t_satData*> it(epoData->satDataGPS);
    334332  int iObsBanc = 0;
    335333  while (it.hasNext()) {
     334    ++iObsBanc;
    336335    it.next();
     336    QString    prn     = it.key();
    337337    t_satData* satData = it.value();
    338     if (satData->system() == 'G') {
    339       ++iObsBanc;
    340       QString    prn     = it.key();
    341       BB(iObsBanc, 1) = satData->xx(1);
    342       BB(iObsBanc, 2) = satData->xx(2);
    343       BB(iObsBanc, 3) = satData->xx(3);
    344       BB(iObsBanc, 4) = satData->P3 + satData->clk;
    345     }
     338    BB(iObsBanc, 1) = satData->xx(1);
     339    BB(iObsBanc, 2) = satData->xx(2);
     340    BB(iObsBanc, 3) = satData->xx(3);
     341    BB(iObsBanc, 4) = satData->P3 + satData->clk;
    346342  }
    347343
    … …  
    354350  // Compute Satellite Elevations
    355351  // ----------------------------
    356   QMutableMapIterator<QString, t_satData*> im(epoData->satData);
    357   while (im.hasNext()) {
    358     im.next();
    359     t_satData* satData = im.value();
     352  QMutableMapIterator<QString, t_satData*> iGPS(epoData->satDataGPS);
     353  while (iGPS.hasNext()) {
     354    iGPS.next();
     355    t_satData* satData = iGPS.value();
    360356    cmpEle(satData);
    361     if (satData->eleSat < MINELE) {
     357    if (satData->eleSat < MINELE_GPS) {
    362358      delete satData;
    363       im.remove();
     359      iGPS.remove();
     360    }
     361  }
     362
     363  QMutableMapIterator<QString, t_satData*> iGlo(epoData->satDataGlo);
     364  while (iGlo.hasNext()) {
     365    iGlo.next();
     366    t_satData* satData = iGlo.value();
     367    cmpEle(satData);
     368    if (satData->eleSat < MINELE_GLO) {
     369      delete satData;
     370      iGlo.remove();
     371    }
     372  }
     373
     374  QMutableMapIterator<QString, t_satData*> iGal(epoData->satDataGal);
     375  while (iGal.hasNext()) {
     376    iGal.next();
     377    t_satData* satData = iGal.value();
     378    cmpEle(satData);
     379    if (satData->eleSat < MINELE_GAL) {
     380      delete satData;
     381      iGal.remove();
    364382    }
    365383  }
    … …  
    582600    // ------------------------------------------------
    583601    int iPar = 0;
    584     QMutableVectorIterator<bncParam*> im(_params);
    585     while (im.hasNext()) {
    586       bncParam* par = im.next();
     602    QMutableVectorIterator<bncParam*> it(_params);
     603    while (it.hasNext()) {
     604      bncParam* par = it.next();
    587605      bool removed = false;
    588606      if (par->type == bncParam::AMB_L3) {
    589         if (epoData->satData.find(par->prn) == epoData->satData.end()) {
     607        if (epoData->satDataGPS.find(par->prn) == epoData->satDataGPS.end() &&
     608            epoData->satDataGlo.find(par->prn) == epoData->satDataGlo.end() &&
     609            epoData->satDataGal.find(par->prn) == epoData->satDataGal.end() ) {
    590610          removed = true;
    591611          delete par;
    592           im.remove();
     612          it.remove();
    593613        }
    594614      }
    … …  
    601621    // Add new ambiguity parameters
    602622    // ----------------------------
    603     QMapIterator<QString, t_satData*> it(epoData->satData);
    604     while (it.hasNext()) {
    605       it.next();
    606       t_satData* satData = it.value();
     623    QMapIterator<QString, t_satData*> iGPS(epoData->satDataGPS);
     624    while (iGPS.hasNext()) {
     625      iGPS.next();
     626      t_satData* satData = iGPS.value();
     627      addAmb(satData);
     628    }
     629
     630    QMapIterator<QString, t_satData*> iGlo(epoData->satDataGlo);
     631    while (iGlo.hasNext()) {
     632      iGlo.next();
     633      t_satData* satData = iGlo.value();
     634      addAmb(satData);
     635    }
     636
     637    QMapIterator<QString, t_satData*> iGal(epoData->satDataGal);
     638    while (iGal.hasNext()) {
     639      iGal.next();
     640      t_satData* satData = iGal.value();
    607641      addAmb(satData);
    608642    }
    … …  
    657691  // ----------------------
    658692  if (update_p(epoData) != success) {
    659     emit newMessage(_log, false);
    660693    return failure;
    661694  }
    … …  
    903936// Outlier Detection
    904937////////////////////////////////////////////////////////////////////////////
    905 bool bncModel::outlierDetection(int iPhase, const ColumnVector& vv,
    906                                 QMap<QString, t_satData*>& satData) {
     938int bncModel::outlierDetection(int iPhase, const SymmetricMatrix& QQsav,
     939                               const ColumnVector& vv,
     940                               QMap<QString, t_satData*>& satDataGPS,
     941                               QMap<QString, t_satData*>& satDataGlo,
     942                               QMap<QString, t_satData*>& satDataGal) {
    907943
    908944  Tracer tracer("bncModel::outlierDetection");
    909945
    910   QString prn;
    911   double  maxRes  = 0.0;
    912   findMaxRes(vv, satData, prn, maxRes);
    913 
    914   if      (iPhase == 1 && maxRes > MAXRES_PHASE) {
    915     _log += "Outlier Phase " + prn + " "
     946  QString prnCode;
     947  QString prnPhase;
     948  double  maxResCode  = 0.0;
     949  double  maxResPhase = 0.0;
     950
     951  QString prnRemoved;
     952  double  maxRes;
     953
     954  int irc = 0;
     955
     956  if (iPhase == 0) {
     957
     958    // Check GPS Code
     959    // --------------
     960    if (irc == 0) {
     961      findMaxRes(iPhase, vv,satDataGPS, prnCode, maxResCode, prnPhase, maxResPhase);
     962      if (maxResCode > MAXRES_CODE_GPS) {
     963        satDataGPS.remove(prnCode);
     964        prnRemoved = prnCode;
     965        maxRes     = maxResCode;
     966        irc        = 1;
     967      }
     968    }
     969   
     970    // Check Galileo Code
     971    // ------------------
     972    if (irc == 0) {
     973      findMaxRes(iPhase, vv,satDataGal, prnCode, maxResCode, prnPhase, maxResPhase);
     974      if (maxResCode > MAXRES_CODE_GAL) {
     975        satDataGal.remove(prnCode);
     976        prnRemoved = prnCode;
     977        maxRes     = maxResCode;
     978        irc        = 1;
     979      }
     980    }
     981  }
     982
     983  else {
     984
     985    // Check Glonass Phase
     986    // -------------------
     987    if (irc == 0) {
     988      findMaxRes(iPhase, vv,satDataGlo, prnCode, maxResCode, prnPhase, maxResPhase);
     989      if (maxResPhase > MAXRES_PHASE_GLO) {
     990        satDataGlo.remove(prnPhase);
     991        prnRemoved = prnPhase;
     992        maxRes     = maxResPhase;
     993        irc        = 1;
     994      }
     995    }
     996   
     997    // Check Galileo Phase
     998    // -------------------
     999    if (irc == 0) {
     1000      findMaxRes(iPhase, vv,satDataGal, prnCode, maxResCode, prnPhase, maxResPhase);
     1001      if      (maxResPhase > MAXRES_PHASE_GAL) {
     1002        satDataGal.remove(prnPhase);
     1003        prnRemoved = prnPhase;
     1004        maxRes     = maxResPhase;
     1005        irc        = 1;
     1006      }
     1007    }
     1008   
     1009    // Check GPS Phase
     1010    // ---------------
     1011    if (irc == 0) {
     1012      findMaxRes(iPhase, vv,satDataGPS, prnCode, maxResCode, prnPhase, maxResPhase);
     1013      if      (maxResPhase > MAXRES_PHASE_GPS) {
     1014        satDataGPS.remove(prnPhase);
     1015        prnRemoved = prnPhase;
     1016        maxRes     = maxResPhase;
     1017        irc        = 1;
     1018      }
     1019    }
     1020  }
     1021 
     1022  if (irc != 0) {
     1023    _log += "Outlier " + prnRemoved.toAscii() + " "
    9161024          + QByteArray::number(maxRes, 'f', 3) + "\n";
    917     return true;
    918   }
    919   else if (iPhase == 0 && maxRes > MAXRES_CODE) {
    920     _log += "Outlier Code  " + prn + " "
    921           + QByteArray::number(maxRes, 'f', 3) + "\n";
    922     return true;
    923   }
    924   else {
    925     return false;
    926   }
     1025    _QQ = QQsav;
     1026    QVectorIterator<bncParam*> itPar(_params);
     1027    while (itPar.hasNext()) {
     1028      bncParam* par = itPar.next();
     1029      if (par->type == bncParam::AMB_L3 && par->prn == prnRemoved) {
     1030        par->numEpo = 0;
     1031      }
     1032    }
     1033  }
     1034
     1035  return irc;
    9271036}
    9281037
    … …  
    9491058}
    9501059
    951 //
     1060////
    9521061//////////////////////////////////////////////////////////////////////////////
    9531062void bncModel::kalman(const Matrix& AA, const ColumnVector& ll,
    … …  
    11211230  }
    11221231
    1123   // Remember Observation Index
    1124   // --------------------------
    1125   ++iObs;
    1126   satData->obsIndex = iObs;
    1127 
    11281232  // Phase Observations
    11291233  // ------------------
    11301234  if (iPhase == 1) {
     1235    ++iObs;
    11311236    ll(iObs)      = satData->L3 - cmpValue(satData, true);
    1132     double sigL3 = _sigL3;
     1237    PP(iObs,iObs) = 1.0 / (_sigL3 * _sigL3) / (ellWgtCoef * ellWgtCoef);
    11331238    if (satData->system() == 'R') {
    1134       sigL3 *= GLONASS_WEIGHT_FACTOR;
    1135     }
    1136     PP(iObs,iObs) = 1.0 / (sigL3 * sigL3) / (ellWgtCoef * ellWgtCoef);
     1239      PP(iObs,iObs) /= 25.0;
     1240    }
    11371241    for (int iPar = 1; iPar <= _params.size(); iPar++) {
    11381242      if (_params[iPar-1]->type == bncParam::AMB_L3 &&
    … …  
    11421246      AA(iObs, iPar) = _params[iPar-1]->partial(satData, true);
    11431247    }
     1248    satData->indexPhase = iObs;
    11441249  }
    11451250
    … …  
    11471252  // -----------------
    11481253  else {
     1254    ++iObs;
    11491255    ll(iObs)      = satData->P3 - cmpValue(satData, false);
    11501256    PP(iObs,iObs) = 1.0 / (_sigP3 * _sigP3) / (ellWgtCoef * ellWgtCoef);
    … …  
    11521258      AA(iObs, iPar) = _params[iPar-1]->partial(satData, false);
    11531259    }
     1260    satData->indexCode = iObs;
    11541261  }
    11551262}
    … …  
    11571264//
    11581265///////////////////////////////////////////////////////////////////////////
    1159 QByteArray bncModel::printRes(int iPhase, const ColumnVector& vv,
    1160                               const QMap<QString, t_satData*>& satDataMap) {
    1161 
     1266void bncModel::printRes(int iPhase, const ColumnVector& vv,
     1267                        ostringstream& str, t_satData* satData) {
    11621268  Tracer tracer("bncModel::printRes");
    1163 
    1164   ostringstream str;
    1165   str.setf(ios::fixed);
    1166        
    1167   QMapIterator<QString, t_satData*> it(satDataMap);
    1168   while (it.hasNext()) {
    1169     it.next();
    1170     t_satData* satData = it.value();
    1171     if (satData->obsIndex != 0) {
    1172       str << _time.timestr(1)
    1173           << " RES " << satData->prn.toAscii().data()
    1174           << (iPhase ? "   L3 " : "   P3 ")
    1175           << setw(9) << setprecision(4) << vv(satData->obsIndex) << endl;
    1176     }
    1177   }
    1178 
    1179   return QByteArray(str.str().c_str());
     1269  if (iPhase == 1) {
     1270    str << _time.timestr(1)
     1271        << " RES " << satData->prn.toAscii().data() << "   L3 "
     1272        << setw(9) << setprecision(4) << vv(satData->indexPhase) << endl;
     1273  }
     1274  else {
     1275    str << _time.timestr(1)
     1276        << " RES " << satData->prn.toAscii().data() << "   P3 "
     1277        << setw(9) << setprecision(4) << vv(satData->indexCode) << endl;
     1278  }
    11801279}
    11811280
    11821281//
    11831282///////////////////////////////////////////////////////////////////////////
    1184 void bncModel::findMaxRes(const ColumnVector& vv,
     1283void bncModel::findMaxRes(int iPhase, const ColumnVector& vv,
    11851284                          const QMap<QString, t_satData*>& satData,
    1186                           QString& prn,  double& maxRes) {
    1187 
     1285                          QString& prnCode,  double& maxResCode,
     1286                          QString& prnPhase, double& maxResPhase) {
    11881287  Tracer tracer("bncModel::findMaxRes");
    1189 
    1190   maxRes  = 0.0;
     1288  maxResCode  = 0.0;
     1289  maxResPhase = 0.0;
    11911290
    11921291  QMapIterator<QString, t_satData*> it(satData);
    … …  
    11941293    it.next();
    11951294    t_satData* satData = it.value();
    1196     if (satData->obsIndex != 0 && fabs(vv(satData->obsIndex)) > maxRes) {
    1197       maxRes = fabs(vv(satData->obsIndex));
    1198       prn    = satData->prn;
     1295    if (iPhase == 0) {
     1296      if (satData->indexCode) {
     1297        if (fabs(vv(satData->indexCode)) > maxResCode) {
     1298          maxResCode = fabs(vv(satData->indexCode));
     1299          prnCode    = satData->prn;
     1300        }
     1301      }
     1302    }
     1303    else {
     1304      if (satData->indexPhase) {
     1305        if (fabs(vv(satData->indexPhase)) > maxResPhase) {
     1306          maxResPhase = fabs(vv(satData->indexPhase));
     1307          prnPhase    = satData->prn;
     1308        }
     1309      }
    11991310    }
    12001311  }
    … …  
    12071318  Tracer tracer("bncModel::update_p");
    12081319
    1209   // Save Variance-Covariance Matrix, and Status Vector
    1210   // --------------------------------------------------
    1211   rememberState(epoData);
    1212 
    1213   ColumnVector dx;
    1214 
    1215   std::vector<QString> allPrns;
    1216   QMapIterator<QString, t_satData*> it(epoData->satData);
    1217   while (it.hasNext()) {
    1218     it.next();
    1219     t_satData* satData = it.value();
    1220     allPrns.push_back(satData->prn);
    1221   }
    1222   std::vector<QString> usedPrns;
    1223 
    1224   // Try with all satellites, then with all minus one, etc.
    1225   // ------------------------------------------------------
    1226   const unsigned MAX_NEGLECT = 1;
    1227   for (unsigned nNeglected = 0; nNeglected <= MAX_NEGLECT; nNeglected++) {
    1228     usedPrns = allPrns;
    1229 
    1230     for (unsigned ii = 0; ii < nNeglected && usedPrns.size() > 0; ii++) {
    1231       usedPrns.pop_back();
    1232     }
    1233 
    1234     // Loop over all Combinations of used Satellites
    1235     // ---------------------------------------------
     1320  rememberState();
     1321
     1322  for (int iPhase = 0; iPhase <= (_usePhase ? 1 : 0); iPhase++) {
     1323
     1324    SymmetricMatrix QQsav;
     1325    ColumnVector    vv;
     1326    ColumnVector    dx;
     1327
    12361328    do {
    1237 
    1238       QByteArray strResCode;
    1239       QByteArray strResPhase;
    1240       QString    strNeglected;
    1241 
    1242       // Remove Neglected Satellites from epoData
    1243       // ----------------------------------------
    1244       for (unsigned ip = 0; ip < allPrns.size(); ip++) {
    1245         QString prn = allPrns[ip];
    1246         if ( !findInVector(usedPrns, prn) ) {
    1247           epoData->satData.remove(prn);
    1248           strNeglected += prn + " ";
    1249         }
    1250       }
    1251       if (epoData->sizeSys('G') < MINOBS) {
    1252         continue;
    1253       }
    12541329
    12551330      // Bancroft Solution
    12561331      // -----------------
    1257       if (cmpBancroft(epoData) != success) {
    1258         continue;
    1259       }
    1260 
    1261       // First update using code observations, then phase observations
    1262       // -------------------------------------------------------------     
    1263       for (int iPhase = 0; iPhase <= (_usePhase ? 1 : 0); iPhase++) {
     1332      if (iPhase == 0) {     
     1333        if (cmpBancroft(epoData) != success) {
     1334          restoreState();
     1335          emit newMessage(_log, false);
     1336          return failure;
     1337        }
     1338      }
     1339      else {
     1340        if (epoData->sizeGPS() < MINOBS) {
     1341          restoreState();
     1342          _log += "bncModel::update_p: not enough data\n";
     1343          emit newMessage(_log, false);
     1344          return failure;
     1345        }
     1346        unsigned numSatNoSlip = 0;
     1347        QVectorIterator<bncParam*> itPar(_params);
     1348        while (itPar.hasNext()) {
     1349          bncParam* par = itPar.next();
     1350          if (par->type == bncParam::AMB_L3 && par->prn[0] == 'G') {
     1351            if (par->numEpo >= 1) {
     1352              ++numSatNoSlip;
     1353            }
     1354          }
     1355        }
     1356        if (numSatNoSlip > 0 && numSatNoSlip < MINOBS) {
     1357          restoreState();
     1358          _log += "bncModel::update_p: not enough GPS satellites without cycle-slips\n";
     1359          emit newMessage(_log, false);
     1360          return failure;
     1361        }
     1362      }
    12641363     
    1265         // Status Prediction
    1266         // -----------------
    1267         predict(iPhase, epoData);
    1268        
    1269         // Create First-Design Matrix
    1270         // --------------------------
    1271         unsigned nPar = _params.size();
    1272         unsigned nObs = 0;
    1273         if (iPhase == 0) {
    1274           nObs = epoData->sizeAll() - epoData->sizeSys('R'); // Glonass code not used
    1275         }
    1276         else {
    1277           nObs = epoData->sizeAll();
    1278         }
    1279        
    1280         // Prepare first-design Matrix, vector observed-computed
    1281         // -----------------------------------------------------
    1282         Matrix          AA(nObs, nPar);  // first design matrix
    1283         ColumnVector    ll(nObs);        // tems observed-computed
    1284         DiagonalMatrix  PP(nObs); PP = 0.0;
    1285        
    1286         unsigned iObs = 0;
    1287         QMapIterator<QString, t_satData*> it(epoData->satData);
    1288         while (it.hasNext()) {
    1289           it.next();
    1290           t_satData* satData = it.value();
    1291           if (iPhase == 1 || satData->system() != 'R') {
    1292             QString prn = satData->prn;
    1293             addObs(iPhase, iObs, satData, AA, ll, PP);
     1364      // Status Prediction
     1365      // -----------------
     1366      predict(iPhase, epoData);
     1367     
     1368      // Create First-Design Matrix
     1369      // --------------------------
     1370      unsigned nPar = _params.size();
     1371      unsigned nObs = 0;
     1372      if (iPhase == 0) {
     1373        nObs = epoData->sizeGPS() + epoData->sizeGal(); // Glonass code not used
     1374      }
     1375      else {
     1376        nObs = epoData->sizeGPS() + epoData->sizeGal() + epoData->sizeGlo();
     1377      }
     1378   
     1379      Matrix          AA(nObs, nPar);  // first design matrix
     1380      ColumnVector    ll(nObs);        // tems observed-computed
     1381      DiagonalMatrix  PP(nObs); PP = 0.0;
     1382     
     1383      unsigned iObs = 0;
     1384   
     1385      // GPS
     1386      // ---
     1387      QMapIterator<QString, t_satData*> itGPS(epoData->satDataGPS);
     1388      while (itGPS.hasNext()) {
     1389        itGPS.next();
     1390        t_satData* satData = itGPS.value();
     1391        addObs(iPhase, iObs, satData, AA, ll, PP);
     1392      }
     1393   
     1394      // Glonass
     1395      // -------
     1396      if (iPhase == 1) {
     1397        QMapIterator<QString, t_satData*> itGlo(epoData->satDataGlo);
     1398        while (itGlo.hasNext()) {
     1399          itGlo.next();
     1400          t_satData* satData = itGlo.value();
     1401          addObs(iPhase, iObs, satData, AA, ll, PP);
     1402        }
     1403      }
     1404
     1405      // Galileo
     1406      // -------
     1407      QMapIterator<QString, t_satData*> itGal(epoData->satDataGal);
     1408      while (itGal.hasNext()) {
     1409        itGal.next();
     1410        t_satData* satData = itGal.value();
     1411        addObs(iPhase, iObs, satData, AA, ll, PP);
     1412      }
     1413   
     1414      // Compute Filter Update
     1415      // ---------------------
     1416      QQsav = _QQ;
     1417   
     1418      kalman(AA, ll, PP, _QQ, dx);
     1419   
     1420      vv = ll - AA * dx;
     1421
     1422      // Print Residuals
     1423      // ---------------
     1424      if (true) {
     1425        ostringstream str;
     1426        str.setf(ios::fixed);
     1427   
     1428        QMapIterator<QString, t_satData*> itGPS(epoData->satDataGPS);
     1429        while (itGPS.hasNext()) {
     1430          itGPS.next();
     1431          t_satData* satData = itGPS.value();
     1432          printRes(iPhase, vv, str, satData);
     1433        }
     1434        if (iPhase == 1) {
     1435          QMapIterator<QString, t_satData*> itGlo(epoData->satDataGlo);
     1436          while (itGlo.hasNext()) {
     1437            itGlo.next();
     1438            t_satData* satData = itGlo.value();
     1439            printRes(iPhase, vv, str, satData);
    12941440          }
    12951441        }
    1296 
    1297         // Compute Filter Update
    1298         // ---------------------
    1299         kalman(AA, ll, PP, _QQ, dx);
    1300        
    1301         ColumnVector vv = ll - AA * dx;
    1302        
    1303         // Print Residuals
    1304         // ---------------
    1305         if (iPhase == 0) {
    1306           strResCode  = printRes(iPhase, vv, epoData->satData);
    1307         }
    1308         else {
    1309           strResPhase = printRes(iPhase, vv, epoData->satData);
    1310         }
    1311 
    1312         // Check the residuals
    1313         // -------------------
    1314         if ( outlierDetection(iPhase, vv, epoData->satData) ) {
    1315           restoreState(epoData);
    1316           break;
    1317         }
    1318 
    1319         // Set estimated values
    1320         // --------------------
    1321         else if (!_usePhase || iPhase == 1) {
    1322           QVectorIterator<bncParam*> itPar(_params);
    1323           while (itPar.hasNext()) {
    1324             bncParam* par = itPar.next();
    1325             par->xx += dx(par->index);
    1326           }
    1327 
    1328           // Print Neglected PRNs
    1329           // --------------------
    1330           if (nNeglected > 0) {
    1331             _log += "Neglected PRNs: " + strNeglected + '\n';
    1332           }
    1333 
    1334           _log += strResCode + strResPhase;
    1335        
    1336           return success;
    1337         }
    1338 
    1339       } // for (int iPhase = 0; iPhase <= (_usePhase ? 1 : 0); iPhase++)
    1340 
    1341     } while ( next_combination(allPrns.begin(), allPrns.end(),
    1342                                usedPrns.begin(), usedPrns.end()) );
    1343 
    1344   } // for (unsigned nNeglected
    1345 
    1346   return failure;
     1442        QMapIterator<QString, t_satData*> itGal(epoData->satDataGal);
     1443        while (itGal.hasNext()) {
     1444          itGal.next();
     1445          t_satData* satData = itGal.value();
     1446          printRes(iPhase, vv, str, satData);
     1447        }
     1448        _log += str.str().c_str();
     1449      }
     1450   
     1451    } while (outlierDetection(iPhase, QQsav, vv, epoData->satDataGPS,
     1452                              epoData->satDataGlo, epoData->satDataGal) != 0);
     1453
     1454    // Update Parameters
     1455    // -----------------
     1456    QVectorIterator<bncParam*> itPar(_params);
     1457    while (itPar.hasNext()) {
     1458      bncParam* par = itPar.next();
     1459      par->xx += dx(par->index);
     1460    }
     1461  }
     1462
     1463  return success;
    13471464}
    13481465
    13491466// Remeber Original State Vector and Variance-Covariance Matrix
    13501467////////////////////////////////////////////////////////////////////////////
    1351 void bncModel::rememberState(t_epoData* epoData) {
     1468void bncModel::rememberState() {
    13521469
    13531470  _QQ_sav = _QQ;
    … …  
    13651482    _params_sav.push_back(new bncParam(*par));
    13661483  }
    1367 
    1368   _epoData_sav->deepCopy(epoData);
    13691484}
    13701485
    13711486// Restore Original State Vector and Variance-Covariance Matrix
    13721487////////////////////////////////////////////////////////////////////////////
    1373 void bncModel::restoreState(t_epoData* epoData) {
     1488void bncModel::restoreState() {
    13741489
    13751490  _QQ = _QQ_sav;
    … …  
    13871502    _params.push_back(new bncParam(*par));
    13881503  }
    1389 
    1390   epoData->deepCopy(_epoData_sav);
    1391 }
     1504}
     1505

  • trunk/BNC/bncmodel.h

    r3400 r3405  
    2929#include <QtNetwork>
    3030#include <newmat.h>
    31 #include <vector>
    3231
    3332#include "bncconst.h"
    … …  
    6059  bncModel(QByteArray staID);
    6160  ~bncModel();
     61  t_irc cmpBancroft(t_epoData* epoData);
    6262  t_irc update(t_epoData* epoData);
    6363  bncTime time()  const {return _time;}
    … …  
    9494
    9595 private:
    96   t_irc cmpBancroft(t_epoData* epoData);
    9796  void   reset();
    9897  void   cmpEle(t_satData* satData);
    … …  
    10099  void   addObs(int iPhase, unsigned& iObs, t_satData* satData,
    101100                Matrix& AA, ColumnVector& ll, DiagonalMatrix& PP);
    102   QByteArray printRes(int iPhase, const ColumnVector& vv,
    103                       const QMap<QString, t_satData*>& satDataMap);
    104   void   findMaxRes(const ColumnVector& vv,
     101  void  printRes(int iPhase, const ColumnVector& vv,
     102                  std::ostringstream& str, t_satData* satData);
     103  void   findMaxRes(int iPhase, const ColumnVector& vv,
    105104                    const QMap<QString, t_satData*>& satData,
    106                     QString& prn,  double& maxRes);
     105                    QString& prnCode,  double& maxResCode,
     106                    QString& prnPhase, double& maxResPhase);
    107107  double cmpValue(t_satData* satData, bool phase);
    108108  double delay_saast(double Ele);
    109109  void   predict(int iPhase, t_epoData* epoData);
    110110  t_irc  update_p(t_epoData* epoData);
    111   bool   outlierDetection(int iPhase, const ColumnVector& vv,
    112                           QMap<QString, t_satData*>& satData);
     111  int    outlierDetection(int iPhase, const SymmetricMatrix& QQsav,
     112                          const ColumnVector& vv,
     113                          QMap<QString, t_satData*>& satDataGPS,
     114                          QMap<QString, t_satData*>& satDataGlo,
     115                          QMap<QString, t_satData*>& satDataGal);
    113116  void writeNMEAstr(const QString& nmStr);
    114117
    … …  
    118121  bncTime  _startTime;
    119122
    120   void rememberState(t_epoData* epoData);
    121   void restoreState(t_epoData* epoData);
     123  void rememberState();
     124  void restoreState();
    122125
    123126  class pppPos {
    … …  
    139142  QVector<bncParam*>    _params_sav;
    140143  SymmetricMatrix       _QQ_sav;
    141   t_epoData*            _epoData_sav;
    142144  ColumnVector          _xcBanc;
    143145  ColumnVector          _ellBanc;

  • trunk/BNC/bncpppclient.cpp

    r3383 r3405  
    194194      satData->lambda3 = c1 * t_CST::c / f1 + c2 * t_CST::c / f2;
    195195
    196       _epoData.back()->satData[satData->prn] = satData;
     196      _epoData.back()->satDataGPS[satData->prn] = satData;
    197197    }
    198198    else {
    … …  
    227227      satData->lambda3 = c1 * t_CST::c / f1 + c2 * t_CST::c / f2;
    228228
    229       _epoData.back()->satData[satData->prn] = satData;
     229      _epoData.back()->satDataGlo[satData->prn] = satData;
    230230    }
    231231    else {
    … …  
    250250      satData->L3      = c1 * satData->L1 + c5 * satData->L5;
    251251      satData->lambda3 = c1 * t_CST::c / f1 + c5 * t_CST::c / f5;
    252       _epoData.back()->satData[satData->prn] = satData;
     252      _epoData.back()->satDataGal[satData->prn] = satData;
    253253    }
    254254    else {
    … …  
    459459  // Data Pre-Processing
    460460  // -------------------
    461   QMutableMapIterator<QString, t_satData*> it(_epoData.front()->satData);
    462   while (it.hasNext()) {
    463     it.next();
    464     QString    prn     = it.key();
    465     t_satData* satData = it.value();
     461  QMutableMapIterator<QString, t_satData*> iGPS(_epoData.front()->satDataGPS);
     462  while (iGPS.hasNext()) {
     463    iGPS.next();
     464    QString    prn     = iGPS.key();
     465    t_satData* satData = iGPS.value();
    466466
    467467    if (cmpToT(satData) != success) {
    468468      delete satData;
    469       it.remove();
     469      iGPS.remove();
     470      continue;
     471    }
     472  }
     473
     474  QMutableMapIterator<QString, t_satData*> iGlo(_epoData.front()->satDataGlo);
     475  while (iGlo.hasNext()) {
     476    iGlo.next();
     477    QString    prn     = iGlo.key();
     478    t_satData* satData = iGlo.value();
     479
     480    if (cmpToT(satData) != success) {
     481      delete satData;
     482      iGlo.remove();
     483      continue;
     484    }
     485  }
     486
     487  QMutableMapIterator<QString, t_satData*> iGal(_epoData.front()->satDataGal);
     488  while (iGal.hasNext()) {
     489    iGal.next();
     490    QString    prn     = iGal.key();
     491    t_satData* satData = iGal.value();
     492
     493    if (cmpToT(satData) != success) {
     494      delete satData;
     495      iGal.remove();
    470496      continue;
    471497    }

  • trunk/BNC/bncpppclient.h

    r3392 r3405  
    3535 public:
    3636  t_satData() {
    37     obsIndex = 0;
     37    indexCode  = 0;
     38    indexPhase = 0;
    3839  }
    3940  ~t_satData() {}
    … …  
    5657  bool         slipFlag;
    5758  double       lambda3;
    58   unsigned     obsIndex;
     59  unsigned     indexCode;
     60  unsigned     indexPhase;
    5961  char system() const {return prn.toAscii()[0];}
    6062};
    … …  
    6365 public:
    6466  t_epoData() {}
    65 
    6667  ~t_epoData() {
    67     clear();
    68   }
    69 
    70   void clear() {
    71     QMapIterator<QString, t_satData*> it(satData);
    72     while (it.hasNext()) {
    73       it.next();
    74       delete it.value();
     68    QMapIterator<QString, t_satData*> itGPS(satDataGPS);
     69    while (itGPS.hasNext()) {
     70      itGPS.next();
     71      delete itGPS.value();
    7572    }
    76     satData.clear();
    77   }
    78 
    79   void deepCopy(const t_epoData* from) {
    80     clear();
    81     tt = from->tt;
    82     QMapIterator<QString, t_satData*> it(from->satData);
    83     while (it.hasNext()) {
    84       it.next();
    85       satData[it.key()] = new t_satData(*it.value());
     73    QMapIterator<QString, t_satData*> itGlo(satDataGlo);
     74    while (itGlo.hasNext()) {
     75      itGlo.next();
     76      delete itGlo.value();
     77    }
     78    QMapIterator<QString, t_satData*> itGal(satDataGal);
     79    while (itGal.hasNext()) {
     80      itGal.next();
     81      delete itGal.value();
    8682    }
    8783  }
    88 
    89   unsigned sizeSys(char system) const {
    90     unsigned ans = 0;
    91     QMapIterator<QString, t_satData*> it(satData);
    92     while (it.hasNext()) {
    93       it.next();
    94       if (it.value()->system() == system) {
    95         ++ans;
    96       }
    97     }
    98     return ans;
    99   }
    100   unsigned sizeAll() const {return satData.size();}
    101 
    102   bncTime                   tt;
    103   QMap<QString, t_satData*> satData;
     84  unsigned sizeGPS() const {return satDataGPS.size();}
     85  unsigned sizeGlo() const {return satDataGlo.size();}
     86  unsigned sizeGal() const {return satDataGal.size();}
     87  unsigned sizeAll() const {return satDataGPS.size() + satDataGlo.size() +
     88                                   satDataGal.size();}
     89  bncTime                    tt;
     90  QMap<QString, t_satData*> satDataGPS;
     91  QMap<QString, t_satData*> satDataGlo;
     92  QMap<QString, t_satData*> satDataGal;
    10493};
    10594

  • trunk/BNC/bncutils.cpp

    r3394 r3405  
    374374               dateTime.time().msec() / 1000.0) / 60.0) / 60.0) / 24.0;
    375375}
    376 
    377 //
    378 ////////////////////////////////////////////////////////////////////////////
    379 bool findInVector(const vector<QString>& vv, const QString& str) {
    380   std::vector<QString>::const_iterator it;
    381   for (it = vv.begin(); it != vv.end(); ++it) {
    382     if ( (*it) == str) {
    383       return true;
    384     }
    385   }
    386   return false;
    387 }
    388 

  • trunk/BNC/bncutils.h

    r3394 r3405  
    2525#ifndef BNCUTILS_H
    2626#define BNCUTILS_H
    27 
    28 #include <vector>
    2927
    3028#include <QString>
    … …  
    7068void mjdFromDateAndTime(const QDateTime& dateTime, int& mjd, double& dayfrac);
    7169
    72 bool findInVector(const std::vector<QString>& vv, const QString& str);
    73 
    7470#endif
Note: See TracChangeset for help on using the changeset viewer.