1 /* Copyright 2002-2025 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.semianalytical.dsst.utilities;
18
19 import org.hipparchus.geometry.euclidean.threed.Vector3D;
20 import org.hipparchus.util.FastMath;
21 import org.hipparchus.util.MathUtils;
22 import org.orekit.frames.Frame;
23 import org.orekit.orbits.Orbit;
24 import org.orekit.time.AbsoluteDate;
25
26
27 /** Container class for common parameters used by all DSST forces.
28 * <p>
29 * Most of them are defined in Danielson paper at § 2.1.
30 * </p>
31 * @author Pascal Parraud
32 */
33 public class AuxiliaryElements {
34
35 /** Orbit date. */
36 private final AbsoluteDate date;
37
38 /** Orbit frame. */
39 private final Frame frame;
40
41 /** Eccentricity. */
42 private final double ecc;
43
44 /** Keplerian mean motion. */
45 private final double n;
46
47 /** Keplerian period. */
48 private final double period;
49
50 /** Semi-major axis. */
51 private final double sma;
52
53 /** x component of eccentricity vector. */
54 private final double k;
55
56 /** y component of eccentricity vector. */
57 private final double h;
58
59 /** x component of inclination vector. */
60 private final double q;
61
62 /** y component of inclination vector. */
63 private final double p;
64
65 /** Mean longitude. */
66 private final double lm;
67
68 /** True longitude. */
69 private final double lv;
70
71 /** Eccentric longitude. */
72 private final double le;
73
74 /** Retrograde factor I.
75 * <p>
76 * DSST model needs equinoctial orbit as internal representation.
77 * Classical equinoctial elements have discontinuities when inclination
78 * is close to zero. In this representation, I = +1. <br>
79 * To avoid this discontinuity, another representation exists and equinoctial
80 * elements can be expressed in a different way, called "retrograde" orbit.
81 * This implies I = -1. <br>
82 * As Orekit doesn't implement the retrograde orbit, I is always set to +1.
83 * But for the sake of consistency with the theory, the retrograde factor
84 * has been kept in the formulas.
85 * </p>
86 */
87 private final int I;
88
89 /** Orbit. */
90 private Orbit orbit;
91
92 /** B = sqrt(1 - h² - k²). */
93 private final double B;
94
95 /** C = 1 + p² + q². */
96 private final double C;
97
98 /** Equinoctial frame f vector. */
99 private final Vector3D f;
100
101 /** Equinoctial frame g vector. */
102 private final Vector3D g;
103
104 /** Equinoctial frame w vector. */
105 private final Vector3D w;
106
107 /** Simple constructor.
108 * @param orbit related mean orbit for auxiliary elements
109 * @param retrogradeFactor retrograde factor I [Eq. 2.1.2-(2)]
110 */
111 public AuxiliaryElements(final Orbit orbit, final int retrogradeFactor) {
112
113 // Orbit
114 this.orbit = orbit;
115
116 // Date of the orbit
117 date = orbit.getDate();
118
119 // Orbit definition frame
120 frame = orbit.getFrame();
121
122 // Eccentricity
123 ecc = orbit.getE();
124
125 // Keplerian mean motion
126 n = orbit.getKeplerianMeanMotion();
127
128 // Keplerian period
129 period = orbit.getKeplerianPeriod();
130
131 // Equinoctial elements [Eq. 2.1.2-(1)]
132 sma = orbit.getA();
133 k = orbit.getEquinoctialEx();
134 h = orbit.getEquinoctialEy();
135 q = orbit.getHx();
136 p = orbit.getHy();
137 lm = MathUtils.normalizeAngle(orbit.getLM(), FastMath.PI);
138 lv = MathUtils.normalizeAngle(orbit.getLv(), FastMath.PI);
139 le = MathUtils.normalizeAngle(orbit.getLE(), FastMath.PI);
140
141 // Retrograde factor [Eq. 2.1.2-(2)]
142 I = retrogradeFactor;
143
144 final double k2 = k * k;
145 final double h2 = h * h;
146 final double q2 = q * q;
147 final double p2 = p * p;
148
149 // A, B, C parameters [Eq. 2.1.6-(1)]
150 B = FastMath.sqrt(1 - k2 - h2);
151 C = 1 + q2 + p2;
152
153 // Equinoctial reference frame [Eq. 2.1.4-(1)]
154 final double ooC = 1. / C;
155 final double px2 = 2. * p;
156 final double qx2 = 2. * q;
157 final double pq2 = px2 * q;
158 f = new Vector3D(ooC, new Vector3D(1. - p2 + q2, pq2, -px2 * I));
159 g = new Vector3D(ooC, new Vector3D(pq2 * I, (1. + p2 - q2) * I, qx2));
160 w = new Vector3D(ooC, new Vector3D(px2, -qx2, (1. - p2 - q2) * I));
161 }
162
163 /** Get the orbit.
164 * @return the orbit
165 */
166 public Orbit getOrbit() {
167 return orbit;
168 }
169
170 /** Get the date of the orbit.
171 * @return the date
172 */
173 public AbsoluteDate getDate() {
174 return date;
175 }
176
177 /** Get the definition frame of the orbit.
178 * @return the definition frame
179 */
180 public Frame getFrame() {
181 return frame;
182 }
183
184 /** Get the eccentricity.
185 * @return ecc
186 */
187 public double getEcc() {
188 return ecc;
189 }
190
191 /** Get the Keplerian mean motion.
192 * @return n
193 */
194 public double getMeanMotion() {
195 return n;
196 }
197
198 /** Get the Keplerian period.
199 * @return period
200 */
201 public double getKeplerianPeriod() {
202 return period;
203 }
204
205 /** Get the semi-major axis.
206 * @return the semi-major axis a
207 */
208 public double getSma() {
209 return sma;
210 }
211
212 /** Get the x component of eccentricity vector.
213 * <p>
214 * This element called k in DSST corresponds to ex
215 * for the {@link org.orekit.orbits.EquinoctialOrbit}
216 * </p>
217 * @return k
218 */
219 public double getK() {
220 return k;
221 }
222
223 /** Get the y component of eccentricity vector.
224 * <p>
225 * This element called h in DSST corresponds to ey
226 * for the {@link org.orekit.orbits.EquinoctialOrbit}
227 * </p>
228 * @return h
229 */
230 public double getH() {
231 return h;
232 }
233
234 /** Get the x component of inclination vector.
235 * <p>
236 * This element called q in DSST corresponds to hx
237 * for the {@link org.orekit.orbits.EquinoctialOrbit}
238 * </p>
239 * @return q
240 */
241 public double getQ() {
242 return q;
243 }
244
245 /** Get the y component of inclination vector.
246 * <p>
247 * This element called p in DSST corresponds to hy
248 * for the {@link org.orekit.orbits.EquinoctialOrbit}
249 * </p>
250 * @return p
251 */
252 public double getP() {
253 return p;
254 }
255
256 /** Get the mean longitude.
257 * @return lm
258 */
259 public double getLM() {
260 return lm;
261 }
262
263 /** Get the true longitude.
264 * @return lv
265 */
266 public double getLv() {
267 return lv;
268 }
269
270 /** Get the eccentric longitude.
271 * @return lf
272 */
273 public double getLf() {
274 return le;
275 }
276
277 /** Get the retrograde factor.
278 * @return the retrograde factor I
279 */
280 public int getRetrogradeFactor() {
281 return I;
282 }
283
284 /** Get B = sqrt(1 - e²).
285 * @return B
286 */
287 public double getB() {
288 return B;
289 }
290
291 /** Get C = 1 + p² + q².
292 * @return C
293 */
294 public double getC() {
295 return C;
296 }
297
298 /** Get equinoctial frame vector f.
299 * @return f vector
300 */
301 public Vector3D getVectorF() {
302 return f;
303 }
304
305 /** Get equinoctial frame vector g.
306 * @return g vector
307 */
308 public Vector3D getVectorG() {
309 return g;
310 }
311
312 /** Get equinoctial frame vector w.
313 * @return w vector
314 */
315 public Vector3D getVectorW() {
316 return w;
317 }
318
319 }