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.data;
18
19 import org.hipparchus.RealFieldElement;
20
21 /** Class for general terms.
22 * @param <T> the type of the field elements
23 * @author Luc Maisonobe
24 */
25 class GeneralTerm extends SeriesTerm {
26
27 /** Coefficient for mean anomaly of the Moon. */
28 private final int cL;
29
30 /** Coefficient for mean anomaly of the Sun. */
31 private final int cLPrime;
32
33 /** Coefficient for L - Ω where L is the mean longitude of the Moon. */
34 private final int cF;
35
36 /** Coefficient for mean elongation of the Moon from the Sun. */
37 private final int cD;
38
39 /** Coefficient for mean longitude of the ascending node of the Moon. */
40 private final int cOmega;
41
42 /** Coefficient for mean Mercury longitude. */
43 private final int cMe;
44
45 /** Coefficient for mean Venus longitude. */
46 private final int cVe;
47
48 /** Coefficient for mean Earth longitude. */
49 private final int cE;
50
51 /** Coefficient for mean Mars longitude. */
52 private final int cMa;
53
54 /** Coefficient for mean Jupiter longitude. */
55 private final int cJu;
56
57 /** Coefficient for mean Saturn longitude. */
58 private final int cSa;
59
60 /** Coefficient for mean Uranus longitude. */
61 private final int cUr;
62
63 /** Coefficient for mean Neptune longitude. */
64 private final int cNe;
65
66 /** Coefficient for general accumulated precession in longitude. */
67 private final int cPa;
68
69 /** Build a general term for nutation series.
70 * @param cL coefficient for mean anomaly of the Moon
71 * @param cLPrime coefficient for mean anomaly of the Sun
72 * @param cF coefficient for L - Ω where L is the mean longitude of the Moon
73 * @param cD coefficient for mean elongation of the Moon from the Sun
74 * @param cOmega coefficient for mean longitude of the ascending node of the Moon
75 * @param cMe coefficient for mean Mercury longitude
76 * @param cVe coefficient for mean Venus longitude
77 * @param cE coefficient for mean Earth longitude
78 * @param cMa coefficient for mean Mars longitude
79 * @param cJu coefficient for mean Jupiter longitude
80 * @param cSa coefficient for mean Saturn longitude
81 * @param cUr coefficient for mean Uranus longitude
82 * @param cNe coefficient for mean Neptune longitude
83 * @param cPa coefficient for general accumulated precession in longitude
84 */
85 GeneralTerm(final int cL, final int cLPrime, final int cF, final int cD, final int cOmega,
86 final int cMe, final int cVe, final int cE, final int cMa, final int cJu,
87 final int cSa, final int cUr, final int cNe, final int cPa) {
88 this.cL = cL;
89 this.cLPrime = cLPrime;
90 this.cF = cF;
91 this.cD = cD;
92 this.cOmega = cOmega;
93 this.cMe = cMe;
94 this.cVe = cVe;
95 this.cE = cE;
96 this.cMa = cMa;
97 this.cJu = cJu;
98 this.cSa = cSa;
99 this.cUr = cUr;
100 this.cNe = cNe;
101 this.cPa = cPa;
102 }
103
104 /** {@inheritDoc} */
105 protected double argument(final BodiesElements elements) {
106 return cL * elements.getL() + cLPrime * elements.getLPrime() + cF * elements.getF() +
107 cD * elements.getD() + cOmega * elements.getOmega() +
108 cMe * elements.getLMe() + cVe * elements.getLVe() + cE * elements.getLE() +
109 cMa * elements.getLMa() + cJu * elements.getLJu() +
110 cSa * elements.getLSa() + cUr * elements.getLUr() +
111 cNe * elements.getLNe() + cPa * elements.getPa();
112 }
113
114 /** {@inheritDoc} */
115 protected double argumentDerivative(final BodiesElements elements) {
116 return cL * elements.getLDot() + cLPrime * elements.getLPrimeDot() + cF * elements.getFDot() +
117 cD * elements.getDDot() + cOmega * elements.getOmegaDot() +
118 cMe * elements.getLMeDot() + cVe * elements.getLVeDot() + cE * elements.getLEDot() +
119 cMa * elements.getLMaDot() + cJu * elements.getLJuDot() +
120 cSa * elements.getLSaDot() + cUr * elements.getLUrDot() +
121 cNe * elements.getLNeDot() + cPa * elements.getPaDot();
122 }
123
124 /** {@inheritDoc} */
125 protected <T extends RealFieldElement<T>> T argument(final FieldBodiesElements<T> elements) {
126 return elements.getL().multiply(cL).
127 add(elements.getLPrime().multiply(cLPrime)).
128 add(elements.getF().multiply(cF)).
129 add(elements.getD().multiply(cD)).
130 add(elements.getOmega().multiply(cOmega)).
131 add(elements.getLMe().multiply(cMe)).
132 add(elements.getLVe().multiply(cVe)).
133 add(elements.getLE().multiply(cE)).
134 add(elements.getLMa().multiply(cMa)).
135 add(elements.getLJu().multiply(cJu)).
136 add(elements.getLSa().multiply(cSa)).
137 add(elements.getLUr().multiply(cUr)).
138 add(elements.getLNe().multiply(cNe)).
139 add(elements.getPa().multiply(cPa));
140 }
141
142 /** {@inheritDoc} */
143 protected <T extends RealFieldElement<T>> T argumentDerivative(final FieldBodiesElements<T> elements) {
144 return elements.getLDot().multiply(cL).
145 add(elements.getLPrimeDot().multiply(cLPrime)).
146 add(elements.getFDot().multiply(cF)).
147 add(elements.getDDot().multiply(cD)).
148 add(elements.getOmegaDot().multiply(cOmega)).
149 add(elements.getLMeDot().multiply(cMe)).
150 add(elements.getLVeDot().multiply(cVe)).
151 add(elements.getLEDot().multiply(cE)).
152 add(elements.getLMaDot().multiply(cMa)).
153 add(elements.getLJuDot().multiply(cJu)).
154 add(elements.getLSaDot().multiply(cSa)).
155 add(elements.getLUrDot().multiply(cUr)).
156 add(elements.getLNeDot().multiply(cNe)).
157 add(elements.getPaDot().multiply(cPa));
158 }
159
160 }