Examples with OrekitException org.orekit.errors.OrekitException used on opensource projects

Example 6 with OrekitException

use of org.orekit.errors.OrekitException in project Orekit by CS-SI.

the class HarrisPriesterTest method gradientComponent.

private double gradientComponent(final Atmosphere atm, final Vector3D position, final Vector3D direction, final double dP, final Frame frame, final ComponentGetter getter) {
    FiniteDifferencesDifferentiator differentiator = new FiniteDifferencesDifferentiator(5, dP);
    UnivariateFunction f = delta -> {
        try {
            return getter.get(atm.getVelocity(date, new Vector3D(1, position, delta, direction), frame));
        } catch (OrekitException oe) {
            return Double.NaN;
        }
    };
    return differentiator.differentiate(f).value(new DSFactory(1, 1).variable(0, 0.0)).getPartialDerivative(1);
}

Example 7 with OrekitException

use of org.orekit.errors.OrekitException in project Orekit by CS-SI.

the class TurnAroundRangeAnalyticTest method genericTestParameterDerivatives.

/**
 * Generic test function for derivatives with respect to parameters (station's position in station's topocentric frame)
 * @param isModifier Use of atmospheric modifiers
 * @param isFiniteDifferences Finite differences reference calculation if true, TurnAroundRange class otherwise
 * @param printResults Print the results ?
 * @throws OrekitException
 */
