BrouwerLyddanePropagatorBuilder.java

  1. /* Copyright 2002-2022 CS GROUP
  2.  * Licensed to CS GROUP (CS) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * CS licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *   http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.orekit.propagation.conversion;


  20. import java.util.List;

  22. import org.hipparchus.util.FastMath;
  23. import org.orekit.attitudes.AttitudeProvider;
  24. import org.orekit.attitudes.InertialProvider;
  25. import org.orekit.estimation.leastsquares.AbstractBatchLSModel;
  26. import org.orekit.estimation.leastsquares.BatchLSModel;
  27. import org.orekit.estimation.leastsquares.ModelObserver;
  28. import org.orekit.estimation.measurements.ObservedMeasurement;
  29. import org.orekit.estimation.sequential.AbstractKalmanModel;
  30. import org.orekit.estimation.sequential.CovarianceMatrixProvider;
  31. import org.orekit.estimation.sequential.KalmanModel;
  32. import org.orekit.forces.gravity.potential.GravityFieldFactory;
  33. import org.orekit.forces.gravity.potential.TideSystem;
  34. import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
  35. import org.orekit.orbits.Orbit;
  36. import org.orekit.orbits.OrbitType;
  37. import org.orekit.orbits.PositionAngle;
  38. import org.orekit.propagation.analytical.BrouwerLyddanePropagator;
  39. import org.orekit.propagation.analytical.tle.TLE;
  40. import org.orekit.utils.ParameterDriver;
  41. import org.orekit.utils.ParameterDriversList;

  43. /** Builder for Brouwer-Lyddane propagator.
  44.  * <p>
  45.  * By default, Brouwer-Lyddane model considers only the perturbations due to zonal harmonics.
  46.  * However, for low Earth orbits, the magnitude of the perturbative acceleration due to
  47.  * atmospheric drag can be significant. Warren Phipps' 1992 thesis considered the atmospheric
  48.  * drag by time derivatives of the <i>mean</i> mean anomaly using the catch-all coefficient M2.
  49.  *
  50.  * Usually, M2 is adjusted during an orbit determination process and it represents the
  51.  * combination of all unmodeled secular along-track effects (i.e. not just the atmospheric drag).
  52.  * The behavior of M2 is closed to the {@link TLE#getBStar()} parameter for the TLE.
  53.  *
  54.  * If the value of M2 is equal to {@link BrouwerLyddanePropagator#M2 0.0}, the along-track
  55.  * secular effects are not considered in the dynamical model. Typical values for M2 are not known.
  56.  * It depends on the orbit type. However, the value of M2 must be very small (e.g. between 1.0e-14 and 1.0e-15).
  57.  * The unit of M2 is rad/s².
  58.  * <p>
  59.  * To estimate the M2 parameter, it is necessary to call the {@link #getPropagationParametersDrivers()} method
  60.  * as follow:
  61.  * <pre>
  62.  *  for (ParameterDriver driver : builder.getPropagationParametersDrivers().getDrivers()) {
  63.  *     if (BrouwerLyddanePropagator.M2_NAME.equals(driver.getName())) {
  64.  *        driver.setSelected(true);
  65.  *     }
  66.  *  }
  67.  * </pre>
  68.  * @author Melina Vanel
  69.  * @author Bryan Cazabonne
  70.  * @since 11.1
  71.  */
  72. public class BrouwerLyddanePropagatorBuilder extends AbstractPropagatorBuilder implements OrbitDeterminationPropagatorBuilder {

  74.     /** Parameters scaling factor.
  75.      * <p>
  76.      * We use a power of 2 to avoid numeric noise introduction
  77.      * in the multiplications/divisions sequences.
  78.      * </p>
  79.      */
  80.     private static final double SCALE = FastMath.scalb(1.0, -32);

  82.     /** Provider for un-normalized coefficients. */
  83.     private final UnnormalizedSphericalHarmonicsProvider provider;

  85.     /** Build a new instance.
  86.      * <p>
  87.      * The template orbit is used as a model to {@link
  88.      * #createInitialOrbit() create initial orbit}. It defines the
  89.      * inertial frame, the central attraction coefficient, the orbit type, and is also
  90.      * used together with the {@code positionScale} to convert from the {@link
  91.      * org.orekit.utils.ParameterDriver#setNormalizedValue(double) normalized} parameters used by the
  92.      * callers of this builder to the real orbital parameters.
  93.      * </p>
  94.      *
  95.      * @param templateOrbit reference orbit from which real orbits will be built
  96.      * (note that the mu from this orbit will be overridden with the mu from the
  97.      * {@code provider})
  98.      * @param provider for un-normalized zonal coefficients
  99.      * @param positionAngle position angle type to use
  100.      * @param positionScale scaling factor used for orbital parameters normalization
  101.      * (typically set to the expected standard deviation of the position)
  102.      * @param M2 value of empirical drag coefficient in rad/s².
  103.      *        If equal to {@link BrouwerLyddanePropagator#M2} drag is not computed
  104.      * @see #BrouwerLyddanePropagatorBuilder(Orbit,
  105.      * UnnormalizedSphericalHarmonicsProvider, PositionAngle, double, AttitudeProvider, double)
  106.      */
  107.     public BrouwerLyddanePropagatorBuilder(final Orbit templateOrbit,
  108.                                            final UnnormalizedSphericalHarmonicsProvider provider,
  109.                                            final PositionAngle positionAngle,
  110.                                            final double positionScale,
  111.                                            final double M2) {
  112.         this(templateOrbit, provider, positionAngle, positionScale, InertialProvider.of(templateOrbit.getFrame()), M2);
  113.     }

  115.     /** Build a new instance.
  116.      * <p>
  117.      * The template orbit is used as a model to {@link
  118.      * #createInitialOrbit() create initial orbit}. It defines the
  119.      * inertial frame, the central attraction coefficient, the orbit type, and is also
  120.      * used together with the {@code positionScale} to convert from the {@link
  121.      * org.orekit.utils.ParameterDriver#setNormalizedValue(double) normalized} parameters used by the
  122.      * callers of this builder to the real orbital parameters.
  123.      * </p>
  124.      *
  125.      * @param templateOrbit reference orbit from which real orbits will be built
  126.      * (note that the mu from this orbit will be overridden with the mu from the
  127.      * {@code provider})
  128.      * @param referenceRadius reference radius of the Earth for the potential model (m)
  129.      * @param mu central attraction coefficient (m³/s²)
  130.      * @param tideSystem tide system
  131.      * @param c20 un-normalized zonal coefficient (about -1.08e-3 for Earth)
  132.      * @param c30 un-normalized zonal coefficient (about +2.53e-6 for Earth)
  133.      * @param c40 un-normalized zonal coefficient (about +1.62e-6 for Earth)
  134.      * @param c50 un-normalized zonal coefficient (about +2.28e-7 for Earth)
  135.      * @param orbitType orbit type to use
  136.      * @param positionAngle position angle type to use
  137.      * @param positionScale scaling factor used for orbital parameters normalization
  138.      * (typically set to the expected standard deviation of the position)
  139.      * @param M2 value of empirical drag coefficient in rad/s².
  140.      *        If equal to {@link BrouwerLyddanePropagator#M2} drag is not computed
  141.      * @see #BrouwerLyddanePropagatorBuilder(Orbit,
  142.      * UnnormalizedSphericalHarmonicsProvider, PositionAngle, double, AttitudeProvider, double)
  143.      */
  144.     public BrouwerLyddanePropagatorBuilder(final Orbit templateOrbit,
  145.                                            final double referenceRadius,
  146.                                            final double mu,
  147.                                            final TideSystem tideSystem,
  148.                                            final double c20,
  149.                                            final double c30,
  150.                                            final double c40,
  151.                                            final double c50,
  152.                                            final OrbitType orbitType,
  153.                                            final PositionAngle positionAngle,
  154.                                            final double positionScale,
  155.                                            final double M2) {
  156.         this(templateOrbit,
  157.              GravityFieldFactory.getUnnormalizedProvider(referenceRadius, mu, tideSystem,
  158.                                                          new double[][] {
  159.                                                              {
  160.                                                                  0
  161.                                                              }, {
  162.                                                                  0
  163.                                                              }, {
  164.                                                                  c20
  165.                                                              }, {
  166.                                                                  c30
  167.                                                              }, {
  168.                                                                  c40
  169.                                                              }, {
  170.                                                                  c50
  171.                                                              }
  172.                                                          }, new double[][] {
  173.                                                              {
  174.                                                                  0
  175.                                                              }, {
  176.                                                                  0
  177.                                                              }, {
  178.                                                                  0
  179.                                                              }, {
  180.                                                                  0
  181.                                                              }, {
  182.                                                                  0
  183.                                                              }, {
  184.                                                                  0
  185.                                                              }
  186.                                                          }),
  187.              positionAngle, positionScale, M2);
  188.     }

  190.     /** Build a new instance.
  191.      * <p>
  192.      * The template orbit is used as a model to {@link
  193.      * #createInitialOrbit() create initial orbit}. It defines the
  194.      * inertial frame, the central attraction coefficient, the orbit type, and is also
  195.      * used together with the {@code positionScale} to convert from the {@link
  196.      * org.orekit.utils.ParameterDriver#setNormalizedValue(double) normalized} parameters used by the
  197.      * callers of this builder to the real orbital parameters.
  198.      * </p>
  199.      * @param templateOrbit reference orbit from which real orbits will be built
  200.      * (note that the mu from this orbit will be overridden with the mu from the
  201.      * {@code provider})
  202.      * @param provider for un-normalized zonal coefficients
  203.      * @param positionAngle position angle type to use
  204.      * @param positionScale scaling factor used for orbital parameters normalization
  205.      * (typically set to the expected standard deviation of the position)
  206.      * @param M2 value of empirical drag coefficient in rad/s².
  207.      *        If equal to {@link BrouwerLyddanePropagator#M2} drag is not computed
  208.      * @param attitudeProvider attitude law to use
  209.      */
  210.     public BrouwerLyddanePropagatorBuilder(final Orbit templateOrbit,
  211.                                            final UnnormalizedSphericalHarmonicsProvider provider,
  212.                                            final PositionAngle positionAngle,
  213.                                            final double positionScale,
  214.                                            final AttitudeProvider attitudeProvider,
  215.                                            final double M2) {
  216.         super(overrideMu(templateOrbit, provider, positionAngle), positionAngle, positionScale, true, attitudeProvider);
  217.         this.provider = provider;
  218.         // initialize M2 driver
  219.         final ParameterDriver M2Driver = new ParameterDriver(BrouwerLyddanePropagator.M2_NAME, M2, SCALE,
  220.                                                              Double.NEGATIVE_INFINITY,
  221.                                                              Double.POSITIVE_INFINITY);
  222.         addSupportedParameter(M2Driver);
  223.     }

  225.     /** Override central attraction coefficient.
  226.      * @param templateOrbit template orbit
  227.      * @param provider gravity field provider
  228.      * @param positionAngle position angle type to use
  229.      * @return orbit with overridden central attraction coefficient
  230.      */
  231.     private static Orbit overrideMu(final Orbit templateOrbit,
  232.                                     final UnnormalizedSphericalHarmonicsProvider provider,
  233.                                     final PositionAngle positionAngle) {
  234.         final double[] parameters    = new double[6];
  235.         final double[] parametersDot = templateOrbit.hasDerivatives() ? new double[6] : null;
  236.         templateOrbit.getType().mapOrbitToArray(templateOrbit, positionAngle, parameters, parametersDot);
  237.         return templateOrbit.getType().mapArrayToOrbit(parameters, parametersDot, positionAngle,
  238.                                                        templateOrbit.getDate(),
  239.                                                        provider.getMu(),
  240.                                                        templateOrbit.getFrame());
  241.     }

  243.     /** {@inheritDoc} */
  244.     public BrouwerLyddanePropagator buildPropagator(final double[] normalizedParameters) {
  245.         setParameters(normalizedParameters);

  247.         // Update M2 value and selection
  248.         double  newM2      = 0.0;
  249.         boolean isSelected = false;
  250.         for (final ParameterDriver driver : getPropagationParametersDrivers().getDrivers()) {
  251.             if (BrouwerLyddanePropagator.M2_NAME.equals(driver.getName())) {
  252.                 newM2      = driver.getValue();
  253.                 isSelected = driver.isSelected();
  254.             }
  255.         }

  257.         // Initialize propagator
  258.         final BrouwerLyddanePropagator propagator = new BrouwerLyddanePropagator(createInitialOrbit(), getAttitudeProvider(), provider, newM2);
  259.         propagator.getParametersDrivers().get(0).setSelected(isSelected);

  261.         // Return
  262.         return propagator;

  264.     }

  266.     /** {@inheritDoc} */
  267.     @Override
  268.     public AbstractBatchLSModel buildLSModel(final OrbitDeterminationPropagatorBuilder[] builders,
  269.                                              final List<ObservedMeasurement<?>> measurements,
  270.                                              final ParameterDriversList estimatedMeasurementsParameters,
  271.                                              final ModelObserver observer) {
  272.         return new BatchLSModel(builders, measurements, estimatedMeasurementsParameters, observer);
  273.     }

  275.     /** {@inheritDoc} */
  276.     @Override
  277.     public AbstractKalmanModel buildKalmanModel(final List<OrbitDeterminationPropagatorBuilder> propagatorBuilders,
  278.                                                 final List<CovarianceMatrixProvider> covarianceMatricesProviders,
  279.                                                 final ParameterDriversList estimatedMeasurementsParameters,
  280.                                                 final CovarianceMatrixProvider measurementProcessNoiseMatrix) {
  281.         return new KalmanModel(propagatorBuilders, covarianceMatricesProviders, estimatedMeasurementsParameters, measurementProcessNoiseMatrix);
  282.     }

  284. }
Created with JaCoCo 0.8.7.202105040129