1   /* Copyright 2002-2026 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.estimation.measurements;
18  
19  import java.util.Map;
20  
21  import org.hipparchus.analysis.differentiation.Gradient;
22  import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
23  import org.hipparchus.geometry.euclidean.threed.Vector3D;
24  import org.orekit.estimation.measurements.model.RaDecModel;
25  import org.orekit.frames.Frame;
26  import org.orekit.propagation.SpacecraftState;
27  import org.orekit.signal.FieldSignalReceptionCondition;
28  import org.orekit.signal.SignalReceptionCondition;
29  import org.orekit.signal.SignalTravelTimeModel;
30  import org.orekit.time.AbsoluteDate;
31  import org.orekit.time.FieldAbsoluteDate;
32  import org.orekit.utils.FieldPVCoordinatesProvider;
33  import org.orekit.utils.PVCoordinatesProvider;
34  import org.orekit.utils.TimeStampedFieldPVCoordinates;
35  import org.orekit.utils.TimeStampedPVCoordinates;
36  
37  /** Class modeling a Right Ascension - Declination measurement from an optical sensor, typically via astrometry.
38   * The angles are given using the axes of an inertial reference frame.
39   * The date of the measurement corresponds to the reception of the reflected signal.
40   *
41   * @author Thierry Ceolin
42   * @author Maxime Journot
43   * @since 9.0
44   */
45  public class AngularRaDec extends AngularMeasurement<AngularRaDec> {
46  
47      /** Type of the measurement. */
48      public static final String MEASUREMENT_TYPE = "AngularRaDec";
49  
50      /** Reference frame in which the right ascension - declination angles are given. */
51      private final Frame referenceFrame;
52  
53      /** Observer that receives signal from satellite. */
54      private final Observer observer;
55  
56      /** Perfect measurement model. */
57      private final RaDecModel measurementModel;
58  
59      /** Simple constructor using default light time delay.
60       * @param observer sensor from which measurement is performed
61       * @param referenceFrame Reference frame in which the right ascension - declination angles are given
62       * @param date date of the measurement
63       * @param angular observed value
64       * @param sigma theoretical standard deviation
65       * @param baseWeight base weight
66       * @param satellite satellite related to this measurement
67       * @since 9.3
68       */
69      public AngularRaDec(final Observer observer, final Frame referenceFrame, final AbsoluteDate date,
70                          final double[] angular, final double[] sigma, final double[] baseWeight,
71                          final ObservableSatellite satellite) {
72          this(observer, referenceFrame, date, angular, new MeasurementQuality(sigma, baseWeight), new SignalTravelTimeModel(), satellite);
73      }
74  
75      /** Constructor.
76       * @param observer sensor from which measurement is performed
77       * @param referenceFrame Reference frame in which the right ascension - declination angles are given
78       * @param date date of the measurement
79       * @param angular observed value
80       * @param measurementQuality measurement quality as used in estimation (in Orekit, the crossed-terms
81       *                           of the covariance matrix are only used by Kalman filters, not least squares)
82       * @param signalTravelTimeModel signal travel time model
83       * @param satellite satellite related to this measurement
84       * @since 14.0
85       */
86      public AngularRaDec(final Observer observer, final Frame referenceFrame, final AbsoluteDate date,
87                          final double[] angular, final MeasurementQuality measurementQuality,
88                          final SignalTravelTimeModel signalTravelTimeModel, final ObservableSatellite satellite) {
89          super(date, angular, measurementQuality, signalTravelTimeModel, satellite);
90          this.referenceFrame = referenceFrame;
91          this.measurementModel = new RaDecModel(referenceFrame, getSignalTravelTimeModel().getWarmedUpModel());
92          this.observer = observer;
93          observer.getParametersDrivers().forEach(this::addParameterDriver);
94      }
95  
96      /**
97       * Getter for the observer.
98       * @return observer
99       * @since 14.0
100      */
101     public Observer getObserver() {
102         return observer;
103     }
104 
105     /**
106      * Getter for the ground station.
107      * @return station
108      * @deprecated as of 14.0
109      */
110     @Deprecated
111     public GroundStation getStation() {
112         if (observer instanceof GroundStation station) {
113             return station;
114         } else {
115             return null;
116         }
117     }
118 
119     /** Get the reference frame in which the right ascension - declination angles are given.
120      * @return reference frame in which the right ascension - declination angles are given
121      */
122     public Frame getReferenceFrame() {
123         return referenceFrame;
124     }
125 
126     /** {@inheritDoc} */
127     @Override
128     protected EstimatedMeasurementBase<AngularRaDec> theoreticalEvaluationWithoutDerivatives(final int iteration,
129                                                                                              final int evaluation,
130                                                                                              final SpacecraftState[] states) {
131         // Compute emission date
132         final AbsoluteDate receptionDate = observer.getCorrectedReceptionDate(getDate());
133         final PVCoordinatesProvider receiver = observer.getPVCoordinatesProvider();
134         final SpacecraftState state = states[0];
135         final PVCoordinatesProvider emitter = AbstractParticipant.extractPVCoordinatesProvider(state, state.getPVCoordinates());
136         final Frame frame = state.getFrame();
137         final TimeStampedPVCoordinates receiverPV = receiver.getPVCoordinates(receptionDate, frame);
138         final SignalReceptionCondition receptionCondition = new SignalReceptionCondition(receptionDate,
139                 receiverPV.getPosition(), frame);
140         final AbsoluteDate emissionDate = computeEmissionDate(receptionCondition, emitter);
141 
142         // Evaluate angular measurement model (use state frame to avoid rounding error in case reference one is not Earth-centered)
143         final double[] raDec = measurementModel.value(receptionCondition, emitter, emissionDate);
144 
145         // Prepare the estimation
146         final double shift = emissionDate.durationFrom(state);
147         final SpacecraftState shiftedState = state.shiftedBy(shift);
148         final EstimatedMeasurementBase<AngularRaDec> estimated = new EstimatedMeasurementBase<>(this, iteration, evaluation,
149                 new SpacecraftState[] { shiftedState },
150                 new TimeStampedPVCoordinates[] { shiftedState.getPVCoordinates(), receiverPV });
151         estimated.setEstimatedValue(wrapFirstAngle(raDec[0]), raDec[1]);
152         return estimated;
153     }
154 
155     /** {@inheritDoc} */
156     @Override
157     protected EstimatedMeasurement<AngularRaDec> theoreticalEvaluation(final int iteration, final int evaluation,
158                                                                        final SpacecraftState[] states) {
159         // Right Ascension/declination (in reference frame) derivatives are computed with respect to spacecraft state in inertial frame
160         // and station parameters
161         // ----------------------
162         //
163         // Parameters:
164         //  - 0..2 - Position of the spacecraft in inertial frame
165         //  - 3..5 - Velocity of the spacecraft in inertial frame
166         //  - 6..n - station parameters (clock offset, station offsets, pole, prime meridian...)
167 
168         // Create the parameter indices map
169         final Map<String, Integer> paramIndices = getParameterIndices(states);
170         final int                  nbParams     = 6 * states.length + paramIndices.size();
171         final SpacecraftState state = states[0];
172         final TimeStampedFieldPVCoordinates<Gradient> pva = AbstractMeasurement.getCoordinates(state, 0, nbParams);
173 
174         // Compute emission date
175         final FieldAbsoluteDate<Gradient> receptionDate = observer.getCorrectedReceptionDateField(getDate(), nbParams, paramIndices);
176         final FieldPVCoordinatesProvider<Gradient> receiver = observer.getFieldPVCoordinatesProvider(nbParams, paramIndices);
177         final FieldPVCoordinatesProvider<Gradient> emitter = AbstractParticipant.extractFieldPVCoordinatesProvider(state, pva);
178         final Frame                frame        = states[0].getFrame();
179         final FieldVector3D<Gradient> receiverPosition = receiver.getPosition(receptionDate, frame);
180         final FieldAbsoluteDate<Gradient> emissionDate = computeEmissionDateField(referenceFrame, receiverPosition, receptionDate, emitter);
181 
182         // Evaluate angular measurement model (use state frame to avoid rounding error in case reference one is not Earth-centered)
183         final FieldSignalReceptionCondition<Gradient> receptionCondition = new FieldSignalReceptionCondition<>(receptionDate,
184                 receiverPosition, frame);
185         final Gradient[] raDec = measurementModel.value(receptionCondition, emitter, emissionDate);
186 
187         // Prepare the estimation
188         final double shift = emissionDate.toAbsoluteDate().durationFrom(state);
189         final SpacecraftState shiftedState = state.shiftedBy(shift);
190         final EstimatedMeasurement<AngularRaDec> estimated = new EstimatedMeasurement<>(this, iteration, evaluation,
191                 new SpacecraftState[] { shiftedState },
192                 new TimeStampedPVCoordinates[] { shiftedState.getPVCoordinates(),
193                 getObserver().getPVCoordinatesProvider().getPVCoordinates(receptionDate.toAbsoluteDate(), frame)});
194         fillEstimatedMeasurement(raDec[0], raDec[1], paramIndices, estimated);
195         return estimated;
196     }
197 
198     /** Calculate the Line Of Sight of the given measurement.
199      * @param outputFrame output frame of the line of sight vector
200      * @return Vector3D the line of Sight of the measurement
201      * @since 12.0
202      */
203     public Vector3D getObservedLineOfSight(final Frame outputFrame) {
204         return referenceFrame.getStaticTransformTo(outputFrame, getDate())
205             .transformVector(new Vector3D(getObservedValue()[0], getObservedValue()[1]));
206     }
207 }