void genericTestParameterDerivatives(final boolean isModifier, final boolean isFiniteDifferences, final boolean printResults, final double refErrorQMMedian, final double refErrorQMMean, final double refErrorQMMax, final double refErrorQSMedian, final double refErrorQSMean, final double refErrorQSMax) 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 TAR measurements
    for (Map.Entry<GroundStation, GroundStation> entry : context.TARstations.entrySet()) {
        final GroundStation masterStation = entry.getKey();
        final GroundStation slaveStation = entry.getValue();
        masterStation.getEastOffsetDriver().setSelected(true);
        masterStation.getNorthOffsetDriver().setSelected(true);
        masterStation.getZenithOffsetDriver().setSelected(true);
        slaveStation.getEastOffsetDriver().setSelected(true);
        slaveStation.getNorthOffsetDriver().setSelected(true);
        slaveStation.getZenithOffsetDriver().setSelected(true);
    }
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new TurnAroundRangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    propagator.setSlaveMode();
    // Print results on console ? Header
    if (printResults) {
        System.out.format(Locale.US, "%-15s %-15s %-23s  %-23s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s%n", "Master Station", "Slave Station", "Measurement Date", "State Date", "ΔdQMx", "rel ΔdQMx", "ΔdQMy", "rel ΔdQMy", "ΔdQMz", "rel ΔdQMz", "ΔdQSx", "rel ΔdQSx", "ΔdQSy", "rel ΔdQSy", "ΔdQSz", "rel ΔdQSz");
    }
    // List to store the results for master and slave station
    final List<Double> relErrorQMList = new ArrayList<Double>();
    final List<Double> relErrorQSList = new ArrayList<Double>();
    // Loop on the measurements
    for (final ObservedMeasurement<?> measurement : measurements) {
        // Add modifiers if test implies it
        final TurnAroundRangeTroposphericDelayModifier modifier = new TurnAroundRangeTroposphericDelayModifier(SaastamoinenModel.getStandardModel());
        if (isModifier) {
            ((TurnAroundRange) measurement).addModifier(modifier);
        }
        // parameter corresponding to station position offset
        final GroundStation masterStationParameter = ((TurnAroundRange) measurement).getMasterStation();
        final GroundStation slaveStationParameter = ((TurnAroundRange) measurement).getSlaveStation();
        // 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 = propagator.propagate(date);
        final ParameterDriver[] drivers = new ParameterDriver[] { masterStationParameter.getEastOffsetDriver(), masterStationParameter.getNorthOffsetDriver(), masterStationParameter.getZenithOffsetDriver(), slaveStationParameter.getEastOffsetDriver(), slaveStationParameter.getNorthOffsetDriver(), slaveStationParameter.getZenithOffsetDriver() };
        // Print results on console ? Stations' names
        if (printResults) {
            String masterStationName = masterStationParameter.getBaseFrame().getName();
            String slaveStationName = slaveStationParameter.getBaseFrame().getName();
            System.out.format(Locale.US, "%-15s %-15s %-23s  %-23s  ", masterStationName, slaveStationName, measurement.getDate(), date);
        }
        // Loop on the parameters
        for (int i = 0; i < 6; ++i) {
            // Analytical computation of the parameters derivatives
            final EstimatedMeasurement<TurnAroundRange> TAR = new TurnAroundRangeAnalytic((TurnAroundRange) measurement).theoreticalEvaluationAnalytic(0, 0, propagator.getInitialState(), state);
            // Optional modifier addition
            if (isModifier) {
                modifier.modify(TAR);
            }
            final double[] gradient = TAR.getParameterDerivatives(drivers[i]);
            Assert.assertEquals(1, measurement.getDimension());
            Assert.assertEquals(1, gradient.length);
            // Reference value
            double ref;
            if (isFiniteDifferences) {
                // Compute a reference value using finite differences
                final ParameterFunction dMkdP = Differentiation.differentiate(new ParameterFunction() {

                    /**
                     * {@inheritDoc}
                     */
                    @Override
                    public double value(final ParameterDriver parameterDriver) throws OrekitException {
                        return measurement.estimate(0, 0, new SpacecraftState[] { state }).getEstimatedValue()[0];
                    }
                }, drivers[i], 3, 20.0);
                ref = dMkdP.value(drivers[i]);
            } else {
                // Compute a reference value using TurnAroundRange function
                ref = measurement.estimate(0, 0, new SpacecraftState[] { state }).getParameterDerivatives(drivers[i])[0];
            }
            // Deltas
            double dGradient = gradient[0] - ref;
            double dGradientRelative = FastMath.abs(dGradient / ref);
            // Print results on console ? Gradient difference
            if (printResults) {
                System.out.format(Locale.US, "%10.3e  %10.3e  ", dGradient, dGradientRelative);
            }
            // Add relative error to the list
            if (i < 3) {
                relErrorQMList.add(dGradientRelative);
            } else {
                relErrorQSList.add(dGradientRelative);
            }
        }
        // End for loop on the parameters
        if (printResults) {
            System.out.format(Locale.US, "%n");
        }
    }
    // End for loop on the measurements
    // Convert error list to double[]
    final double[] relErrorQM = relErrorQMList.stream().mapToDouble(Double::doubleValue).toArray();
    final double[] relErrorQS = relErrorQSList.stream().mapToDouble(Double::doubleValue).toArray();
    // Compute statistics
    final double relErrorsQMMedian = new Median().evaluate(relErrorQM);
    final double relErrorsQMMean = new Mean().evaluate(relErrorQM);
    final double relErrorsQMMax = new Max().evaluate(relErrorQM);
    final double relErrorsQSMedian = new Median().evaluate(relErrorQS);
    final double relErrorsQSMean = new Mean().evaluate(relErrorQS);
    final double relErrorsQSMax = new Max().evaluate(relErrorQS);
    // Print the results on console ?
    if (printResults) {
        System.out.println();
        System.out.format(Locale.US, "Relative errors dR/dQ master station -> Median: %6.3e / Mean: %6.3e / Max: %6.3e%n", relErrorsQMMedian, relErrorsQMMean, relErrorsQMMax);
        System.out.format(Locale.US, "Relative errors dR/dQ slave station  -> Median: %6.3e / Mean: %6.3e / Max: %6.3e%n", relErrorsQSMedian, relErrorsQSMean, relErrorsQSMax);
    }
    // Check values
    Assert.assertEquals(0.0, relErrorsQMMedian, refErrorQMMedian);
    Assert.assertEquals(0.0, relErrorsQMMean, refErrorQMMean);
    Assert.assertEquals(0.0, relErrorsQMMax, refErrorQMMax);
    Assert.assertEquals(0.0, relErrorsQSMedian, refErrorQSMedian);
    Assert.assertEquals(0.0, relErrorsQSMean, refErrorQSMean);
    Assert.assertEquals(0.0, relErrorsQSMax, refErrorQSMax);
}

