Search in sources :

Example 6 with RangeTroposphericDelayModifier

use of org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier in project Orekit by CS-SI.

the class OrbitDeterminationTest method createStationsData.

/**
 * Set up stations.
 * @param parser input file parser
 * @param body central body
 * @return name to station data map
 * @exception OrekitException if some frame transforms cannot be computed
 * @throws NoSuchElementException if input parameters are missing
 */
private Map<String, StationData> createStationsData(final KeyValueFileParser<ParameterKey> parser, final OneAxisEllipsoid body) throws OrekitException, NoSuchElementException {
    final Map<String, StationData> stations = new HashMap<String, StationData>();
    final String[] stationNames = parser.getStringArray(ParameterKey.GROUND_STATION_NAME);
    final double[] stationLatitudes = parser.getAngleArray(ParameterKey.GROUND_STATION_LATITUDE);
    final double[] stationLongitudes = parser.getAngleArray(ParameterKey.GROUND_STATION_LONGITUDE);
    final double[] stationAltitudes = parser.getDoubleArray(ParameterKey.GROUND_STATION_ALTITUDE);
    final boolean[] stationPositionEstimated = parser.getBooleanArray(ParameterKey.GROUND_STATION_POSITION_ESTIMATED);
    final double[] stationRangeSigma = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_SIGMA);
    final double[] stationRangeBias = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_BIAS);
    final double[] stationRangeBiasMin = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_BIAS_MIN);
    final double[] stationRangeBiasMax = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_BIAS_MAX);
    final boolean[] stationRangeBiasEstimated = parser.getBooleanArray(ParameterKey.GROUND_STATION_RANGE_BIAS_ESTIMATED);
    final double[] stationRangeRateSigma = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_RATE_SIGMA);
    final double[] stationRangeRateBias = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_RATE_BIAS);
    final double[] stationRangeRateBiasMin = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_RATE_BIAS_MIN);
    final double[] stationRangeRateBiasMax = parser.getDoubleArray(ParameterKey.GROUND_STATION_RANGE_RATE_BIAS_MAX);
    final boolean[] stationRangeRateBiasEstimated = parser.getBooleanArray(ParameterKey.GROUND_STATION_RANGE_RATE_BIAS_ESTIMATED);
    final double[] stationAzimuthSigma = parser.getAngleArray(ParameterKey.GROUND_STATION_AZIMUTH_SIGMA);
    final double[] stationAzimuthBias = parser.getAngleArray(ParameterKey.GROUND_STATION_AZIMUTH_BIAS);
    final double[] stationAzimuthBiasMin = parser.getAngleArray(ParameterKey.GROUND_STATION_AZIMUTH_BIAS_MIN);
    final double[] stationAzimuthBiasMax = parser.getAngleArray(ParameterKey.GROUND_STATION_AZIMUTH_BIAS_MAX);
    final double[] stationElevationSigma = parser.getAngleArray(ParameterKey.GROUND_STATION_ELEVATION_SIGMA);
    final double[] stationElevationBias = parser.getAngleArray(ParameterKey.GROUND_STATION_ELEVATION_BIAS);
    final double[] stationElevationBiasMin = parser.getAngleArray(ParameterKey.GROUND_STATION_ELEVATION_BIAS_MIN);
    final double[] stationElevationBiasMax = parser.getAngleArray(ParameterKey.GROUND_STATION_ELEVATION_BIAS_MAX);
    final boolean[] stationAzElBiasesEstimated = parser.getBooleanArray(ParameterKey.GROUND_STATION_AZ_EL_BIASES_ESTIMATED);
    final boolean[] stationElevationRefraction = parser.getBooleanArray(ParameterKey.GROUND_STATION_ELEVATION_REFRACTION_CORRECTION);
    final boolean[] stationRangeTropospheric = parser.getBooleanArray(ParameterKey.GROUND_STATION_RANGE_TROPOSPHERIC_CORRECTION);
    for (int i = 0; i < stationNames.length; ++i) {
        // the station itself
        final GeodeticPoint position = new GeodeticPoint(stationLatitudes[i], stationLongitudes[i], stationAltitudes[i]);
        final TopocentricFrame topo = new TopocentricFrame(body, position, stationNames[i]);
        final GroundStation station = new GroundStation(topo);
        station.getEastOffsetDriver().setSelected(stationPositionEstimated[i]);
        station.getNorthOffsetDriver().setSelected(stationPositionEstimated[i]);
        station.getZenithOffsetDriver().setSelected(stationPositionEstimated[i]);
        // range
        final double rangeSigma = stationRangeSigma[i];
        final Bias<Range> rangeBias;
        if (FastMath.abs(stationRangeBias[i]) >= Precision.SAFE_MIN || stationRangeBiasEstimated[i]) {
            rangeBias = new Bias<Range>(new String[] { stationNames[i] + "/range bias" }, new double[] { stationRangeBias[i] }, new double[] { rangeSigma }, new double[] { stationRangeBiasMin[i] }, new double[] { stationRangeBiasMax[i] });
            rangeBias.getParametersDrivers().get(0).setSelected(stationRangeBiasEstimated[i]);
        } else {
            // bias fixed to zero, we don't need to create a modifier for this
            rangeBias = null;
        }
        // range rate
        final double rangeRateSigma = stationRangeRateSigma[i];
        final Bias<RangeRate> rangeRateBias;
        if (FastMath.abs(stationRangeRateBias[i]) >= Precision.SAFE_MIN || stationRangeRateBiasEstimated[i]) {
            rangeRateBias = new Bias<RangeRate>(new String[] { stationNames[i] + "/range rate bias" }, new double[] { stationRangeRateBias[i] }, new double[] { rangeRateSigma }, new double[] { stationRangeRateBiasMin[i] }, new double[] { stationRangeRateBiasMax[i] });
            rangeRateBias.getParametersDrivers().get(0).setSelected(stationRangeRateBiasEstimated[i]);
        } else {
            // bias fixed to zero, we don't need to create a modifier for this
            rangeRateBias = null;
        }
        // angular biases
        final double[] azELSigma = new double[] { stationAzimuthSigma[i], stationElevationSigma[i] };
        final Bias<AngularAzEl> azELBias;
        if (FastMath.abs(stationAzimuthBias[i]) >= Precision.SAFE_MIN || FastMath.abs(stationElevationBias[i]) >= Precision.SAFE_MIN || stationAzElBiasesEstimated[i]) {
            azELBias = new Bias<AngularAzEl>(new String[] { stationNames[i] + "/az bias", stationNames[i] + "/el bias" }, new double[] { stationAzimuthBias[i], stationElevationBias[i] }, azELSigma, new double[] { stationAzimuthBiasMin[i], stationElevationBiasMin[i] }, new double[] { stationAzimuthBiasMax[i], stationElevationBiasMax[i] });
            azELBias.getParametersDrivers().get(0).setSelected(stationAzElBiasesEstimated[i]);
            azELBias.getParametersDrivers().get(1).setSelected(stationAzElBiasesEstimated[i]);
        } else {
            // bias fixed to zero, we don't need to create a modifier for this
            azELBias = null;
        }
        // Refraction correction
        final AngularRadioRefractionModifier refractionCorrection;
        if (stationElevationRefraction[i]) {
            final double altitude = station.getBaseFrame().getPoint().getAltitude();
            final AtmosphericRefractionModel refractionModel = new EarthITU453AtmosphereRefraction(altitude);
            refractionCorrection = new AngularRadioRefractionModifier(refractionModel);
        } else {
            refractionCorrection = null;
        }
        // Tropospheric correction
        final RangeTroposphericDelayModifier rangeTroposphericCorrection;
        if (stationRangeTropospheric[i]) {
            final SaastamoinenModel troposphericModel = SaastamoinenModel.getStandardModel();
            rangeTroposphericCorrection = new RangeTroposphericDelayModifier(troposphericModel);
        } else {
            rangeTroposphericCorrection = null;
        }
        stations.put(stationNames[i], new StationData(station, rangeSigma, rangeBias, rangeRateSigma, rangeRateBias, azELSigma, azELBias, refractionCorrection, rangeTroposphericCorrection));
    }
    return stations;
}
Also used : GroundStation(org.orekit.estimation.measurements.GroundStation) HashMap(java.util.HashMap) EarthITU453AtmosphereRefraction(org.orekit.models.earth.EarthITU453AtmosphereRefraction) TopocentricFrame(org.orekit.frames.TopocentricFrame) AngularRadioRefractionModifier(org.orekit.estimation.measurements.modifiers.AngularRadioRefractionModifier) AtmosphericRefractionModel(org.orekit.models.AtmosphericRefractionModel) Range(org.orekit.estimation.measurements.Range) RangeTroposphericDelayModifier(org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier) GeodeticPoint(org.orekit.bodies.GeodeticPoint) RangeRate(org.orekit.estimation.measurements.RangeRate) GeodeticPoint(org.orekit.bodies.GeodeticPoint) AngularAzEl(org.orekit.estimation.measurements.AngularAzEl) SaastamoinenModel(org.orekit.models.earth.SaastamoinenModel)

