1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17 package org.orekit.propagation.semianalytical.dsst;
18
19 import java.io.IOException;
20 import java.text.ParseException;
21 import java.util.ArrayList;
22 import java.util.Arrays;
23 import java.util.List;
24
25 import org.hipparchus.geometry.euclidean.threed.Rotation;
26 import org.hipparchus.geometry.euclidean.threed.RotationOrder;
27 import org.hipparchus.geometry.euclidean.threed.Vector3D;
28 import org.junit.Assert;
29 import org.junit.Before;
30 import org.junit.Test;
31 import org.orekit.Utils;
32 import org.orekit.attitudes.Attitude;
33 import org.orekit.attitudes.AttitudeProvider;
34 import org.orekit.attitudes.InertialProvider;
35 import org.orekit.attitudes.LofOffset;
36 import org.orekit.bodies.CelestialBody;
37 import org.orekit.bodies.CelestialBodyFactory;
38 import org.orekit.bodies.OneAxisEllipsoid;
39 import org.orekit.errors.OrekitException;
40 import org.orekit.forces.BoxAndSolarArraySpacecraft;
41 import org.orekit.forces.gravity.potential.GravityFieldFactory;
42 import org.orekit.forces.gravity.potential.UnnormalizedSphericalHarmonicsProvider;
43 import org.orekit.frames.Frame;
44 import org.orekit.frames.FramesFactory;
45 import org.orekit.frames.LOFType;
46 import org.orekit.models.earth.atmosphere.Atmosphere;
47 import org.orekit.models.earth.atmosphere.HarrisPriester;
48 import org.orekit.orbits.EquinoctialOrbit;
49 import org.orekit.orbits.Orbit;
50 import org.orekit.orbits.PositionAngle;
51 import org.orekit.propagation.PropagationType;
52 import org.orekit.propagation.SpacecraftState;
53 import org.orekit.propagation.semianalytical.dsst.forces.DSSTAtmosphericDrag;
54 import org.orekit.propagation.semianalytical.dsst.forces.DSSTForceModel;
55 import org.orekit.propagation.semianalytical.dsst.forces.ShortPeriodTerms;
56 import org.orekit.propagation.semianalytical.dsst.utilities.AuxiliaryElements;
57 import org.orekit.time.AbsoluteDate;
58 import org.orekit.time.DateComponents;
59 import org.orekit.time.TimeComponents;
60 import org.orekit.time.TimeScalesFactory;
61 import org.orekit.utils.Constants;
62 import org.orekit.utils.IERSConventions;
63 import org.orekit.utils.TimeStampedAngularCoordinates;
64
65 public class DSSTAtmosphericDragTest {
66
67 @Test
68 public void testGetMeanElementRate() throws IllegalArgumentException, OrekitException {
69
70
71 final UnnormalizedSphericalHarmonicsProvider provider =
72 GravityFieldFactory.getUnnormalizedProvider(2, 0);
73
74 final Frame earthFrame = FramesFactory.getEME2000();
75 final AbsoluteDate initDate = new AbsoluteDate(2003, 07, 01, 0, 0, 0, TimeScalesFactory.getUTC());
76 final double mu = 3.986004415E14;
77
78
79
80
81
82
83 final Orbit orbit = new EquinoctialOrbit(7204535.84810944,
84 -0.001119677138261611,
85 5.333650671984143E-4,
86 0.847841707880348,
87 0.7998014061193262,
88 3.897842092486239,
89 PositionAngle.TRUE,
90 earthFrame,
91 initDate,
92 mu);
93
94
95 final OneAxisEllipsoid earth = new OneAxisEllipsoid(provider.getAe(),
96 Constants.WGS84_EARTH_FLATTENING,
97 CelestialBodyFactory.getEarth().getBodyOrientedFrame());
98 final Atmosphere atm = new HarrisPriester(CelestialBodyFactory.getSun(), earth, 6);
99 final double cd = 2.0;
100 final double area = 25.0;
101 DSSTForceModel drag = new DSSTAtmosphericDrag(atm, cd, area, mu);
102
103
104 Rotation rotation = new Rotation(1.0, 0.0, 0.0, 0.0, false);
105 Vector3D rotationRate = new Vector3D(0., 0., 0.);
106 Vector3D rotationAcceleration = new Vector3D(0., 0., 0.);
107 TimeStampedAngularCoordinates orientation = new TimeStampedAngularCoordinates(initDate, rotation, rotationRate, rotationAcceleration);
108 final Attitude att = new Attitude(earthFrame, orientation);
109
110 final SpacecraftState state = new SpacecraftState(orbit, att, 1000.0);
111 final AuxiliaryElements auxiliaryElements = new AuxiliaryElements(state.getOrbit(), 1);
112
113
114 final double[] parameters = drag.getParameters();
115
116 drag.initializeShortPeriodTerms(auxiliaryElements, PropagationType.MEAN, parameters);
117
118
119 AttitudeProvider attitudeProvider = new InertialProvider(rotation);
120 drag.registerAttitudeProvider(attitudeProvider );
121
122
123 final double[] elements = new double[7];
124 Arrays.fill(elements, 0.0);
125 final double[] daidt = drag.getMeanElementRate(state, auxiliaryElements, parameters);
126 for (int i = 0; i < daidt.length; i++) {
127 elements[i] = daidt[i];
128 }
129
130 Assert.assertEquals(-3.415320567871035E-5, elements[0], 1.e-20);
131 Assert.assertEquals(6.276312897745139E-13, elements[1], 1.9e-27);
132 Assert.assertEquals(-9.303357008691404E-13, elements[2], 0.7e-27);
133 Assert.assertEquals(-7.052316604063199E-14, elements[3], 1.e-28);
134 Assert.assertEquals(-6.793277250493389E-14, elements[4], 3.e-29);
135 Assert.assertEquals(-1.3565284454826392E-15, elements[5], 1.e-27);
136
137 }
138
139 @Test
140 public void testShortPeriodTerms() throws IllegalArgumentException, OrekitException {
141
142 final AbsoluteDate initDate = new AbsoluteDate(new DateComponents(2003, 03, 21), new TimeComponents(1, 0, 0.), TimeScalesFactory.getUTC());
143
144 final Orbit orbit = new EquinoctialOrbit(7069219.9806427825,
145 -4.5941811292223825E-4,
146 1.309932339472599E-4,
147 -1.002996107003202,
148 0.570979900577994,
149 2.62038786211518,
150 PositionAngle.TRUE,
151 FramesFactory.getEME2000(),
152 initDate,
153 3.986004415E14);
154
155 final SpacecraftState meanState = new SpacecraftState(orbit);
156
157 final CelestialBody sun = CelestialBodyFactory.getSun();
158 final OneAxisEllipsoid earth = new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
159 Constants.WGS84_EARTH_FLATTENING,
160 FramesFactory.getITRF(IERSConventions.IERS_2010,
161 true));
162
163 final BoxAndSolarArraySpacecraft boxAndWing = new BoxAndSolarArraySpacecraft(5.0, 2.0, 2.0,
164 sun,
165 50.0, Vector3D.PLUS_J,
166 2.0, 0.1,
167 0.2, 0.6);
168
169 final Atmosphere atmosphere = new HarrisPriester(CelestialBodyFactory.getSun(), earth, 6);
170 final AttitudeProvider attitudeProvider = new LofOffset(meanState.getFrame(),
171 LOFType.LVLH_CCSDS, RotationOrder.XYZ,
172 0.0, 0.0, 0.0);
173
174 final DSSTForceModel drag = new DSSTAtmosphericDrag(atmosphere, boxAndWing, meanState.getMu());
175
176
177 final AuxiliaryElements aux = new AuxiliaryElements(meanState.getOrbit(), 1);
178
179
180 final List<ShortPeriodTerms> shortPeriodTerms = new ArrayList<ShortPeriodTerms>();
181
182 drag.registerAttitudeProvider(attitudeProvider);
183 shortPeriodTerms.addAll(drag.initializeShortPeriodTerms(aux, PropagationType.OSCULATING, drag.getParameters()));
184 drag.updateShortPeriodTerms(drag.getParameters(), meanState);
185
186 double[] y = new double[6];
187 for (final ShortPeriodTerms spt : shortPeriodTerms) {
188 final double[] shortPeriodic = spt.value(meanState.getOrbit());
189 for (int i = 0; i < shortPeriodic.length; i++) {
190 y[i] += shortPeriodic[i];
191 }
192 }
193
194 Assert.assertEquals(0.03966657233280967, y[0], 1.e-15);
195 Assert.assertEquals(-1.5294381443173415E-8, y[1], 1.e-23);
196 Assert.assertEquals(-2.3614929828516364E-8, y[2], 1.e-23);
197 Assert.assertEquals(-5.901580336558653E-11, y[3], 1.e-26);
198 Assert.assertEquals(1.0287639743124977E-11, y[4], 1.e-26);
199 Assert.assertEquals(2.538427523777691E-8, y[5], 1.e-23);
200 }
201
202 @Before
203 public void setUp() throws IOException, ParseException {
204 Utils.setDataRoot("regular-data:potential/shm-format");
205 }
206 }