Example 8 with OrekitException

use of org.orekit.errors.OrekitException in project Orekit by CS-SI.

the class TurnAroundRangeMeasurementCreator method handleStep.

/**
 * Function handling the steps of the propagator
 * A turn-around measurement needs 2 stations, a master and a slave
 * The measurement is a signal:
 * - Emitted from the master ground station
 * - Reflected on the spacecraft
 * - Reflected on the slave ground station
 * - Reflected on the spacecraft again
 * - Received on the master ground station
 * Its value is the elapsed time between emission and reception
 * divided by 2c were c is the speed of light.
 *
 * The path of the signal is divided into 2 legs:
 *  - The 1st leg goes from emission by the master station to reception by the slave station
 *  - The 2nd leg goes from emission by the slave station to reception by the master station
 *
 * The spacecraft state date should, after a few iterations of the estimation process, be
 * set to the date of arrival/departure of the signal to/from the slave station.
 * It is guaranteed by implementation of the estimated measurement.
 * This is done to avoid big shifts in time to compute the transit states.
 * See TurnAroundRange.java for more
 * Thus the spacecraft date is the date when the 1st leg of the path ends and the 2nd leg begins
 */
public void handleStep(final SpacecraftState currentState, final boolean isLast) throws OrekitException {
    try {
        for (Map.Entry<GroundStation, GroundStation> entry : context.TARstations.entrySet()) {
            final GroundStation masterStation = entry.getKey();
            final GroundStation slaveStation = entry.getValue();
            final AbsoluteDate date = currentState.getDate();
            final Frame inertial = currentState.getFrame();
            final Vector3D position = currentState.toTransform().getInverse().transformPosition(antennaPhaseCenter);
            // Create a TAR measurement only if elevation for both stations is higher than elevationMin°
            if ((masterStation.getBaseFrame().getElevation(position, inertial, date) > FastMath.toRadians(30.0)) && (slaveStation.getBaseFrame().getElevation(position, inertial, date) > FastMath.toRadians(30.0))) {
                // The solver used
                final UnivariateSolver solver = new BracketingNthOrderBrentSolver(1.0e-12, 5);
                // Spacecraft date t = date of arrival/departure of the signal to/from from the slave station
                // Slave station position in inertial frame at t
                final Vector3D slaveStationPosition = slaveStation.getOffsetToInertial(inertial, date).transformPosition(Vector3D.ZERO);
                // Downlink time of flight to slave station
                // The date of arrival/departure of the signal to/from the slave station is known and
                // equal to spacecraft date t.
                // Therefore we can use the function "downlinkTimeOfFlight" from GroundStation class
                // final double slaveTauD = slaveStation.downlinkTimeOfFlight(currentState, date);
                final double slaveTauD = solver.solve(1000, new UnivariateFunction() {

                    public double value(final double x) throws OrekitExceptionWrapper {
                        final SpacecraftState transitState = currentState.shiftedBy(-x);
                        final double d = Vector3D.distance(transitState.toTransform().getInverse().transformPosition(antennaPhaseCenter), slaveStationPosition);
                        return d - x * Constants.SPEED_OF_LIGHT;
                    }
                }, -1.0, 1.0);
                // Uplink time of flight from slave station
                // A solver is used to know where the satellite is when it receives the signal
                // back from the slave station
                final double slaveTauU = solver.solve(1000, new UnivariateFunction() {

                    public double value(final double x) throws OrekitExceptionWrapper {
                        final SpacecraftState transitState = currentState.shiftedBy(+x);
                        final double d = Vector3D.distance(transitState.toTransform().getInverse().transformPosition(antennaPhaseCenter), slaveStationPosition);
                        return d - x * Constants.SPEED_OF_LIGHT;
                    }
                }, -1.0, 1.0);
                // Find the position of the master station at signal departure and arrival
                // ----
                // Transit state position & date for the 1st leg of the signal path
                final SpacecraftState S1 = currentState.shiftedBy(-slaveTauD);
                final Vector3D P1 = S1.toTransform().getInverse().transformPosition(antennaPhaseCenter);
                final AbsoluteDate T1 = date.shiftedBy(-slaveTauD);
                // Transit state position & date for the 2nd leg of the signal path
                final Vector3D P2 = currentState.shiftedBy(+slaveTauU).toTransform().getInverse().transformPosition(antennaPhaseCenter);
                final AbsoluteDate T2 = date.shiftedBy(+slaveTauU);
                // Master station downlink delay - from P2 to master station
                // We use a solver to know where the master station is when it receives
                // the signal back from the satellite on the 2nd leg of the path
                final double masterTauD = solver.solve(1000, new UnivariateFunction() {

                    public double value(final double x) throws OrekitExceptionWrapper {
                        try {
                            final Transform t = masterStation.getOffsetToInertial(inertial, T2.shiftedBy(+x));
                            final double d = Vector3D.distance(P2, t.transformPosition(Vector3D.ZERO));
                            return d - x * Constants.SPEED_OF_LIGHT;
                        } catch (OrekitException oe) {
                            throw new OrekitExceptionWrapper(oe);
                        }
                    }
                }, -1.0, 1.0);
                final AbsoluteDate masterReceptionDate = T2.shiftedBy(+masterTauD);
                final TimeStampedPVCoordinates masterStationAtReception = masterStation.getOffsetToInertial(inertial, masterReceptionDate).transformPVCoordinates(new TimeStampedPVCoordinates(masterReceptionDate, PVCoordinates.ZERO));
                // Master station uplink delay - from master station to P1
                // Here the state date is known. Thus we can use the function "signalTimeOfFlight"
                // of the AbstractMeasurement class
                final double masterTauU = AbstractMeasurement.signalTimeOfFlight(masterStationAtReception, P1, T1);
                final AbsoluteDate masterEmissionDate = T1.shiftedBy(-masterTauU);
                final Vector3D masterStationAtEmission = masterStation.getOffsetToInertial(inertial, masterEmissionDate).transformPosition(Vector3D.ZERO);
                // Uplink/downlink distance from/to slave station
                final double slaveDownLinkDistance = Vector3D.distance(P1, slaveStationPosition);
                final double slaveUpLinkDistance = Vector3D.distance(P2, slaveStationPosition);
                // Uplink/downlink distance from/to master station
                final double masterUpLinkDistance = Vector3D.distance(P1, masterStationAtEmission);
                final double masterDownLinkDistance = Vector3D.distance(P2, masterStationAtReception.getPosition());
                addMeasurement(new TurnAroundRange(masterStation, slaveStation, masterReceptionDate, 0.5 * (masterUpLinkDistance + slaveDownLinkDistance + slaveUpLinkDistance + masterDownLinkDistance), 1.0, 10));
            }
        }
    } catch (OrekitExceptionWrapper oew) {
        throw new OrekitException(oew.getException());
    } catch (OrekitException oe) {
        throw new OrekitException(oe);
    }
}