Example 7 with RangeTroposphericDelayModifier

use of org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier in project Orekit by CS-SI.

the class TropoModifierTest method testRangeTropoModifier.

@Test
public void testRangeTropoModifier() throws OrekitException {
    Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    // create perfect range measurements
    for (final GroundStation station : context.stations) {
        station.getEastOffsetDriver().setSelected(true);
        station.getNorthOffsetDriver().setSelected(true);
        station.getZenithOffsetDriver().setSelected(true);
    }
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new RangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    propagator.setSlaveMode();
    final RangeTroposphericDelayModifier modifier = new RangeTroposphericDelayModifier(SaastamoinenModel.getStandardModel());
    for (final ObservedMeasurement<?> measurement : measurements) {
        final AbsoluteDate date = measurement.getDate();
        final SpacecraftState refState = propagator.propagate(date);
        Range range = (Range) measurement;
        EstimatedMeasurement<Range> evalNoMod = range.estimate(0, 0, new SpacecraftState[] { refState });
        // add modifier
        range.addModifier(modifier);
        EstimatedMeasurement<Range> eval = range.estimate(0, 0, new SpacecraftState[] { refState });
        final double diffMeters = eval.getEstimatedValue()[0] - evalNoMod.getEstimatedValue()[0];
        final double epsilon = 1e-6;
        Assert.assertTrue(Precision.compareTo(diffMeters, 12., epsilon) < 0);
        Assert.assertTrue(Precision.compareTo(diffMeters, 0., epsilon) > 0);
    }
}
Also used : Context(org.orekit.estimation.Context) GroundStation(org.orekit.estimation.measurements.GroundStation) TurnAroundRange(org.orekit.estimation.measurements.TurnAroundRange) Range(org.orekit.estimation.measurements.Range) RangeTroposphericDelayModifier(org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier) TurnAroundRangeTroposphericDelayModifier(org.orekit.estimation.measurements.modifiers.TurnAroundRangeTroposphericDelayModifier) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) Propagator(org.orekit.propagation.Propagator) RangeMeasurementCreator(org.orekit.estimation.measurements.RangeMeasurementCreator) TurnAroundRangeMeasurementCreator(org.orekit.estimation.measurements.TurnAroundRangeMeasurementCreator) ObservedMeasurement(org.orekit.estimation.measurements.ObservedMeasurement) Test(org.junit.Test)

