Index: /trunk/BNC/Example_Configs/10_PPPQuickStart.bnc =================================================================== --- /trunk/BNC/Example_Configs/10_PPPQuickStart.bnc (revision 11057) +++ /trunk/BNC/Example_Configs/10_PPPQuickStart.bnc (revision 11058) @@ -33,5 +33,5 @@ miscPort= miscScanRTCM=0 -mountPoints=//Example:Configs@euref-ip.net:2101/WTZR00DEU0 RTCM_3.3 DEU 49.14 12.88 no 2, //Example:Configs@products.igs-ip.net:2101/BCEP00BKG0 RTCM_3 DEU 50.09 8.66 no 2, //Example:Configs@products.igs-ip.net:2101/SSRA01CAS1 RTCM_3.1 CHN 40.07 116.28 no 2 +mountPoints=//Example:Configs@euref-ip.net:2101/WTZR00DEU0 RTCM_3.3 DEU 49.14 12.88 no 2, //Example:Configs@products.igs-ip.net:2101/BCEP00BKG0 RTCM_3 DEU 50.09 8.66 no 2, //Example:Configs@products.igs-ip.net:2101/SSRA00CNE1 RTCM_3.1 FRA 43.56 1.48 no 2 ntripVersion=2 onTheFlyInterval=no @@ -151,5 +151,5 @@ constraints=Ionosphere: pseudo-obs corrFile= -corrMount=SSRA00CAS1 +corrMount=SSRA00CNE1 corrWaitTime=5 crdFile=Input/BNC_IGS20.CRD Index: /trunk/BNC/Example_Configs/17_Combi.bnc =================================================================== --- /trunk/BNC/Example_Configs/17_Combi.bnc (revision 11057) +++ /trunk/BNC/Example_Configs/17_Combi.bnc (revision 11058) @@ -17,5 +17,5 @@ cmbNavic=0 cmbQzss=0 -cmbSampl=1 sec +cmbSampl=5 sec cmbSbas=0 cmbStreams="SSRA00BKG1 BKG 1.0 Computed ", "SSRA00CNE1 CNE 1.0 Computed " Index: /trunk/BNC/Example_Configs/18_CombiPPP.bnc =================================================================== --- /trunk/BNC/Example_Configs/18_CombiPPP.bnc (revision 11057) +++ /trunk/BNC/Example_Configs/18_CombiPPP.bnc (revision 11058) @@ -17,5 +17,5 @@ cmbNavic=0 cmbQzss=0 -cmbSampl=1 sec +cmbSampl=5 sec cmbSbas=0 cmbStreams="SSRA00BKG1 BKG 1.0 Computed ", "SSRA00CNE1 CNE 1.0 Computed " Index: /trunk/BNC/src/RTCM3/RTCM3coDecoder.cpp =================================================================== --- /trunk/BNC/src/RTCM3/RTCM3coDecoder.cpp (revision 11057) +++ /trunk/BNC/src/RTCM3/RTCM3coDecoder.cpp (revision 11058) @@ -1132,4 +1132,11 @@ // Set _lastTime close to currentTime // ---------------------------------- + if (_type == IGSssr) { + if (epoSecGPS != -1) _igsRefEpoSec = epoSecGPS; + else if (epoSecGal != -1) _igsRefEpoSec = epoSecGal; + else if (epoSecQzss != -1) _igsRefEpoSec = epoSecQzss; + else if (epoSecSbas != -1) _igsRefEpoSec = epoSecSbas; + } + if (epoSecGPS != -1) { _lastTime.set(currentWeek, epoSecGPS); @@ -1137,10 +1144,28 @@ else if (epoSecGlo != -1) { QDate date = dateAndTimeFromGPSweek(currentTime.gpsw(), currentTime.gpssec()).date(); + double leap = gnumleap(date.year(), date.month(), date.day()); // IGS-SSR: the SSR Epoch Time (IDF003, 20 bit) is given in one continuous // time scale for all GNSS, i.e. GPS seconds of week also for GLONASS - // (IGS SSR v1.00) - no conversion + // (IGS SSR v1.00) - no conversion. Old encoders (e.g. BNC before + // 2026-09-30) wrote UTC seconds of day instead: recognised by comparison + // with the GPS/Galileo epoch of the same stream and converted, with a + // warning, during a transition period. + if (_type == IGSssr && _igsRefEpoSec >= 0.0 && + igsLegacyEpoch(epoSecGlo, fmod(_igsRefEpoSec - leap + 7 * 86400.0, 86400.0), 86400.0)) { + double sow = floor(_igsRefEpoSec / 86400.0) * 86400.0 + epoSecGlo + leap; + if (sow - _igsRefEpoSec > 43200.0) sow -= 86400.0; + else if (_igsRefEpoSec - sow > 43200.0) sow += 86400.0; + epoSecGlo = sow; + bncTime& warned = _legacyEpochWarned['R']; + if (!warned.valid() || currentTime - warned >= 3600.0) { + warned = currentTime; + emit newMessage(_staID.toLatin1() + ": IGS-SSR GLONASS epoch time in outdated convention " + "(UTC seconds of day instead of GPS seconds of week, IGS SSR v1.00 IDF003) " + "detected and converted - the SSR encoder should be updated", true); + } + } if (_type == RTCMssr || _type == RTCMnewssr) { // Same GLONASS epoch wire convention (UTC+3h) as RTCMssr - epoSecGlo = epoSecGlo - 3 * 3600 + gnumleap(date.year(), date.month(), date.day()); + epoSecGlo = epoSecGlo - 3 * 3600 + leap; } _lastTime.set(currentWeek, epoSecGlo); @@ -1156,6 +1181,19 @@ } else if (epoSecBds != -1) { - // RTCM-SSR: BDT; IGS-SSR: GPS time scale for all GNSS (IDF003) - if (_type != IGSssr) { + // RTCM-SSR: BDT; IGS-SSR: GPS time scale for all GNSS (IDF003). Old + // IGS-SSR encoders wrote BDT: recognised and converted as for GLONASS. + bool bdt = (_type != IGSssr); + if (_type == IGSssr && _igsRefEpoSec >= 0.0 && + igsLegacyEpoch(epoSecBds, fmod(_igsRefEpoSec - 14.0 + 604800.0, 604800.0), 604800.0)) { + bdt = true; + bncTime& warned = _legacyEpochWarned['C']; + if (!warned.valid() || currentTime - warned >= 3600.0) { + warned = currentTime; + emit newMessage(_staID.toLatin1() + ": IGS-SSR BDS epoch time in outdated convention " + "(BDT instead of GPS seconds of week, IGS SSR v1.00 IDF003) " + "detected and converted - the SSR encoder should be updated", true); + } + } + if (bdt) { epoSecBds += 14.0; if (epoSecBds > 604800.0) { @@ -1182,4 +1220,18 @@ } } +} + +// IGS-SSR: true if value (seconds, periodic with period) matches the epoch +// expected in the outdated convention (legacyRef) better than the GPS time +// reference of the stream, within the spread of epochs of one update cycle +//////////////////////////////////////////////////////////////////////////// +bool RTCM3coDecoder::igsLegacyEpoch(double value, double legacyRef, double period) const { + auto dist = [period](double a, double b) { + double d = fmod(fabs(a - b), period); + return std::min(d, period - d); + }; + double dSpec = dist(value, fmod(_igsRefEpoSec, period)); + double dLegacy = dist(value, legacyRef); + return dLegacy < dSpec && dLegacy <= 30.0; } Index: /trunk/BNC/src/RTCM3/RTCM3coDecoder.h =================================================================== --- /trunk/BNC/src/RTCM3/RTCM3coDecoder.h (revision 11057) +++ /trunk/BNC/src/RTCM3/RTCM3coDecoder.h (revision 11058) @@ -82,4 +82,10 @@ void checkPrnRange(const t_prn& prn); QSet _prnOutOfRangeLogged; + // IGS-SSR: latest GPS/Galileo/QZSS/SBAS epoch (GPS seconds of week) of + // this stream, reference to detect GLONASS/BDS epochs in the outdated + // convention of old encoders (GLONASS: UTC seconds of day, BDS: BDT) + double _igsRefEpoSec = -1.0; + QMap _legacyEpochWarned; // warned at most once per hour + bool igsLegacyEpoch(double value, double legacyRef, double period) const; Index: /trunk/BNC/src/bnchelp.html =================================================================== --- /trunk/BNC/src/bnchelp.html (revision 11057) +++ /trunk/BNC/src/bnchelp.html (revision 11058) @@ -8105,4 +8105,10 @@ Specify a combination sampling interval. Orbit and clock corrections will be produced following that interval. A value of 10 sec may be an appropriate choice.

