1 /* Copyright 2002-2019 CS Systèmes d'Information
2 * Licensed to CS Systèmes d'Information (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.analytical.tle;
18
19 import org.hipparchus.util.FastMath;
20 import org.orekit.attitudes.AttitudeProvider;
21
22 /** This class contains methods to compute propagated coordinates with the SGP4 model.
23 * <p>
24 * The user should not bother in this class since it is handled internaly by the
25 * {@link TLEPropagator}.
26 * </p>
27 * <p>This implementation is largely inspired from the paper and source code <a
28 * href="http://www.celestrak.com/publications/AIAA/2006-6753/">Revisiting Spacetrack
29 * Report #3</a> and is fully compliant with its results and tests cases.</p>
30 * @author Felix R. Hoots, Ronald L. Roehrich, December 1980 (original fortran)
31 * @author David A. Vallado, Paul Crawford, Richard Hujsak, T.S. Kelso (C++ translation and improvements)
32 * @author Fabien Maussion (java translation)
33 */
34 public class SGP4 extends TLEPropagator {
35
36 /** If perige is less than 220 km, some calculus are avoided. */
37 private boolean lessThan220;
38
39 /** (1 + eta * cos(M0))³. */
40 private double delM0;
41
42 // CHECKSTYLE: stop JavadocVariable check
43 private double d2;
44 private double d3;
45 private double d4;
46 private double t3cof;
47 private double t4cof;
48 private double t5cof;
49 private double sinM0;
50 private double omgcof;
51 private double xmcof;
52 private double c5;
53 // CHECKSTYLE: resume JavadocVariable check
54
55 /** Constructor for a unique initial TLE.
56 * @param initialTLE the TLE to propagate.
57 * @param attitudeProvider provider for attitude computation
58 * @param mass spacecraft mass (kg)
59 */
60 public SGP4(final TLE initialTLE, final AttitudeProvider attitudeProvider,
61 final double mass) {
62 super(initialTLE, attitudeProvider, mass);
63 }
64
65 /** Initialization proper to each propagator (SGP or SDP).
66 */
67 protected void sxpInitialize() {
68
69 // For perigee less than 220 kilometers, the equations are truncated to
70 // linear variation in sqrt a and quadratic variation in mean anomaly.
71 // Also, the c3 term, the delta omega term, and the delta m term are dropped.
72 lessThan220 = perige < 220;
73 if (!lessThan220) {
74 final double c1sq = c1 * c1;
75 delM0 = 1.0 + eta * FastMath.cos(tle.getMeanAnomaly());
76 delM0 *= delM0 * delM0;
77 d2 = 4 * a0dp * tsi * c1sq;
78 final double temp = d2 * tsi * c1 / 3.0;
79 d3 = (17 * a0dp + s4) * temp;
80 d4 = 0.5 * temp * a0dp * tsi * (221 * a0dp + 31 * s4) * c1;
81 t3cof = d2 + 2 * c1sq;
82 t4cof = 0.25 * (3 * d3 + c1 * (12 * d2 + 10 * c1sq));
83 t5cof = 0.2 * (3 * d4 + 12 * c1 * d3 + 6 * d2 * d2 + 15 * c1sq * (2 * d2 + c1sq));
84 sinM0 = FastMath.sin(tle.getMeanAnomaly());
85 if (tle.getE() < 1e-4) {
86 omgcof = 0.;
87 xmcof = 0.;
88 } else {
89 final double c3 = coef * tsi * TLEConstants.A3OVK2 * xn0dp *
90 TLEConstants.NORMALIZED_EQUATORIAL_RADIUS * sini0 / tle.getE();
91 xmcof = -TLEConstants.TWO_THIRD * coef * tle.getBStar() *
92 TLEConstants.NORMALIZED_EQUATORIAL_RADIUS / eeta;
93 omgcof = tle.getBStar() * c3 * FastMath.cos(tle.getPerigeeArgument());
94 }
95 }
96
97 c5 = 2 * coef1 * a0dp * beta02 * (1 + 2.75 * (etasq + eeta) + eeta * etasq);
98 // initialized
99 }
100
101 /** Propagation proper to each propagator (SGP or SDP).
102 * @param tSince the offset from initial epoch (min)
103 */
104 protected void sxpPropagate(final double tSince) {
105
106 // Update for secular gravity and atmospheric drag.
107 final double xmdf = tle.getMeanAnomaly() + xmdot * tSince;
108 final double omgadf = tle.getPerigeeArgument() + omgdot * tSince;
109 final double xn0ddf = tle.getRaan() + xnodot * tSince;
110 omega = omgadf;
111 double xmp = xmdf;
112 final double tsq = tSince * tSince;
113 xnode = xn0ddf + xnodcf * tsq;
114 double tempa = 1 - c1 * tSince;
115 double tempe = tle.getBStar() * c4 * tSince;
116 double templ = t2cof * tsq;
117
118 if (!lessThan220) {
119 final double delomg = omgcof * tSince;
120 double delm = 1. + eta * FastMath.cos(xmdf);
121 delm = xmcof * (delm * delm * delm - delM0);
122 final double temp = delomg + delm;
123 xmp = xmdf + temp;
124 omega = omgadf - temp;
125 final double tcube = tsq * tSince;
126 final double tfour = tSince * tcube;
127 tempa = tempa - d2 * tsq - d3 * tcube - d4 * tfour;
128 tempe = tempe + tle.getBStar() * c5 * (FastMath.sin(xmp) - sinM0);
129 templ = templ + t3cof * tcube + tfour * (t4cof + tSince * t5cof);
130 }
131
132 a = a0dp * tempa * tempa;
133 e = tle.getE() - tempe;
134
135 // A highly arbitrary lower limit on e, of 1e-6:
136 if (e < 1e-6) {
137 e = 1e-6;
138 }
139
140 xl = xmp + omega + xnode + xn0dp * templ;
141
142 i = tle.getI();
143
144 }
145
146 }