Example 8 with RangeTroposphericDelayModifier

use of org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier in project Orekit by CS-SI.

the class RangeTest method genericTestStateDerivatives.

void genericTestStateDerivatives(final boolean isModifier, final boolean printResults, final double refErrorsPMedian, final double refErrorsPMean, final double refErrorsPMax, final double refErrorsVMedian, final double refErrorsVMean, final double refErrorsVMax) throws OrekitException {
    Context context = EstimationTestUtils.eccentricContext("regular-data:potential:tides");
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    // Create perfect range measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new RangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    // Lists for results' storage - Used only for derivatives with respect to state
    // "final" value to be seen by "handleStep" function of the propagator
    final List<Double> errorsP = new ArrayList<Double>();
    final List<Double> errorsV = new ArrayList<Double>();
    // Set master mode
    // Use a lambda function to implement "handleStep" function
    propagator.setMasterMode((OrekitStepInterpolator interpolator, boolean isLast) -> {
        for (final ObservedMeasurement<?> measurement : measurements) {
            // Play test if the measurement date is between interpolator previous and current date
            if ((measurement.getDate().durationFrom(interpolator.getPreviousState().getDate()) > 0.) && (measurement.getDate().durationFrom(interpolator.getCurrentState().getDate()) <= 0.)) {
                // Add modifiers if test implies it
                final RangeTroposphericDelayModifier modifier = new RangeTroposphericDelayModifier(SaastamoinenModel.getStandardModel());
                if (isModifier) {
                    ((Range) measurement).addModifier(modifier);
                }
                // We intentionally propagate to a date which is close to the
                // real spacecraft state but is *not* the accurate date, by
                // compensating only part of the downlink delay. This is done
                // in order to validate the partial derivatives with respect
                // to velocity. If we had chosen the proper state date, the
                // range would have depended only on the current position but
                // not on the current velocity.
                final double meanDelay = measurement.getObservedValue()[0] / Constants.SPEED_OF_LIGHT;
                final AbsoluteDate date = measurement.getDate().shiftedBy(-0.75 * meanDelay);
                final SpacecraftState state = interpolator.getInterpolatedState(date);
                final double[][] jacobian = measurement.estimate(0, 0, new SpacecraftState[] { state }).getStateDerivatives(0);
                // Jacobian reference value
                final double[][] jacobianRef;
                // Compute a reference value using finite differences
                jacobianRef = Differentiation.differentiate(new StateFunction() {

                    public double[] value(final SpacecraftState state) throws OrekitException {
                        return measurement.estimate(0, 0, new SpacecraftState[] { state }).getEstimatedValue();
                    }
                }, measurement.getDimension(), propagator.getAttitudeProvider(), OrbitType.CARTESIAN, PositionAngle.TRUE, 2.0, 3).value(state);
                Assert.assertEquals(jacobianRef.length, jacobian.length);
                Assert.assertEquals(jacobianRef[0].length, jacobian[0].length);
                // Errors & relative errors on the Jacobian
                double[][] dJacobian = new double[jacobian.length][jacobian[0].length];
                double[][] dJacobianRelative = new double[jacobian.length][jacobian[0].length];
                for (int i = 0; i < jacobian.length; ++i) {
                    for (int j = 0; j < jacobian[i].length; ++j) {
                        dJacobian[i][j] = jacobian[i][j] - jacobianRef[i][j];
                        dJacobianRelative[i][j] = FastMath.abs(dJacobian[i][j] / jacobianRef[i][j]);
                        if (j < 3) {
                            errorsP.add(dJacobianRelative[i][j]);
                        } else {
                            errorsV.add(dJacobianRelative[i][j]);
                        }
                    }
                }
                // Print values in console ?
                if (printResults) {
                    String stationName = ((Range) measurement).getStation().getBaseFrame().getName();
                    System.out.format(Locale.US, "%-15s  %-23s  %-23s  " + "%10.3e  %10.3e  %10.3e  " + "%10.3e  %10.3e  %10.3e  " + "%10.3e  %10.3e  %10.3e  " + "%10.3e  %10.3e  %10.3e%n", stationName, measurement.getDate(), date, dJacobian[0][0], dJacobian[0][1], dJacobian[0][2], dJacobian[0][3], dJacobian[0][4], dJacobian[0][5], dJacobianRelative[0][0], dJacobianRelative[0][1], dJacobianRelative[0][2], dJacobianRelative[0][3], dJacobianRelative[0][4], dJacobianRelative[0][5]);
                }
            }
        // End if measurement date between previous and current interpolator step
        }
    // End for loop on the measurements
    });
    // Print results on console ?
    if (printResults) {
        System.out.format(Locale.US, "%-15s  %-23s  %-23s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s%n", "Station", "Measurement Date", "State Date", "ΔdPx", "ΔdPy", "ΔdPz", "ΔdVx", "ΔdVy", "ΔdVz", "rel ΔdPx", "rel ΔdPy", "rel ΔdPz", "rel ΔdVx", "rel ΔdVy", "rel ΔdVz");
    }
    // Rewind the propagator to initial date
    propagator.propagate(context.initialOrbit.getDate());
    // Sort measurements chronologically
    measurements.sort(new ChronologicalComparator());
    // Propagate to final measurement's date
    propagator.propagate(measurements.get(measurements.size() - 1).getDate());
    // Convert lists to double[] and evaluate some statistics
    final double[] relErrorsP = errorsP.stream().mapToDouble(Double::doubleValue).toArray();
    final double[] relErrorsV = errorsV.stream().mapToDouble(Double::doubleValue).toArray();
    final double errorsPMedian = new Median().evaluate(relErrorsP);
    final double errorsPMean = new Mean().evaluate(relErrorsP);
    final double errorsPMax = new Max().evaluate(relErrorsP);
    final double errorsVMedian = new Median().evaluate(relErrorsV);
    final double errorsVMean = new Mean().evaluate(relErrorsV);
    final double errorsVMax = new Max().evaluate(relErrorsV);
    // Print the results on console ?
    if (printResults) {
        System.out.println();
        System.out.format(Locale.US, "Relative errors dR/dP -> Median: %6.3e / Mean: %6.3e / Max: %6.3e%n", errorsPMedian, errorsPMean, errorsPMax);
        System.out.format(Locale.US, "Relative errors dR/dV -> Median: %6.3e / Mean: %6.3e / Max: %6.3e%n", errorsVMedian, errorsVMean, errorsVMax);
    }
    Assert.assertEquals(0.0, errorsPMedian, refErrorsPMedian);
    Assert.assertEquals(0.0, errorsPMean, refErrorsPMean);
    Assert.assertEquals(0.0, errorsPMax, refErrorsPMax);
    Assert.assertEquals(0.0, errorsVMedian, refErrorsVMedian);
    Assert.assertEquals(0.0, errorsVMean, refErrorsVMean);
    Assert.assertEquals(0.0, errorsVMax, refErrorsVMax);
}
Also used : Mean(org.hipparchus.stat.descriptive.moment.Mean) Max(org.hipparchus.stat.descriptive.rank.Max) ArrayList(java.util.ArrayList) Median(org.hipparchus.stat.descriptive.rank.Median) AbsoluteDate(org.orekit.time.AbsoluteDate) SpacecraftState(org.orekit.propagation.SpacecraftState) Propagator(org.orekit.propagation.Propagator) OrekitException(org.orekit.errors.OrekitException) Context(org.orekit.estimation.Context) RangeTroposphericDelayModifier(org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier) OrekitStepInterpolator(org.orekit.propagation.sampling.OrekitStepInterpolator) NumericalPropagatorBuilder(org.orekit.propagation.conversion.NumericalPropagatorBuilder) StateFunction(org.orekit.utils.StateFunction) ChronologicalComparator(org.orekit.time.ChronologicalComparator)