Example 9 with OrekitException

use of org.orekit.errors.OrekitException in project Orekit by CS-SI.

the class TurnAroundRangeTest 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");
    // Context context = EstimationTestUtils.geoStationnaryContext();
    final NumericalPropagatorBuilder propagatorBuilder = context.createBuilder(OrbitType.KEPLERIAN, PositionAngle.TRUE, true, 1.0e-6, 60.0, 0.001);
    // create perfect range2 measurements
    final Propagator propagator = EstimationTestUtils.createPropagator(context.initialOrbit, propagatorBuilder);
    final List<ObservedMeasurement<?>> measurements = EstimationTestUtils.createMeasurements(propagator, new TurnAroundRangeMeasurementCreator(context), 1.0, 3.0, 300.0);
    propagator.setSlaveMode();
    double[] errorsP = new double[3 * measurements.size()];
    double[] errorsV = new double[3 * measurements.size()];
    int indexP = 0;
    int indexV = 0;
    // Print the results ? Header
    if (printResults) {
        System.out.format(Locale.US, "%-15s  %-15s  %-23s  %-23s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s  " + "%10s  %10s  %10s%n", "Master Station", "Slave Station", "Measurement Date", "State Date", "ΔdPx", "ΔdPy", "ΔdPz", "ΔdVx", "ΔdVy", "ΔdVz", "rel ΔdPx", "rel ΔdPy", "rel ΔdPz", "rel ΔdVx", "rel ΔdVy", "rel ΔdVz");
    }
    // Loop on the measurements
    for (final ObservedMeasurement<?> measurement : measurements) {
        // Add modifiers if test implies it
        final TurnAroundRangeTroposphericDelayModifier modifier = new TurnAroundRangeTroposphericDelayModifier(SaastamoinenModel.getStandardModel());
        if (isModifier) {
            ((TurnAroundRange) 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 = propagator.propagate(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);
        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[indexP++] = dJacobianRelative[i][j];
                } else {
                    errorsV[indexV++] = dJacobianRelative[i][j];
                }
            }
        }
        // Print results on the console ? Print the Jacobian
        if (printResults) {
            String masterStationName = ((TurnAroundRange) measurement).getMasterStation().getBaseFrame().getName();
            String slaveStationName = ((TurnAroundRange) measurement).getSlaveStation().getBaseFrame().getName();
            System.out.format(Locale.US, "%-15s  %-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", masterStationName, slaveStationName, 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 loop on the measurements
    // Compute some statistics
    final double errorsPMedian = new Median().evaluate(errorsP);
    final double errorsPMean = new Mean().evaluate(errorsP);
    final double errorsPMax = new Max().evaluate(errorsP);
    final double errorsVMedian = new Median().evaluate(errorsV);
    final double errorsVMean = new Mean().evaluate(errorsV);
    final double errorsVMax = new Max().evaluate(errorsV);
    // Print the results on console ? Final results
    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);
}

