PhaseMinusCodeCycleSlipDetector.java
- /* Copyright 2002-2024 CS GROUP
- * Licensed to CS GROUP (CS) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * CS licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
- package org.orekit.estimation.measurements.gnss;
- import java.util.ArrayList;
- import java.util.List;
- import java.util.Map;
- import org.hipparchus.fitting.PolynomialCurveFitter;
- import org.hipparchus.fitting.WeightedObservedPoint;
- import org.hipparchus.util.FastMath;
- import org.orekit.files.rinex.observation.ObservationData;
- import org.orekit.files.rinex.observation.ObservationDataSet;
- import org.orekit.gnss.Frequency;
- import org.orekit.gnss.MeasurementType;
- import org.orekit.gnss.SatelliteSystem;
- import org.orekit.time.AbsoluteDate;
- import org.orekit.utils.Constants;
- /**
- * Phase minus code cycle slip detectors.
- * The detector is based the algorithm given in <a
- * href="https://gssc.esa.int/navipedia/index.php/Examples_of_single_frequency_Cycle-Slip_Detectors">
- * Examples of single frequency Cycle-Slip Detectors</a> by Zornoza and M. Hernández-Pajares. Within this class
- * a polynomial is used to smooth the data. We consider a cycle_slip occurring if the current measurement is too
- * far from the one predicted with the polynomial (algorithm 1 on Navipedia).
- * <p>
- * For building the detector, one should give a threshold and a gap time limit.
- * After construction of the detectors, one can have access to a List of CycleData. Each CycleDate represents
- * a link between the station (define by the RINEX file) and a satellite at a specific frequency. For each cycle data,
- * one has access to the begin and end of availability, and a sorted set which contains all the date at which
- * cycle-slip have been detected
- * </p>
- * @author David Soulard
- * @since 10.2
- */
- public class PhaseMinusCodeCycleSlipDetector extends AbstractCycleSlipDetector {
- /** Mega Hertz to Hertz conversion factor. */
- private static final double MHZ_TO_HZ = 1.0e6;
- /** Order of the polynomial used for fitting. */
- private final int order;
- /** Threshold above which cycle-slip occurs. */
- private final double threshold;
- /** Polynomial single frequency cycle-slip detector Constructor.
- * @param dt time gap threshold between two consecutive measurement (if time between two consecutive measurement is greater than dt, a cycle slip is declared)
- * @param threshold threshold above which cycle-slip occurs
- * @param n number of measurement before starting
- * @param order polynomial order
- */
- public PhaseMinusCodeCycleSlipDetector(final double dt, final double threshold,
- final int n, final int order) {
- super(dt, n);
- this.threshold = threshold;
- this.order = order;
- }
- /** {@inheritDoc} */
- @Override
- protected void manageData(final ObservationDataSet observation) {
- // Extract observation data
- final SatelliteSystem system = observation.getSatellite().getSystem();
- final int prn = observation.getSatellite().getPRN();
- final AbsoluteDate date = observation.getDate();
- // Initialize list of measurements
- final List<ObservationData> pseudoRanges = new ArrayList<>();
- final List<ObservationData> phases = new ArrayList<>();
- // Loop on observation data to fill lists
- for (final ObservationData od : observation.getObservationData()) {
- if (!Double.isNaN(od.getValue())) {
- if (od.getObservationType().getMeasurementType() == MeasurementType.PSEUDO_RANGE) {
- pseudoRanges.add(od);
- } else if (od.getObservationType().getMeasurementType() == MeasurementType.CARRIER_PHASE) {
- phases.add(od);
- }
- }
- }
- // Loop on phase measurements
- for (final ObservationData phase : phases) {
- // Loop on range measurement
- for (final ObservationData pseudoRange : pseudoRanges) {
- // Change unit of phase measurement
- final double frequency = phase.getObservationType().getFrequency(system).getMHzFrequency() * MHZ_TO_HZ;
- final double cOverF = Constants.SPEED_OF_LIGHT / frequency;
- final ObservationData phaseInMeters = new ObservationData(phase.getObservationType(),
- cOverF * phase.getValue(),
- phase.getLossOfLockIndicator(),
- phase.getSignalStrength());
- // Check if measurement frequencies are the same
- if (phase.getObservationType().getFrequency(system) == pseudoRange.getObservationType().getFrequency(system)) {
- // Phase minus Code combination
- final PhaseMinusCodeCombination phaseMinusCode = MeasurementCombinationFactory.getPhaseMinusCodeCombination(system);
- final CombinedObservationData cod = phaseMinusCode.combine(phaseInMeters, pseudoRange);
- final String nameSat = setName(prn, observation.getSatellite().getSystem());
- // Check for cycle-slip detection
- final boolean slip = cycleSlipDetection(nameSat, date, cod.getValue(), phase.getObservationType().getFrequency(system));
- if (!slip) {
- // Update cycle slip data
- cycleSlipDataSet(nameSat, date, cod.getValue(), phase.getObservationType().getFrequency(system));
- }
- }
- }
- }
- }
- /**
- * Compute if there is a cycle slip at a specific date.
- * @param nameSat name of the satellite, on the predefined format (e.g. GPS - 07 for satellite 7 of GPS constellation)
- * @param currentDate the date at which we check if a cycle-slip occurs
- * @param phaseMinusCode phase measurement minus code measurement
- * @param frequency frequency used
- * @return true if a cycle slip has been detected.
- */
- private boolean cycleSlipDetection(final String nameSat, final AbsoluteDate currentDate,
- final double phaseMinusCode, final Frequency frequency) {
- // Access the cycle slip results to know if a cycle-slip already occurred
- final List<CycleSlipDetectorResults> data = getResults();
- final List<Map<Frequency, DataForDetection>> stuff = getStuffReference();
- // If a cycle-slip already occurred
- if (data != null) {
- // Loop on cycle-slip results
- for (CycleSlipDetectorResults resultPmC : data) {
- // Found the right cycle data
- if (resultPmC.getSatelliteName().compareTo(nameSat) == 0 && resultPmC.getCycleSlipMap().containsKey(frequency)) {
- final Map<Frequency, DataForDetection> values = stuff.get(data.indexOf(resultPmC));
- final DataForDetection v = values.get(frequency);
- // Check the time gap condition
- if (FastMath.abs(currentDate.durationFrom(v.getFiguresReference()[v.getWrite()].getDate())) > getMaxTimeBeetween2Measurement()) {
- resultPmC.addCycleSlipDate(frequency, currentDate);
- v.resetFigures( new SlipComputationData[getMinMeasurementNumber()], phaseMinusCode, currentDate);
- resultPmC.setDate(frequency, currentDate);
- return true;
- }
- // Compute the fitting polynomial if there are enough measurement since last cycle-slip
- if (v.getCanBeComputed() >= getMinMeasurementNumber()) {
- final List<WeightedObservedPoint> xy = new ArrayList<>();
- for (int i = 0; i < getMinMeasurementNumber(); i++) {
- final SlipComputationData current = v.getFiguresReference()[i];
- xy.add(new WeightedObservedPoint(1.0, current.getDate().durationFrom(currentDate),
- current.getValue()));
- }
- final PolynomialCurveFitter fitting = PolynomialCurveFitter.create(order);
- // Check if there is a cycle_slip
- if (FastMath.abs(fitting.fit(xy)[0] - phaseMinusCode) > threshold) {
- resultPmC.addCycleSlipDate(frequency, currentDate);
- v.resetFigures( new SlipComputationData[getMinMeasurementNumber()], phaseMinusCode, currentDate);
- resultPmC.setDate(frequency, currentDate);
- return true;
- }
- } else {
- break;
- }
- }
- }
- }
- // No cycle-slip
- return false;
- }
- }