Index: /trunk/BNC/scripts/Bnc.pm =================================================================== --- /trunk/BNC/scripts/Bnc.pm (revision 10387) +++ /trunk/BNC/scripts/Bnc.pm (revision 10388) @@ -493,6 +493,6 @@ my $epochSec = 0; my $epochDiff = 0; - my ( @hlp, @N, @E, @U, %SATNUM, @TRP, @CLK, @OFF_GLO, @OFF_GAL, @OFF_BDS); - my ( @EPOCH, @EPOCH_CLK, @EPOCH_OFF_GLO, @EPOCH_OFF_GAL, @EPOCH_OFF_BDS ); + my ( @hlp, @N, @E, @U, %SATNUM, @TRP, @CLK, @OFF_GPS, @OFF_GLO, @OFF_GAL, @OFF_BDS); + my ( @EPOCH, @EPOCH_CLK, @EPOCH_OFF_GPS, @EPOCH_OFF_GLO, @EPOCH_OFF_GAL, @EPOCH_OFF_BDS ); my ( %AMB, %RES, %ELE, %ION, %BIA ); my ( $station, $lki, $sys, $sat, $amb ); @@ -588,14 +588,18 @@ push ( @CLK, $hlp[2] + $hlp[3] ); } + elsif ( $ln =~ /\bOFF_GPS\b/ ) { # 2015-08... OFF_GPS 52.6806 -3.8042 +- 9.0077 + push ( @EPOCH_OFF_GPS, $epochSec ); + push ( @OFF_GPS, $hlp[2] + $hlp[3] ); + } elsif ( $ln =~ /\bOFF_GLO\b/ ) { # 2015-08... OFF_GLO 52.6806 -3.8042 +- 9.0077 -# push ( @EPOCH_OFF_GLO, $epochSec ); + push ( @EPOCH_OFF_GLO, $epochSec ); push ( @OFF_GLO, $hlp[2] + $hlp[3] ); } elsif ( $ln =~ /\bOFF_GAL\b/ ) { # 2015-08... OFF_GAL 52.6806 -3.8042 +- 9.0077 -# push ( @EPOCH_OFF_GAL, $epochSec ); + push ( @EPOCH_OFF_GAL, $epochSec ); push ( @OFF_GAL, $hlp[2] + $hlp[3] ); } elsif ( $ln =~ /\bOFF_BDS\b/ ) { # 2015-08... OFF_BDS 52.6806 -3.8042 +- 9.0077 -# push ( @EPOCH_OFF_BDS, $epochSec ); + push ( @EPOCH_OFF_BDS, $epochSec ); push ( @OFF_BDS, $hlp[2] + $hlp[3] ); } @@ -624,7 +628,8 @@ if ( $nof_epochs != scalar @TRP ) { LOGDIE "number of epochs and TRP not equal\n" } if ( @CLK && scalar @EPOCH_CLK != scalar @CLK ) { LOGDIE "number of epochs and CLK not equal\n" } -# if ( @OFF_GLO && scalar @EPOCH_OFF_GLO != scalar @OFF_GLO ) { LOGDIE "number of epochs and OFF_GLO not equal\n" } -# if ( @OFF_GAL && scalar @EPOCH_OFF_GAL != scalar @OFF_GAL ) { LOGDIE "number of epochs and OFF_GAL not equal\n" } -# if ( @OFF_BDS && scalar @EPOCH_OFF_BDS != scalar @OFF_BDS ) { LOGDIE "number of epochs and OFF_BDS not equal\n" } + if ( @OFF_GPS && scalar @EPOCH_OFF_GPS != scalar @OFF_GPS ) { LOGDIE "number of epochs and OFF_GPS not equal\n" } + if ( @OFF_GLO && scalar @EPOCH_OFF_GLO != scalar @OFF_GLO ) { LOGDIE "number of epochs and OFF_GLO not equal\n" } + if ( @OFF_GAL && scalar @EPOCH_OFF_GAL != scalar @OFF_GAL ) { LOGDIE "number of epochs and OFF_GAL not equal\n" } + if ( @OFF_BDS && scalar @EPOCH_OFF_BDS != scalar @OFF_BDS ) { LOGDIE "number of epochs and OFF_BDS not equal\n" } if ( !$station ) { WARN "could not grep stationname from file: $file\n" } @@ -638,4 +643,5 @@ TRP => \@TRP, CLK => \@CLK, + OFF_GPS => \@OFF_GPS, OFF_GLO => \@OFF_GLO, OFF_GAL => \@OFF_GAL, Index: /trunk/BNC/scripts/pppPlot.pl =================================================================== --- /trunk/BNC/scripts/pppPlot.pl (revision 10387) +++ /trunk/BNC/scripts/pppPlot.pl (revision 10388) @@ -107,4 +107,5 @@ my ( $station, $file ) = Bnc::parsePPPLogfile( $file, $sampling ); my $EPOCH = $file->{'EPOCH'}; + my $EPOCH_OFF_GPS = $file->{'EPOCH_OFF_GPS'}; my $EPOCH_OFF_GLO = $file->{'EPOCH_OFF_GLO'}; my $EPOCH_OFF_GAL = $file->{'EPOCH_OFF_GAL'}; @@ -243,4 +244,37 @@ } + ######### OFF_GPS ##################### + if ( grep ( $_ eq "ALL", @plotTypes ) ) { + if ( scalar @{ $file->{'OFF_GPS'} } < 1 ) { + DEBUG "No OFF_GPS found"; + } + else { + DEBUG "Plot OFF_GPS"; + my $pngName = sprintf ( "%s_OFF_G.png", $station ); + my $chartOFF_GPS = newChart($station); + $chartOFF_GPS->set( output => $pngName ); + $chartOFF_GPS->set( ylabel => "Offset GPS [m]" ); + + $dataset = Chart::Gnuplot::DataSet->new( + xdata => $EPOCH_OFF_GPS, + ydata => $file->{'OFF_GPS'}, + title => "Receiver Offset GPS", + timefmt => '%s', + style => "linespoints", + ); + $chartOFF_GPS->plot2d($dataset); #system ("display $pngName&"); + $y = $y - $dy; + + if ( $y < 30 / mm ) { + $page = $pdf->page(); + $page->mediabox('A4'); + $y = $y0; + } + $png = $page->gfx(); + LOGDIE("could not find image file: $!\n") unless -e $pngName; + $png->image( $pdf->image_png($pngName), $x, $y, $width, $height ); + } + } + ######### OFF_GLO ##################### @@ -257,6 +291,5 @@ $dataset = Chart::Gnuplot::DataSet->new( - #xdata => $EPOCH_OFF_GLO, - xdata => $EPOCH, + xdata => $EPOCH_OFF_GLO, ydata => $file->{'OFF_GLO'}, title => "Receiver Offset GLONASS", @@ -291,6 +324,5 @@ $dataset = Chart::Gnuplot::DataSet->new( - #xdata => $EPOCH_OFF_GAL, - xdata => $EPOCH, + xdata => $EPOCH_OFF_GAL, ydata => $file->{'OFF_GAL'}, title => "Receiver Offset Galileo", @@ -325,6 +357,5 @@ $dataset = Chart::Gnuplot::DataSet->new( - #xdata => $EPOCH_OFF_BDS, - xdata => $EPOCH, + xdata => $EPOCH_OFF_BDS, ydata => $file->{'OFF_BDS'}, title => "Receiver Offset Beidou", Index: /trunk/BNC/src/PPP/pppClient.cpp =================================================================== --- /trunk/BNC/src/PPP/pppClient.cpp (revision 10387) +++ /trunk/BNC/src/PPP/pppClient.cpp (revision 10388) @@ -67,4 +67,5 @@ } } + _offGps = 0.0; _offGlo = 0.0; _offGal = 0.0; @@ -301,5 +302,5 @@ const t_pppSatObs* satObs = obsVector.at(ii); char sys = satObs->prn().system(); - if (satObs->isValid() && + if (satObs->isValid() && (!satObs->modelSet() || satObs->eleSat() >= _opt->_minEle) ) { ColumnVector rr = satObs->xc().Rows(1,3) - xyzc.Rows(1,3); @@ -337,4 +338,58 @@ return success; } + +// Compute A Priori Gps Clock Offset +////////////////////////////////////////////////////////////////////////////// +double t_pppClient::cmpOffGps(vector& obsVector) { + + t_lc::type tLC = t_lc::dummy; + double offGps = 0.0; + + if (_opt->useSystem('G')) { + while (obsVector.size() > 0) { + offGps = 0.0; + double maxRes = 0.0; + int maxResIndex = -1; + unsigned nObs = 0; + t_prn maxResPrn; + for (unsigned ii = 0; ii < obsVector.size(); ii++) { + const t_pppSatObs* satObs = obsVector.at(ii); + if (satObs->prn().system() == 'G') { + if (tLC == t_lc::dummy) { + tLC = satObs->isValid(t_lc::cIF) ? t_lc::cIF : t_lc::c1; + } + if (satObs->isValid(tLC) && (!satObs->modelSet() || satObs->eleSat() >= _opt->_minEle)) { + double ll = satObs->obsValue(tLC) - satObs->cmpValue(tLC); + ++nObs; + offGps += ll; + if (fabs(ll) > fabs(maxRes)) { + maxRes = ll; + maxResIndex = ii; + maxResPrn = satObs->prn(); + } + } + } + } + + if (nObs > 0) { + offGps = offGps / nObs; + } + else { + offGps = 0.0; + } + + if (fabs(maxRes) > 100.0) { + LOG << "t_pppClient::cmpOffGps outlier " << maxResPrn.toString() << " " << maxRes << endl; + delete obsVector.at(maxResIndex); + obsVector.erase(obsVector.begin() + maxResIndex); + } + else { + break; + } + } + } + return offGps; +} + // Compute A Priori Glonass Clock Offset @@ -443,5 +498,4 @@ return offGal; } - // Compute A Priori BDS Clock Offset ////////////////////////////////////////////////////////////////////////////// @@ -638,4 +692,5 @@ } + _offGps = cmpOffGps(_obsRover); _offGlo = cmpOffGlo(_obsRover); _offGal = cmpOffGal(_obsRover); Index: /trunk/BNC/src/PPP/pppClient.h =================================================================== --- /trunk/BNC/src/PPP/pppClient.h (revision 10387) +++ /trunk/BNC/src/PPP/pppClient.h (revision 10388) @@ -36,7 +36,12 @@ const bncAntex* antex() const {return _antex;} const t_pppStation* staRover() const {return _staRover;} + double offGps() const {return _offGps;} double offGlo() const {return _offGlo;} double offGal() const {return _offGal;} double offBds() const {return _offBds;} + void resetOffGps() {_offGps = 0.0;} + void resetOffGlo() {_offGlo = 0.0;} + void resetOffGal() {_offGal = 0.0;} + void resetOffBds() {_offBds = 0.0;} @@ -61,4 +66,5 @@ t_irc cmpBancroft(const bncTime& epoTime, std::vector& obsVector, ColumnVector& xyzc, bool print); + double cmpOffGps(std::vector& obsVector); double cmpOffGlo(std::vector& obsVector); double cmpOffGal(std::vector& obsVector); @@ -73,4 +79,5 @@ bncAntex* _antex; t_pppFilter* _filter; + double _offGps; double _offGlo; double _offGal; Index: /trunk/BNC/src/PPP/pppFilter.cpp =================================================================== --- /trunk/BNC/src/PPP/pppFilter.cpp (revision 10387) +++ /trunk/BNC/src/PPP/pppFilter.cpp (revision 10388) @@ -310,4 +310,6 @@ string epoTimeStr = string(_epoTime); const vector ¶ms = _parlist->params(); + int maxNumberOfReset = (2.0 * obsVector.size()) - 2.0; + int numberOfReset = 0; for (unsigned ii = 0; ii < LCs.size(); ii++) { @@ -347,5 +349,5 @@ // Check Pre-Fit Residuals - /* ----------------------- + // ----------------------- else { ColumnVector AA(params.size()); @@ -358,9 +360,12 @@ if (fabs(vv) > SLIP) { + numberOfReset++; LOG << epoTimeStr << " cycle slip detected " << t_lc::toString(tLC) << ' ' << obs->prn().toString() << ' ' << setw(8) << setprecision(4) << vv << endl; - resetAmb(obs->prn(), obsVector, tLC); - } - }*/ + if (numberOfReset < maxNumberOfReset) { + resetAmb(obs->prn(), obsVector, tLC); + } + } + } } } Index: /trunk/BNC/src/PPP/pppParlist.cpp =================================================================== --- /trunk/BNC/src/PPP/pppParlist.cpp (revision 10387) +++ /trunk/BNC/src/PPP/pppParlist.cpp (revision 10388) @@ -74,4 +74,8 @@ const t_pppSatObs* obs = obsVector->at(ii); if (obs->prn() == _prn) { + double offGps = 0.0; + if (_prn.system() == 'G' && tLC != t_lc::MW) { + offGps = PPP_CLIENT->offGps(); + } double offGlo = 0.0; if (_prn.system() == 'R' && tLC != t_lc::MW) { @@ -86,9 +90,14 @@ offBds = PPP_CLIENT->offBds(); } - _x0 = floor((obs->obsValue(tLC) - offGlo - offGal - offBds - obs->cmpValue(tLC)) / obs->lambda(tLC) + 0.5); + _x0 = floor((obs->obsValue(tLC) - offGps - offGlo - offGal - offBds - obs->cmpValue(tLC)) / obs->lambda(tLC) + 0.5); break; } } } + break; + case offGps: + _epoSpec = true; + _sigma0 = OPT->_aprSigClkOff; + _x0 = PPP_CLIENT->offGps(); break; case offGlo: @@ -171,4 +180,7 @@ if (tLC == t_lc::GIM) {return 0.0;} return 1.0; + case offGps: + if (tLC == t_lc::GIM) {return 0.0;} + return (obs->prn().system() == 'G') ? 1.0 : 0.0; case offGlo: if (tLC == t_lc::GIM) {return 0.0;} @@ -298,4 +310,7 @@ case rClk: ss << "REC_CLK "; + break; + case offGps: + ss << "OFF_GPS "; break; case offGlo: @@ -401,4 +416,13 @@ } + // check which systems have observations + // ------------------------------------- + _usedSystems['G'] = _usedSystems['R'] = _usedSystems['E'] = _usedSystems['C'] = 0; + for (unsigned jj = 0; jj < obsVector.size(); jj++) { + const t_pppSatObs* satObs = obsVector[jj]; + char sys = satObs->prn().system(); + _usedSystems[sys]++; + } + // Check whether parameters have observations // ------------------------------------------ @@ -435,27 +459,44 @@ required.push_back(new t_pppParam(t_pppParam::crdZ, t_prn(), t_lc::dummy)); - // Receiver Clock - // -------------- + // Receiver Clocks + // --------------- required.push_back(new t_pppParam(t_pppParam::rClk, t_prn(), t_lc::dummy)); // GLONASS Clock Offset // -------------------- - if ((OPT->useSystem('G') && OPT->useSystem('R')) || - (OPT->useSystem('E') && OPT->useSystem('R')) || - (OPT->useSystem('C') && OPT->useSystem('R')) ) { + if ( _usedSystems.value('R') && + (_usedSystems.value('G') || _usedSystems.value('E') || _usedSystems.value('C'))) { required.push_back(new t_pppParam(t_pppParam::offGlo, t_prn(), t_lc::dummy)); + } + else { + PPP_CLIENT->resetOffGlo(); } // Galileo Clock Offset // -------------------- - if (OPT->useSystem('G') && OPT->useSystem('E')) { + if (_usedSystems.value('E') && _usedSystems.value('G') && _usedSystems.value('G') >= OPT->_minObs) { required.push_back(new t_pppParam(t_pppParam::offGal, t_prn(), t_lc::dummy)); + } + else { + PPP_CLIENT->resetOffGal(); + } + + // GPS Clock Offset + // -------------------- + if (_usedSystems.value('E') && _usedSystems.value('G') && _usedSystems.value('G') < OPT->_minObs) { + required.push_back(new t_pppParam(t_pppParam::offGps, t_prn(), t_lc::dummy)); + } + else { + PPP_CLIENT->resetOffGps(); } // BDS Clock Offset // ---------------- - if ((OPT->useSystem('G') && OPT->useSystem('C')) || - (OPT->useSystem('E') && OPT->useSystem('C'))) { + if (_usedSystems.contains('C') && + (_usedSystems.contains('G') || _usedSystems.contains('E'))) { required.push_back(new t_pppParam(t_pppParam::offBds, t_prn(), t_lc::dummy)); + } + else { + PPP_CLIENT->resetOffBds(); } Index: /trunk/BNC/src/PPP/pppParlist.h =================================================================== --- /trunk/BNC/src/PPP/pppParlist.h (revision 10387) +++ /trunk/BNC/src/PPP/pppParlist.h (revision 10388) @@ -14,5 +14,5 @@ class t_pppParam { public: - enum e_type {crdX, crdY, crdZ, rClk, offGlo, offGal, offBds, trp, ion, amb, + enum e_type {crdX, crdY, crdZ, rClk, offGps, offGlo, offGal, offBds, trp, ion, amb, cBiasG1, cBiasR1, cBiasE1, cBiasC1, pBiasG1, pBiasR1, pBiasE1, pBiasC1, cBiasG2, cBiasR2, cBiasE2, cBiasC2, pBiasG2, pBiasR2, pBiasE2, pBiasC2}; @@ -104,4 +104,5 @@ const std::vector& params() const {return _params;} std::vector& params() {return _params;} + const QMap& usedSystems() const {return _usedSystems;} void printResult(const bncTime& epoTime, const SymmetricMatrix& QQ, const ColumnVector& xx) const; @@ -110,4 +111,5 @@ private: std::vector _params; + QMap _usedSystems; };