1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.propagation.analytical;
18
19 import java.util.Arrays;
20 import java.util.List;
21
22 import org.hipparchus.analysis.differentiation.Gradient;
23 import org.hipparchus.linear.MatrixUtils;
24 import org.hipparchus.linear.RealMatrix;
25 import org.orekit.orbits.FieldOrbit;
26 import org.orekit.orbits.OrbitParamsType;
27 import org.orekit.orbits.PositionAngleType;
28 import org.orekit.propagation.AbstractMatricesHarvester;
29 import org.orekit.propagation.AdditionalDataProvider;
30 import org.orekit.propagation.FieldSpacecraftState;
31 import org.orekit.propagation.SpacecraftState;
32 import org.orekit.time.AbsoluteDate;
33 import org.orekit.time.FieldAbsoluteDate;
34 import org.orekit.utils.DoubleArrayDictionary;
35 import org.orekit.utils.drivers.ParameterDriver;
36
37
38
39
40
41
42
43
44 public abstract class AbstractAnalyticalMatricesHarvester extends AbstractMatricesHarvester
45 implements AdditionalDataProvider<double[]> {
46
47
48 private List<String> columnsNames;
49
50
51 private AbsoluteDate epoch;
52
53
54 private double[][] analyticalDerivativesStm;
55
56
57 private final DoubleArrayDictionary analyticalDerivativesJacobianColumns;
58
59
60 private final AbstractAnalyticalPropagator propagator;
61
62
63
64
65
66
67
68
69
70 protected AbstractAnalyticalMatricesHarvester(final AbstractAnalyticalPropagator propagator) {
71 this.propagator = propagator;
72 this.epoch = propagator.getInitialState().getDate();
73 this.columnsNames = null;
74 this.analyticalDerivativesJacobianColumns = new DoubleArrayDictionary();
75 }
76
77
78 @Override
79 protected void setInitialStm(final String stmName, final RealMatrix stm) {
80 super.setInitialStm(stmName, stm);
81 this.analyticalDerivativesStm = getInitialStateTransitionMatrix().getData();
82 }
83
84
85 @Override
86 public List<String> getJacobiansColumnsNames() {
87 return columnsNames == null ? propagator.getJacobiansColumnsNames() : columnsNames;
88 }
89
90
91 @Override
92 public void freezeColumnsNames() {
93 columnsNames = getJacobiansColumnsNames();
94 }
95
96
97 @Override
98 public String getName() {
99 return getStmName();
100 }
101
102
103 @Override
104 public double[] getAdditionalData(final SpacecraftState state) {
105
106 updateDerivativesIfNeeded(state);
107
108 return toArray(analyticalDerivativesStm);
109 }
110
111
112 @Override
113 public RealMatrix getStateTransitionMatrix(final SpacecraftState state) {
114
115 if (!state.hasAdditionalData(getName())) {
116 return null;
117 }
118
119 return toSquareMatrix(state.getAdditionalState(getName()));
120 }
121
122
123 @Override
124 public RealMatrix getParametersJacobian(final SpacecraftState state) {
125
126
127 updateDerivativesIfNeeded(state);
128
129
130 final List<String> names = getJacobiansColumnsNames();
131 if (names == null || names.isEmpty()) {
132 return null;
133 }
134
135
136 final RealMatrix dYdP = MatrixUtils.createRealMatrix(getStateDimension(), names.size());
137
138
139 for (int j = 0; j < names.size(); ++j) {
140 final double[] column = analyticalDerivativesJacobianColumns.get(names.get(j));
141 if (column != null) {
142 for (int i = 0; i < getStateDimension(); i++) {
143 dYdP.addToEntry(i, j, column[i]);
144 }
145 }
146 }
147
148
149 return dYdP;
150
151 }
152
153
154 @Override
155 public void setReferenceState(final SpacecraftState reference) {
156
157
158 for (final double[] row : analyticalDerivativesStm) {
159 Arrays.fill(row, 0.0);
160 }
161 analyticalDerivativesJacobianColumns.clear();
162
163 final AbstractAnalyticalGradientConverter converter = getGradientConverter();
164 final FieldAbstractAnalyticalPropagator<Gradient> gPropagator = converter.getPropagator();
165
166
167 final AbsoluteDate target = reference.getDate();
168 final FieldAbsoluteDate<Gradient> start = gPropagator.getInitialState().getDate();
169 final double dt = target.durationFrom(start.toAbsoluteDate());
170 final FieldSpacecraftState<Gradient> state = gPropagator.getInitialState();
171 final Gradient[] parameters = converter.getParameters(state, converter);
172 final FieldOrbit<Gradient> gOrbit = gPropagator.propagateOrbit(start.shiftedBy(dt), parameters);
173 final Gradient[] orbitDerivatives = new Gradient[6];
174 getOrbitParamsType().mapOrbitToArray(gOrbit, getPositionAngleType(), orbitDerivatives, null);
175
176
177 for (int i = 0; i < orbitDerivatives.length; ++i) {
178 addToRow(orbitDerivatives[i].getGradient(), i);
179 }
180
181
182 int paramsIndex = converter.getFreeStateParameters();
183 for (ParameterDriver driver : converter.getParametersDrivers()) {
184 if (driver.isSelected()) {
185
186
187 DoubleArrayDictionary.Entry entry = analyticalDerivativesJacobianColumns.getEntry(driver.getName());
188 if (entry == null) {
189
190 analyticalDerivativesJacobianColumns.put(driver.getName(), new double[getStateDimension()]);
191 entry = analyticalDerivativesJacobianColumns.getEntry(driver.getName());
192 }
193
194
195 entry.increment(new double[] {
196 orbitDerivatives[0].getPartialDerivative(paramsIndex),
197 orbitDerivatives[1].getPartialDerivative(paramsIndex),
198 orbitDerivatives[2].getPartialDerivative(paramsIndex),
199 orbitDerivatives[3].getPartialDerivative(paramsIndex),
200 orbitDerivatives[4].getPartialDerivative(paramsIndex),
201 orbitDerivatives[5].getPartialDerivative(paramsIndex)
202 });
203 ++paramsIndex;
204 }
205 }
206
207
208 epoch = target;
209
210 }
211
212
213
214
215 private void updateDerivativesIfNeeded(final SpacecraftState state) {
216 if (!state.getDate().isEqualTo(epoch)) {
217 setReferenceState(state);
218 }
219 }
220
221
222
223
224
225 private void addToRow(final double[] derivatives, final int index) {
226 for (int i = 0; i < 6; i++) {
227 analyticalDerivativesStm[index][i] += derivatives[i];
228 }
229 }
230
231
232 @Override
233 public OrbitParamsType getOrbitParamsType() {
234
235 return OrbitParamsType.CARTESIAN;
236 }
237
238
239 @Override
240 public PositionAngleType getPositionAngleType() {
241
242 return PositionAngleType.MEAN;
243 }
244
245
246
247
248
249 public abstract AbstractAnalyticalGradientConverter getGradientConverter();
250
251 }