Aggregations

RangeTroposphericDelayModifier (org.orekit.estimation.measurements.modifiers.RangeTroposphericDelayModifier)8 Context (org.orekit.estimation.Context)5 Propagator (org.orekit.propagation.Propagator)5 SpacecraftState (org.orekit.propagation.SpacecraftState)5 NumericalPropagatorBuilder (org.orekit.propagation.conversion.NumericalPropagatorBuilder)5 AbsoluteDate (org.orekit.time.AbsoluteDate)5 ArrayList (java.util.ArrayList)4 Mean (org.hipparchus.stat.descriptive.moment.Mean)4 Max (org.hipparchus.stat.descriptive.rank.Max)4 Median (org.hipparchus.stat.descriptive.rank.Median)4 GroundStation (org.orekit.estimation.measurements.GroundStation)4 Range (org.orekit.estimation.measurements.Range)4 OrekitStepInterpolator (org.orekit.propagation.sampling.OrekitStepInterpolator)4 ChronologicalComparator (org.orekit.time.ChronologicalComparator)4 HashMap (java.util.HashMap)3 GeodeticPoint (org.orekit.bodies.GeodeticPoint)3 OrekitException (org.orekit.errors.OrekitException)3 AngularAzEl (org.orekit.estimation.measurements.AngularAzEl)3 RangeRate (org.orekit.estimation.measurements.RangeRate)3 AngularRadioRefractionModifier (org.orekit.estimation.measurements.modifiers.AngularRadioRefractionModifier)3