Changeset 3412 in ntrip for trunk/BNC


Ignore:
Timestamp:
Sep 3, 2011, 12:31:44 PM (13 years ago)
Author:
mervart
Message:
 
Location:
trunk/BNC
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • trunk/BNC/bncmodel.cpp

    r3411 r3412  
    5252#include "bnctides.h"
    5353#include "bncantex.h"
    54 #include "bnccomb.h"
    5554
    5655using namespace std;
     
    904903// Outlier Detection
    905904////////////////////////////////////////////////////////////////////////////
    906 bool bncModel::outlierDetection(int iPhase, const ColumnVector& vv,
    907                                 QMap<QString, t_satData*>& satData) {
     905QString bncModel::outlierDetection(int iPhase, const ColumnVector& vv,
     906                                   QMap<QString, t_satData*>& satData) {
    908907
    909908  Tracer tracer("bncModel::outlierDetection");
     
    917916    _log += "Outlier Phase " + prn + " "
    918917          + QByteArray::number(maxRes, 'f', 3) + "\n";
    919     return true;
     918    return prn;
    920919  }
    921920  else if (iPhase == 0 && maxRes > MAXRES_CODE) {
    922921    _log += "Outlier Code  " + prn + " "
    923922          + QByteArray::number(maxRes, 'f', 3) + "\n";
    924     return true;
     923    return prn;
    925924  }
    926925  else {
    927     return false;
     926    return QString();
    928927  }
    929928}
     
    12151214  ColumnVector dx;
    12161215
    1217   std::vector<QString> allPrns;
    1218   QMapIterator<QString, t_satData*> it(epoData->satData);
    1219   while (it.hasNext()) {
    1220     it.next();
    1221     t_satData* satData = it.value();
    1222     allPrns.push_back(satData->prn);
    1223   }
    1224   std::vector<QString> usedPrns;
     1216  QString lastOutlierPrn;
    12251217
    12261218  // Try with all satellites, then with all minus one, etc.
    12271219  // ------------------------------------------------------
    1228   for (unsigned nNeglected = 0; nNeglected <= allPrns.size() - MINOBS; nNeglected++) {
    1229     usedPrns = allPrns;
    1230 
    1231     for (unsigned ii = 0; ii < nNeglected && usedPrns.size() > 0; ii++) {
    1232       usedPrns.pop_back();
    1233     }
    1234 
    1235     // Loop over all Combinations of used Satellites
    1236     // ---------------------------------------------
    1237     do {
    1238 
    1239       QByteArray strResCode;
    1240       QByteArray strResPhase;
    1241       QString    strNeglected;
    1242 
    1243       // Remove Neglected Satellites from epoData
    1244       // ----------------------------------------
    1245       unsigned neglectedGPS = 0;
    1246       for (unsigned ip = 0; ip < allPrns.size(); ip++) {
    1247         QString prn = allPrns[ip];
    1248         if ( !findInVector(usedPrns, prn) ) {
    1249           epoData->satData.remove(prn);
    1250           strNeglected += prn + " ";
    1251           if (prn[0] == 'G') {
    1252             ++neglectedGPS;
    1253           }
     1220  while (selectSatellites(lastOutlierPrn, epoData->satData) == success) {
     1221
     1222    QByteArray strResCode;
     1223    QByteArray strResPhase;
     1224
     1225    // Bancroft Solution
     1226    // -----------------
     1227    if (cmpBancroft(epoData) != success) {
     1228      break;
     1229    }
     1230
     1231    // First update using code observations, then phase observations
     1232    // -------------------------------------------------------------     
     1233    for (int iPhase = 0; iPhase <= (_usePhase ? 1 : 0); iPhase++) {
     1234   
     1235      // Status Prediction
     1236      // -----------------
     1237      predict(iPhase, epoData);
     1238     
     1239      // Create First-Design Matrix
     1240      // --------------------------
     1241      unsigned nPar = _params.size();
     1242      unsigned nObs = 0;
     1243      if (iPhase == 0) {
     1244        nObs = epoData->sizeAll() - epoData->sizeSys('R'); // Glonass code not used
     1245      }
     1246      else {
     1247        nObs = epoData->sizeAll();
     1248      }
     1249     
     1250      // Prepare first-design Matrix, vector observed-computed
     1251      // -----------------------------------------------------
     1252      Matrix          AA(nObs, nPar);  // first design matrix
     1253      ColumnVector    ll(nObs);        // tems observed-computed
     1254      DiagonalMatrix  PP(nObs); PP = 0.0;
     1255     
     1256      unsigned iObs = 0;
     1257      QMapIterator<QString, t_satData*> it(epoData->satData);
     1258      while (it.hasNext()) {
     1259        it.next();
     1260        t_satData* satData = it.value();
     1261        if (iPhase == 1 || satData->system() != 'R') {
     1262          QString prn = satData->prn;
     1263          addObs(iPhase, iObs, satData, AA, ll, PP);
    12541264        }
    12551265      }
    1256       if (neglectedGPS > 1) {
    1257         continue;
    1258       }
    1259       if (epoData->sizeSys('G') < MINOBS) {
    1260         continue;
    1261       }
    1262 
    1263       // Bancroft Solution
    1264       // -----------------
    1265       if (cmpBancroft(epoData) != success) {
    1266         continue;
    1267       }
    1268 
    1269       // First update using code observations, then phase observations
    1270       // -------------------------------------------------------------     
    1271       for (int iPhase = 0; iPhase <= (_usePhase ? 1 : 0); iPhase++) {
     1266
     1267      // Compute Filter Update
     1268      // ---------------------
     1269      kalman(AA, ll, PP, _QQ, dx);
    12721270     
    1273         // Status Prediction
    1274         // -----------------
    1275         predict(iPhase, epoData);
    1276        
    1277         // Create First-Design Matrix
    1278         // --------------------------
    1279         unsigned nPar = _params.size();
    1280         unsigned nObs = 0;
    1281         if (iPhase == 0) {
    1282           nObs = epoData->sizeAll() - epoData->sizeSys('R'); // Glonass code not used
    1283         }
    1284         else {
    1285           nObs = epoData->sizeAll();
    1286         }
    1287        
    1288         // Prepare first-design Matrix, vector observed-computed
    1289         // -----------------------------------------------------
    1290         Matrix          AA(nObs, nPar);  // first design matrix
    1291         ColumnVector    ll(nObs);        // tems observed-computed
    1292         DiagonalMatrix  PP(nObs); PP = 0.0;
    1293        
    1294         unsigned iObs = 0;
    1295         QMapIterator<QString, t_satData*> it(epoData->satData);
    1296         while (it.hasNext()) {
    1297           it.next();
    1298           t_satData* satData = it.value();
    1299           if (iPhase == 1 || satData->system() != 'R') {
    1300             QString prn = satData->prn;
    1301             addObs(iPhase, iObs, satData, AA, ll, PP);
    1302           }
    1303         }
    1304 
    1305         // Compute Filter Update
    1306         // ---------------------
    1307         kalman(AA, ll, PP, _QQ, dx);
    1308        
    1309         ColumnVector vv = ll - AA * dx;
    1310        
    1311         // Print Residuals
    1312         // ---------------
    1313         if (iPhase == 0) {
    1314           strResCode  = printRes(iPhase, vv, epoData->satData);
    1315         }
    1316         else {
    1317           strResPhase = printRes(iPhase, vv, epoData->satData);
    1318         }
    1319 
    1320         // Check the residuals
    1321         // -------------------
    1322         if ( outlierDetection(iPhase, vv, epoData->satData) ) {
    1323           restoreState(epoData);
    1324           break;
    1325         }
    1326 
    1327         // Set estimated values
    1328         // --------------------
    1329         else if (!_usePhase || iPhase == 1) {
     1271      ColumnVector vv = ll - AA * dx;
     1272     
     1273      // Print Residuals
     1274      // ---------------
     1275      if (iPhase == 0) {
     1276        strResCode  = printRes(iPhase, vv, epoData->satData);
     1277      }
     1278      else {
     1279        strResPhase = printRes(iPhase, vv, epoData->satData);
     1280      }
     1281
     1282      // Check the residuals
     1283      // -------------------
     1284      lastOutlierPrn = outlierDetection(iPhase, vv, epoData->satData);
     1285
     1286      // No Outlier Detected
     1287      // -------------------
     1288      if (lastOutlierPrn.isEmpty()) {
     1289        if (!_usePhase || iPhase == 1) {
    13301290          QVectorIterator<bncParam*> itPar(_params);
    13311291          while (itPar.hasNext()) {
     
    13331293            par->xx += dx(par->index);
    13341294          }
    1335 
    1336           // Print Neglected PRNs
    1337           // --------------------
    1338           if (nNeglected > 0) {
    1339             _log += "Neglected PRNs: " + strNeglected + '\n';
    1340           }
    1341 
    13421295          _log += strResCode + strResPhase;
    1343        
    13441296          return success;
    13451297        }
    1346 
    1347       } // for (int iPhase = 0; iPhase <= (_usePhase ? 1 : 0); iPhase++)
    1348 
    1349     } while ( next_combination(allPrns.begin(), allPrns.end(),
    1350                                usedPrns.begin(), usedPrns.end()) );
    1351 
    1352   } // for (unsigned nNeglected
    1353 
     1298      }
     1299
     1300      // Outlier Found
     1301      // -------------
     1302      else {
     1303        restoreState(epoData);
     1304        break;
     1305      }
     1306
     1307    } // for iPhase
     1308
     1309  } // while selectSatellites
     1310
     1311  restoreState(epoData);
    13541312  return failure;
    13551313}
     
    13981356  epoData->deepCopy(_epoData_sav);
    13991357}
     1358
     1359//
     1360////////////////////////////////////////////////////////////////////////////
     1361t_irc bncModel::selectSatellites(const QString& lastOutlierPrn,
     1362                                 QMap<QString, t_satData*>& satData) {
     1363
     1364  // First Call
     1365  // ----------
     1366  if (lastOutlierPrn.isEmpty()) {
     1367    _outlierGPS.clear();
     1368    _outlierGlo.clear();
     1369    return success;
     1370  }
     1371
     1372  // Second and next trials
     1373  // ----------------------
     1374  else {
     1375
     1376    if (lastOutlierPrn[0] == 'R') {
     1377      _outlierGlo << lastOutlierPrn;
     1378    }
     1379
     1380    // Remove all Glonass Outliers
     1381    // ---------------------------
     1382    QStringListIterator it(_outlierGlo);
     1383    while (it.hasNext()) {
     1384      QString prn = it.next();
     1385      if (satData.contains(prn)) {
     1386        delete satData.take(prn);
     1387      }
     1388    }
     1389
     1390    if (lastOutlierPrn[0] == 'R') {
     1391      _outlierGPS.clear();
     1392      return success;
     1393    }
     1394
     1395    // GPS Outlier appeared for the first time - try to delete it
     1396    // ----------------------------------------------------------
     1397    if (_outlierGPS.indexOf(lastOutlierPrn) == -1) {
     1398      _outlierGPS << lastOutlierPrn;
     1399      if (satData.contains(lastOutlierPrn)) {
     1400        delete satData.take(lastOutlierPrn);
     1401      }
     1402      return success;
     1403    }
     1404
     1405  }
     1406
     1407  return failure;
     1408}
  • trunk/BNC/bncmodel.h

    r3408 r3412  
    2929#include <QtNetwork>
    3030#include <newmat.h>
    31 #include <vector>
    3231
    3332#include "bncconst.h"
     
    109108  void   predict(int iPhase, t_epoData* epoData);
    110109  t_irc  update_p(t_epoData* epoData);
    111   bool  outlierDetection(int iPhase, const ColumnVector& vv,
    112                           QMap<QString, t_satData*>& satData);
     110  QString outlierDetection(int iPhase, const ColumnVector& vv,
     111                           QMap<QString, t_satData*>& satData);
    113112  void writeNMEAstr(const QString& nmStr);
    114113
     
    120119  void rememberState(t_epoData* epoData);
    121120  void restoreState(t_epoData* epoData);
     121 
     122  t_irc selectSatellites(const QString& lastOutlierPrn,
     123                         QMap<QString, t_satData*>& satData);
    122124
    123125  class pppPos {
     
    168170  bncAntex*             _antex;
    169171  QString               _antennaName;
     172  QStringList           _outlierGPS;
     173  QStringList           _outlierGlo;
    170174};
    171175
Note: See TracChangeset for help on using the changeset viewer.