Added some gpstk based application, to better understand gpstk
authorBenoit Papillault <benoit.papillault@free.fr>
Wed, 21 Dec 2011 21:57:48 +0000 (22:57 +0100)
committerBenoit Papillault <benoit.papillault@free.fr>
Wed, 21 Dec 2011 21:57:48 +0000 (22:57 +0100)
gpstk-prsolution.cpp [new file with mode: 0644]
gpstk-solution2.cpp [new file with mode: 0644]

diff --git a/gpstk-prsolution.cpp b/gpstk-prsolution.cpp
new file mode 100644 (file)
index 0000000..8d3dfda
--- /dev/null
@@ -0,0 +1,90 @@
+#include <iostream>
+#include "PRSolution.hpp"
+#include "SP3EphemerisStore.hpp"
+#include "Position.hpp"
+
+using namespace std;
+using namespace gpstk;
+
+int main() {
+
+       DayTime dayTime(1667, 66448.850000);
+       PRSolution raimSolver;
+
+       // raimSolver.RMSLimit = 10;
+       //raimSolver.Debug = true;
+
+       vector<SatID> satelliteVector;
+       vector<double> pseudoRangeVector;
+
+       // push PRN 
+       satelliteVector.push_back(SatID( 2, SatID::systemGPS));
+       satelliteVector.push_back(SatID(10, SatID::systemGPS));
+       satelliteVector.push_back(SatID( 4, SatID::systemGPS));
+       satelliteVector.push_back(SatID( 5, SatID::systemGPS));
+       satelliteVector.push_back(SatID(13, SatID::systemGPS));
+       satelliteVector.push_back(SatID( 7, SatID::systemGPS));
+       satelliteVector.push_back(SatID(29, SatID::systemGPS));
+
+       // push pseudo range
+       pseudoRangeVector.push_back(21375681.860758);
+       pseudoRangeVector.push_back(20192220.627359);
+       pseudoRangeVector.push_back(21808239.424646);
+       pseudoRangeVector.push_back(23834064.766712);
+       pseudoRangeVector.push_back(20794005.449315);
+       pseudoRangeVector.push_back(21213754.205750);
+       pseudoRangeVector.push_back(25306963.567147);
+
+       SP3EphemerisStore ephemerisStore;
+       // from
+       // ftp://rgpdata.ign.fr/pub/products/ephemerides/1667/esr16670.sp3.Z
+       ephemerisStore.loadFile("esr16670.sp3");
+
+       // no tropo model
+       ZeroTropModel noTropModel;
+
+       raimSolver.RAIMCompute(dayTime, 
+                              satelliteVector,
+                              pseudoRangeVector,
+                              ephemerisStore,
+                              &noTropModel);
+
+       if (raimSolver.isValid()) {
+               cout << fixed << setprecision(3);
+               cout << "x = " << raimSolver.Solution[0] << endl;
+               cout << "y = " << raimSolver.Solution[1] << endl;
+               cout << "z = " << raimSolver.Solution[2] << endl;
+
+               /* we got (the results depends on RMSLimit) : 
+                  x = 4205773.319
+                  y =  166187.065
+                  z = 4776312.107
+
+                  SiRF Star III was computing (sat 2, 10 & 4) :
+                    X-position : 4205766 m
+                    Y-position :  166220 m
+                    Z-position : 4776300 m
+               */
+
+               Position pos(raimSolver.Solution[0],
+                            raimSolver.Solution[1],
+                            raimSolver.Solution[2],
+                            Position::Cartesian);
+
+               cout << fixed << setprecision(6);
+               cout << "longitude = " << pos.asGeodetic().getLongitude()
+                    << endl;
+               cout << "latitude  = " << pos.asGeodetic().geodeticLatitude()
+                    << endl;
+               cout << "altitude  = " << pos.asGeodetic().getAltitude()
+                    << endl;
+
+               /* current position is (as reported by gpsd) : 
+                  longitude =   2.262773
+                  latitude  =  48.803081
+                  altitude  = 156.445000
+               */
+       }
+
+       return 0;
+}
diff --git a/gpstk-solution2.cpp b/gpstk-solution2.cpp
new file mode 100644 (file)
index 0000000..eb3f4f8
--- /dev/null
@@ -0,0 +1,190 @@
+#include <iostream>
+#include "SP3EphemerisStore.hpp"
+#include "Position.hpp"
+
+#include "DataStructures.hpp"
+#include "SimpleFilter.hpp"
+#include "IonoModel.hpp"
+#include "IonoModelStore.hpp"
+#include "ModelObs.hpp"
+#include "SolverLMS.hpp"
+#include "SolverPPP.hpp"
+
+using namespace std;
+using namespace gpstk;
+
+int main() {
+
+       DayTime dayTime(1667, 66448.850000);
+
+       /* satellites list */
+       SatID sat2 ( 2, SatID::systemGPS);
+       SatID sat4 ( 4, SatID::systemGPS);
+       SatID sat5 ( 5, SatID::systemGPS);
+       SatID sat7 ( 7, SatID::systemGPS);
+       SatID sat10(10, SatID::systemGPS);
+       SatID sat13(13, SatID::systemGPS);
+       SatID sat29(29, SatID::systemGPS);
+       
+       /*gnssRinex*/ gnssSatTypeValue gRin;
+       /* header = sourceEpochTypeHeader + inheritance from gnssSatTypeValue
+          = gnssData<sourceEpochHeader, satTypeValueMap
+          satTypeValueMap = std::map<SatID, typeValueMap> 
+          typeValueMap == std::map<TypeID, double>
+       */
+
+       gRin.header.epoch = dayTime;
+
+       typeValueMap TVmap;
+       satTypeValueMap STVmap;
+
+       STVmap.clear();
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 21375681.860758;
+       STVmap[sat2 ] = TVmap;
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 21808239.424646;
+       STVmap[sat4 ] = TVmap;
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 23834064.766712;
+       STVmap[sat5 ] = TVmap;
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 21213754.205750;
+       STVmap[sat7 ] = TVmap;
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 20192220.627359;
+       STVmap[sat10] = TVmap;
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 20794005.449315;
+       STVmap[sat13] = TVmap;
+
+       TVmap.clear();
+       TVmap[TypeID::C1] = 25306963.567147;
+       STVmap[sat29] = TVmap;
+
+       gRin.body = STVmap;
+
+       SP3EphemerisStore ephemerisStore;
+       // from
+       // ftp://rgpdata.ign.fr/pub/products/ephemerides/1667/esr16670.sp3.Z
+       ephemerisStore.loadFile("esr16670.sp3");
+
+       /* we got (the results depends on RMSLimit) : 
+                  x = 4205773.319
+                  y =  166187.065
+                  z = 4776312.107
+
+                  SiRF Star III was computing (sat 2, 10 & 4) :
+                    X-position : 4205766 m
+                    Y-position :  166220 m
+                    Z-position : 4776300 m
+       */
+
+               /* current position is (as reported by gpsd) : 
+                  longitude =   2.262773
+                  latitude  =  48.803081
+                  altitude  = 156.445000
+               */
+
+       try {
+               SimpleFilter filter;
+
+               /* FIXME : what is the role of nominalPos ? */
+               Position nominalPos(4205000.0, 166000.0, 4776000.0);
+
+               double ionAlpha[4] = { 0.0 , 0.0, 0.0, 0.0 };
+               double ionBeta[4]  = { 149500.0, -131100.0, 0.0, -131100.0 };
+
+               IonoModel ioModel;
+               ioModel.setModel(ionAlpha, ionBeta);
+
+               IonoModelStore ionoStore;
+               ionoStore.addIonoModel(DayTime::BEGINNING_OF_TIME, ioModel);
+
+               MOPSTropModel mopsTM( nominalPos.getAltitude(),
+                                     nominalPos.getGeodeticLatitude(),
+                                     30 );
+
+
+               ModelObs model(/* IonoModelStore */ ionoStore,
+                              /* TropModel */ mopsTM, 
+                              /* XvtStore<SatID> */ ephemerisStore,
+                              /* TypeID */ TypeID::C1);
+               model.setModelPrepared(false);
+
+               SolverLMS solver;
+
+               cout << fixed << setprecision(6);
+
+               //model.Prepare(gRin.header.epoch, gRin.body);
+
+               cout << "model.getModelPrepared = " 
+                    << model.getModelPrepared() << endl;
+
+               //cout << "GDS data at start : " << endl;
+               //gRin.body.dump(cout, 1);
+               gRin >> filter >> model >> solver;
+               //cout << "GDS data after 'filter+model+solver' : " <<endl;
+               //gRin.body.dump(cout,1);
+
+               cout << "Solution is : " << solver.isValid() << endl;
+
+               cout << "dx = " << solver.getSolution(TypeID::dx) << endl;
+               cout << "dy = " << solver.getSolution(TypeID::dy) << endl;
+               cout << "dz = " << solver.getSolution(TypeID::dz) << endl;
+               cout << "cdt = " << solver.getSolution(TypeID::cdt) << endl;
+
+               cout << "model.rxPos is : " <<endl;
+               cout << "x = " << model.rxPos.X() << endl;
+               cout << "y = " << model.rxPos.Y() << endl;
+               cout << "z = " << model.rxPos.Z() << endl;
+
+               /*
+                 model.rxPos is : 
+                 x = 4205773.614814
+                 y = 166221.288957
+                 z = 4776329.263813
+               */
+
+               cout << "ECEF result : " << endl;
+               cout << "x = " << model.rxPos.X()+solver.getSolution(TypeID::dx)
+                    << endl;
+               cout << "y = " << model.rxPos.Y()+solver.getSolution(TypeID::dy)
+                    << endl;
+               cout << "z = " << model.rxPos.Z()+solver.getSolution(TypeID::dz)
+                    << endl;
+
+               /*
+                 ECEF result : 
+                 x = 4205741.185280
+                 y = 166208.570230
+                 z = 4776301.430305
+                */
+
+               Position pos(model.rxPos.X()+solver.getSolution(TypeID::dx),
+                            model.rxPos.Y()+solver.getSolution(TypeID::dy),
+                            model.rxPos.Z()+solver.getSolution(TypeID::dz),
+                            Position::Cartesian);
+               cout <<"longitude = " << pos.asGeodetic().getLongitude()<<endl;
+               cout <<"latitude = "<<pos.asGeodetic().geodeticLatitude()<<endl;
+               cout <<"altitude = " << pos.asGeodetic().getAltitude() << endl;
+
+               /*
+                 longitude = 2.263120
+                 latitude = 48.803313
+                 altitude = 160.905461
+               */
+       }
+
+       catch (...) {
+               cout << "An exception occurs" << endl;
+       }
+
+       return 0;
+}