+

+ As the incoming corrections of the individual providers arrive with different latencies, BNC combines an epoch as + soon as each provider has delivered corrections for a later epoch, i.e. has completed the epoch. The combined + product is thus only delayed by the latency of the slowest provider. If a provider is late or missing, the epoch is + combined without it at the latest when the newest incoming corrections are 10 seconds (or the sampling interval, + if larger) younger than the epoch. Corrections arriving after an epoch was combined are ignored.

@@ -9108,4 +9114,11 @@ A value of 10 sec may be appropriate. A value of zero '0' tells BNC to store all available samples into Clock RINEX files. +

+

+ SP3 and Clock RINEX epochs are multiples of the sampling interval. Corrections whose epochs are not on this grid, + e.g. Galileo HAS corrections at ...:02, ...:12 seconds, are written at the preceding grid epoch if they refer to an + epoch at most 5 seconds (and less than half the sampling interval) after it. The satellite clock is then propagated + to the grid epoch using its rate, the orbit using the broadcast orbit at the grid epoch plus the difference between + precise and broadcast orbit at the correction epoch. Corrections on the grid are always used unchanged.

Index: /trunk/BNC/src/bncoutf.h =================================================================== --- /trunk/BNC/src/bncoutf.h (revision 11057) +++ /trunk/BNC/src/bncoutf.h (revision 11058) @@ -16,4 +16,5 @@ QString agencyFromFileName(); bool isProductFile(); + int sampl() const {return _sampl;} protected: virtual t_irc reopen(int GPSweek, double GPSweeks); Index: /trunk/BNC/src/combination/bnccomb.cpp =================================================================== --- /trunk/BNC/src/combination/bnccomb.cpp (revision 11057) +++ /trunk/BNC/src/combination/bnccomb.cpp (revision 11058) @@ -535,4 +535,18 @@ +// True if every AC has delivered clock corrections for an epoch after epoTime +//////////////////////////////////////////////////////////////////////////// +bool bncComb::allACsBeyond(const bncTime& epoTime) const { + QListIterator itAC(_ACs); + while (itAC.hasNext()) { + const cmbAC* AC = itAC.next(); + bncTime acLast = _acLastClkTime.value(AC->name); + if (!acLast.valid() || !(acLast > epoTime)) { + return false; + } + } + return true; +} + // Process clock corrections //////////////////////////////////////////////////////////////////////////// @@ -552,4 +566,19 @@ if (!_cmbSysPrn.contains(sys)){ continue; + } + + // Remember the latest correction epoch of each AC (whether or not the + // correction is used), to know when an AC has finished an epoch + // -------------------------------------------------------------------- + QListIterator itAC(_ACs); + while (itAC.hasNext()) { + cmbAC* AC = itAC.next(); + if (AC->mountPoint == QString(newClk._staID.c_str())) { + bncTime& acLast = _acLastClkTime[AC->name]; + if (!acLast.valid() || newClk._time > acLast) { + acLast = newClk._time; + } + break; + } } @@ -665,7 +694,12 @@ } } - // Process the front epoch - // ----------------------- - if (epoTime < (_lastClkCorrTime - outWait)) { + // Process the front epoch as soon as all ACs have delivered corrections + // for a later epoch (i.e. have finished this one), at the latest when + // the newest correction is maxWait newer (an AC may be late or missing). + // Waiting only outWait after the newest correction lost the corrections + // of ACs with a few seconds more latency ("old correction"). + // ----------------------------------------------------------------------- + const double maxWait = max(10.0, outWait); + if (allACsBeyond(epoTime) || epoTime < (_lastClkCorrTime - maxWait)) { _resTime = epoTime; processEpoch(_resTime, clkCorrVec); @@ -834,4 +868,24 @@ codeBiasesRefSig[frqType] = bias._value; } + } + } + // Code biases present, but not for both reference signals (e.g. + // Galileo HAS): the AC clocks cannot be referred to the reference + // signals and are used as they are - warn (at most once per hour + // per AC and system) + if (codeBiasesRefSig.size() < 2) { + QString key = acName + sys; + bncTime& warned = _refBiasWarned[key]; + if (!warned.valid() || _newCorr->_time - warned >= 3600.0) { + warned = _newCorr->_time; + QString refSig; + for (unsigned ii = 1; ii < cmbRefSig::cIF; ii++) { + t_frequency::type frqType = cmbRefSig::toFreq(sys, static_cast(ii)); + refSig += QString(" C%1%2").arg(t_frequency::toString(frqType)[1]) + .arg(cmbRefSig::toAttrib(sys, static_cast(ii))); + } + emit newMessage("bncComb: " + acName.toLatin1() + " provides no code biases for the " + + QByteArray(1, sys) + " reference signals" + refSig.toLatin1() + + " - its clocks are combined without referring them to these signals", false); } } @@ -1398,5 +1452,10 @@ double dC_apriori_corr = dC_orb;// + dC_att; can be added if useful values are provided by the ACs - ll(iObs) = ((corr->_clkCorr._dClk * t_CST::c) - dC_apriori_corr - corr->_satCodeBiasIF) - DotProduct(AA.Row(iObs), x0); + // Clock referred to the ionosphere-free combination of the reference + // signals: dt - (alpha b1 + beta b2) with biases in Bias-SINEX sign + // (subtracted from observations, Banville et al. 2020, Eq. 1). The code + // biases taken from the SSR streams have the opposite sign (added to + // observations), hence + _satCodeBiasIF. + ll(iObs) = ((corr->_clkCorr._dClk * t_CST::c) - dC_apriori_corr + corr->_satCodeBiasIF) - DotProduct(AA.Row(iObs), x0); PP(iObs, iObs) *= 1.0 / (corr->_weightFactor * corr->_weightFactor); Index: /trunk/BNC/src/combination/bnccomb.h =================================================================== --- /trunk/BNC/src/combination/bnccomb.h (revision 11057) +++ /trunk/BNC/src/combination/bnccomb.h (revision 11058) @@ -267,4 +267,7 @@ std::deque _epoClkData; bncTime _lastClkCorrTime; + QMap _acLastClkTime; // latest correction epoch per AC + QMap _refBiasWarned; // last warning per AC and system + bool allACsBeyond(const bncTime& epoTime) const; bncTime _resTime; cmbCorr* _newCorr; Index: /trunk/BNC/src/upload/bncrtnetuploadcaster.cpp =================================================================== --- /trunk/BNC/src/upload/bncrtnetuploadcaster.cpp (revision 11057) +++ /trunk/BNC/src/upload/bncrtnetuploadcaster.cpp (revision 11058) @@ -548,7 +548,19 @@ // Make sure the clock messages refer to same IOD as orbit messages // ---------------------------------------------------------------- + // The ephemeris is selected at the first correction epoch of each + // sampling interval: the grid epoch itself for corrections on the + // sampling grid, otherwise (e.g. Galileo HAS at ...:02, ...:12) + // the first epoch after it, which previously never selected one if (_usedEph) { - if (fmod(epoTime.gpssec(), _samplRtcmEphCorr) == 0.0) { + bool select = (fmod(epoTime.gpssec(), _samplRtcmEphCorr) == 0.0); + if (!select && _samplRtcmEphCorr > 0.0) { + select = !_usedEphEpoch.valid() || epoTime == _usedEphEpoch || + epoTime.gpsw() != _usedEphEpoch.gpsw() || + floor(epoTime.gpssec() / _samplRtcmEphCorr) != + floor(_usedEphEpoch.gpssec() / _samplRtcmEphCorr); + } + if (select) { (*_usedEph)[prnInternalStr] = eph; + _usedEphEpoch = epoTime; } else { @@ -1170,5 +1182,15 @@ double clkRnxAcc = rtnClk[2] / t_CST::c; // [s/s² = -/s] - if (_rnx) { + // Corrections on epochs off the output sampling grid (e.g. Galileo HAS + // at ...:02, ...:12) are written at the nearest grid epoch, see gridEpoch; + // the clock is propagated with its rate, the orbit via the broadcast + // orbit at the grid epoch plus the (slowly varying) difference + // precise - broadcast at the correction epoch + // ----------------------------------------------------------------------- + bncTime tt(GPSweek, GPSweeks); + bncTime tg; + double clkRate = (clkRnxRate != 0.0) ? clkRnxRate : xB(5); + + if (_rnx && gridEpoch(_rnx->sampl(), tt, tg, _rnxGridEpo, _rnxGridSrc)) { double clkRnxSig, clkRnxRateSig, clkRnxAccSig; int s = rtnClkSig.size(); @@ -1190,11 +1212,50 @@ break; } - _rnx->write(GPSweek, GPSweeks, prn, clkRnx, clkRnxRate, clkRnxAcc, - clkRnxSig, clkRnxRateSig, clkRnxAccSig); - } - if (_sp3) { - _sp3->write(GPSweek, GPSweeks, prn, rtnCoM, clkRnx, rtnVel, clkRnxRate); + _rnx->write(tg.gpsw(), tg.gpssec(), prn, clkRnx + clkRate * (tg - tt), + clkRnxRate, clkRnxAcc, clkRnxSig, clkRnxRateSig, clkRnxAccSig); + } + if (_sp3 && gridEpoch(_sp3->sampl(), tt, tg, _sp3GridEpo, _sp3GridSrc)) { + ColumnVector xyz = rtnCoM; + ColumnVector vel = rtnVel; + if (!(tg == tt)) { + ColumnVector xBg(6); + ColumnVector vBg(3); + if (eph->getCrd(tg, xBg, vBg, false) != success) { + return success; + } + xyz = xBg.Rows(1, 3) + (rtnCoM - xB.Rows(1, 3)); + vel = vBg + (rtnVel - vB); + } + _sp3->write(tg.gpsw(), tg.gpssec(), prn, xyz, clkRnx + clkRate * (tg - tt), + vel, clkRnxRate); } return success; +} + +// Output epoch for SP3/CLK files with sampling interval sampl [s]: the +// correction epoch itself if it is on the sampling grid, otherwise the +// preceding grid epoch if at most 5 s (and less than half the sampling +// interval) earlier. Only preceding grid epochs are used, so a correction +// epoch on the grid always arrives first and takes precedence. Each grid +// epoch is derived from one correction epoch only. +//////////////////////////////////////////////////////////////////////////// +bool bncRtnetUploadCaster::gridEpoch(int sampl, const bncTime& tt, bncTime& tg, + bncTime& lastGrid, bncTime& lastSrc) const { + if (sampl <= 0) { + tg = tt; + return true; + } + double sec = tt.gpssec(); + double back = sec - floor(sec / sampl) * sampl; // [0, sampl) + if (back > 1.e-6 && (back > 5.0 || back >= 0.5 * sampl)) { + return false; + } + tg = (back <= 1.e-6) ? tt : tt - back; + if (lastGrid.valid() && tg == lastGrid && !(tt == lastSrc)) { + return false; + } + lastGrid = tg; + lastSrc = tt; + return true; } Index: /trunk/BNC/src/upload/bncrtnetuploadcaster.h =================================================================== --- /trunk/BNC/src/upload/bncrtnetuploadcaster.h (revision 11057) +++ /trunk/BNC/src/upload/bncrtnetuploadcaster.h (revision 11058) @@ -103,4 +103,9 @@ bncClockRinex* _rnx; bncSP3* _sp3; + // grid epoch last written to SP3/CLK and the epoch it was derived from + bncTime _sp3GridEpo, _sp3GridSrc, _rnxGridEpo, _rnxGridSrc; + bncTime _usedEphEpoch; // epoch the used ephemerides were last selected + bool gridEpoch(int sampl, const bncTime& tt, bncTime& tg, + bncTime& lastGrid, bncTime& lastSrc) const; bncBiasSinex* _bsx; QMap* _usedEph; Index: /trunk/BNC/tests/README.md =================================================================== --- /trunk/BNC/tests/README.md (revision 11057) +++ /trunk/BNC/tests/README.md (revision 11058) @@ -46,9 +46,17 @@ `replay_rinex_eph` (02), `replay_feed_engine` (07, IP port disabled) and `replay_save_sp3` (14, with SSRA00CNE1 instead of SSRA00BKG1). +* `replay_save_sp3_has`: 14 with Galileo HAS (SSRA00EUH0) only; its epochs are + off the SP3/CLK sampling grid and are written at the preceding grid epoch. * Real-time PPP replays (`ppp_rt_*`) from 11 min slices of the CNES and WHU captures (`data/wtzr_cnes_11min.raw`, `data/wtzr_whu_11min.raw`): - `ppp_rt_float` (08), `ppp_rt_quickstart` (10, with SSRA00CNE1 instead of - SSRA01CAS1, NMEA port disabled), `ppp_rt_ar_cnes` (29) and `ppp_rt_ar_whu` + `ppp_rt_float` (08), `ppp_rt_quickstart` (10, NMEA port disabled), `ppp_rt_ar_cnes` (29) and `ppp_rt_ar_whu` (30), compared exactly. +* `ppp_rt_net` (09, two stations) and `ppp_rt_has` (13, Galileo HAS via + SSRA00EUH0) from `data/merged_260930_11min.raw`, a slice of two simultaneous + recordings merged by time (BCEP00BKG0, WTZR00DEU0, FFMJ00DEU0, SSRA00CNE1, + SSRA00EUH0). +* Combination (`combi`, `combi_ppp`; 17 and 18) from `data/combi_260930_11min.raw` + (SSRA00BKG1, SSRA00CNE1, BCEP00BKG0, FFMJ00DEU0), upload to the caster + disabled, SP3/CLK written to files. * SSR decoding (`ssr_*`): 3 min slices of real streams replayed through `03_BrdcCorr.bnc`. They compare the `.ssr` files and the `BNC_DEBUG_SSR` Index: /trunk/BNC/tests/cases/ppp_rt_quickstart/case.json =================================================================== --- /trunk/BNC/tests/cases/ppp_rt_quickstart/case.json (revision 11057) +++ /trunk/BNC/tests/cases/ppp_rt_quickstart/case.json (revision 11058) @@ -1,4 +1,4 @@ { - "description": "10_PPPQuickStart replay: Quick-Start PPP with SSRA00CNE1 instead of SSRA01CAS1, NMEA port disabled, 11 min", + "description": "10_PPPQuickStart replay: Quick-Start PPP with SSRA00CNE1, NMEA port disabled, 11 min", "config": "Example_Configs/10_PPPQuickStart.bnc", "raw": "data/wtzr_cnes_11min.raw", @@ -9,7 +9,4 @@ "rawOutFile", "", - "-key", - "PPP/corrMount", - "SSRA00CNE1", "-key", "PPP/staTable",