Example 10 with OrekitException

use of org.orekit.errors.OrekitException in project Orekit by CS-SI.

the class HolmesFeatherstoneAttractionModelTest method setUp.

@Before
public void setUp() {
    itrf = null;
    propagator = null;
    Utils.setDataRoot("regular-data");
    try {
        // Eigen 6s model truncated to degree 6
        mu = 3.986004415e+14;
        ae = 6378136.460;
        normalizedC20 = -4.84165299820e-04;
        normalizedC30 = 9.57211326674e-07;
        normalizedC40 = 5.39990167207e-07;
        normalizedC50 = 6.86846073356e-08;
        normalizedC60 = -1.49953256913e-07;
        unnormalizedC20 = FastMath.sqrt(5) * normalizedC20;
        unnormalizedC30 = FastMath.sqrt(7) * normalizedC30;
        unnormalizedC40 = FastMath.sqrt(9) * normalizedC40;
        unnormalizedC50 = FastMath.sqrt(11) * normalizedC50;
        unnormalizedC60 = FastMath.sqrt(13) * normalizedC60;
        itrf = FramesFactory.getITRF(IERSConventions.IERS_2010, true);
        double[] absTolerance = { 0.001, 1.0e-9, 1.0e-9, 1.0e-6, 1.0e-6, 1.0e-6, 0.001 };
        double[] relTolerance = { 1.0e-7, 1.0e-4, 1.0e-4, 1.0e-7, 1.0e-7, 1.0e-7, 1.0e-7 };
        AdaptiveStepsizeIntegrator integrator = new DormandPrince853Integrator(0.001, 1000, absTolerance, relTolerance);
        integrator.setInitialStepSize(60);
        propagator = new NumericalPropagator(integrator);
    } catch (OrekitException oe) {
        Assert.fail(oe.getMessage());
    }
}