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  
18  package org.orekit.propagation.events;
19  
20  import org.orekit.annotation.DefaultDataContext;
21  import org.orekit.bodies.GeodeticPoint;
22  import org.orekit.bodies.OneAxisEllipsoid;
23  import org.orekit.data.DataContext;
24  import org.orekit.models.earth.GeoMagneticField;
25  import org.orekit.models.earth.GeoMagneticFieldFactory;
26  import org.orekit.models.earth.GeoMagneticFieldFactory.FieldModel;
27  import org.orekit.propagation.SpacecraftState;
28  import org.orekit.propagation.events.handlers.EventHandler;
29  import org.orekit.propagation.events.handlers.StopOnIncreasing;
30  import org.orekit.time.AbsoluteDate;
31  import org.orekit.time.TimeScale;
32  
33  /** Detector for Earth magnetic field strength.
34   * <p>
35   * The detector is based on the field intensity calculated at the
36   * satellite's latitude and longitude, either at sea level or at
37   * satellite altitude, depending on the value chosen for the
38   * <code>atSeaLevel</code> indicator.<br>
39   * It can detect flyovers of the South-Atlantic anomaly with
40   * a classically accepted limit value of 32,000 nT at sea level.
41   * </p>
42   * @author Romaric Her
43   */
44  public class MagneticFieldDetector extends AbstractDetector<MagneticFieldDetector> {
45  
46      /** Fixed threshold value of Magnetic field to be crossed, in Teslas. */
47      private final double limit;
48  
49      /** Switch for calculating field strength at sea level (true) or satellite altitude (false). */
50      private final boolean atSeaLevel;
51  
52      /** Earth geomagnetic field. */
53      private GeoMagneticField field;
54  
55      /** year of the current state. */
56      private double currentYear;
57  
58      /** Earth geomagnetic field model. */
59      private final FieldModel model;
60  
61      /** Earth body shape. */
62      private final OneAxisEllipsoid body;
63  
64      /** Current data context. */
65      private final DataContext dataContext;
66  
67  
68      /** Build a new detector.
69       *
70       * <p>This constructor uses:
71       * <ul>
72       * <li>the {@link DataContext#getDefault() default data context}</li>
73       * <li>the {@link AbstractDetector#DEFAULT_MAX_CHECK default value} for maximal checking interval</li>
74       * <li>the {@link AbstractDetector#DEFAULT_THRESHOLD default value} for convergence threshold</li>
75       * <li>the <code>atSeaLevel</code> switch set to false</li>
76       * </ul>
77       *
78       * @param limit threshold value for magnetic field detection, in Teslas
79       * @param model magnetic field model
80       * @param body  Earth body shape
81       * @see #MagneticFieldDetector(double, double, double, GeoMagneticFieldFactory.FieldModel, OneAxisEllipsoid, boolean, DataContext)
82       */
83      @DefaultDataContext
84      public MagneticFieldDetector(final double limit, final FieldModel model, final OneAxisEllipsoid body) {
85          this(DEFAULT_MAX_CHECK, DEFAULT_THRESHOLD, limit, model, body, false);
86      }
87  
88      /** Build a new detector.
89       *
90       * <p>This constructor uses:
91       * <ul>
92       * <li>the {@link DataContext#getDefault() default data context}</li>
93       * <li>the {@link AbstractDetector#DEFAULT_MAX_CHECK default value} for maximal checking interval</li>
94       * <li>the {@link AbstractDetector#DEFAULT_THRESHOLD default value} for convergence threshold </li>
95       * </ul>
96       *
97       * @param limit    threshold value for magnetic field detection, in Teslas
98       * @param model    magnetic field model
99       * @param body     Earth body shape
100      * @param atSeaLevel switch for calculating field intensity at sea level (true) or satellite altitude (false)
101      * @see #MagneticFieldDetector(double, double, double, GeoMagneticFieldFactory.FieldModel, OneAxisEllipsoid, boolean, DataContext)
102      */
103     @DefaultDataContext
104     public MagneticFieldDetector(final double limit, final FieldModel model,
105                                  final OneAxisEllipsoid body, final boolean atSeaLevel) {
106         this(DEFAULT_MAX_CHECK, DEFAULT_THRESHOLD, limit, model, body, atSeaLevel);
107     }
108 
109     /** Build a detector.
110      *
111      * <p>This method uses the {@link DataContext#getDefault() default data context}.</p>
112      *
113      * @param maxCheck   maximal checking interval (s)
114      * @param threshold  convergence threshold (s)
115      * @param limit      threshold value for magnetic field detection, in Teslas
116      * @param model      magnetic field model
117      * @param body       Earth body shape
118      * @param atSeaLevel switch for calculating field intensity at sea level (true) or satellite altitude (false)
119      * @see #MagneticFieldDetector(double, double, double, GeoMagneticFieldFactory.FieldModel, OneAxisEllipsoid, boolean, DataContext)
120      */
121     @DefaultDataContext
122     public MagneticFieldDetector(final double maxCheck, final double threshold, final double limit,
123                                  final FieldModel model, final OneAxisEllipsoid body, final boolean atSeaLevel) {
124         this(maxCheck, threshold, limit, model, body, atSeaLevel, DataContext.getDefault());
125     }
126 
127     /**
128      * Build a detector.
129      *
130      * @param maxCheck    maximal checking interval (s)
131      * @param threshold   convergence threshold (s)
132      * @param limit       threshold value for magnetic field detection, in Teslas
133      * @param model       magnetic field model
134      * @param body        Earth body shape
135      * @param atSeaLevel  switch for calculating field intensity at sea level (true) or satellite altitude (false)
136      * @param dataContext used to look up the magnetic field model.
137      * @since 10.1
138      */
139     public MagneticFieldDetector(final double maxCheck,
140                                  final double threshold,
141                                  final double limit,
142                                  final FieldModel model,
143                                  final OneAxisEllipsoid body,
144                                  final boolean atSeaLevel,
145                                  final DataContext dataContext) {
146         this(new EventDetectionSettings(maxCheck, threshold, DEFAULT_MAX_ITER), new StopOnIncreasing(),
147              limit, model, body, atSeaLevel, dataContext);
148     }
149 
150     /** Protected constructor with full parameters.
151      * <p>
152      * This constructor is not public as users are expected to use the builder
153      * API with the various {@code withXxx()} methods to set up the instance
154      * in a readable manner without using a huge amount of parameters.
155      * </p>
156      * @param detectionSettings event detection settings
157      * @param handler     event handler to call at event occurrences
158      * @param limit       threshold value for magnetic field detection, in Teslas
159      * @param model       magnetic field model
160      * @param body        Earth body shape
161      * @param atSeaLevel  switch for calculating field intensity at sea level (true) or satellite altitude (false)
162      * @param dataContext used to look up the magnetic field model.
163      * @since 13.0
164      */
165     protected MagneticFieldDetector(final EventDetectionSettings detectionSettings, final EventHandler handler,
166                                     final double limit, final FieldModel model, final OneAxisEllipsoid body,
167                                     final boolean atSeaLevel, final DataContext dataContext) {
168         super(detectionSettings, handler);
169         this.limit       = limit;
170         this.model       = model;
171         this.body        = body;
172         this.atSeaLevel  = atSeaLevel;
173         this.dataContext = dataContext;
174     }
175 
176     /** {@inheritDoc} */
177     @Override
178     protected MagneticFieldDetector create(final EventDetectionSettings detectionSettings, final EventHandler newHandler) {
179         return new MagneticFieldDetector(detectionSettings, newHandler, limit, model, body, atSeaLevel, dataContext);
180     }
181 
182     /** {@inheritDoc} */
183     @Override
184     public void init(final SpacecraftState s0, final AbsoluteDate t) {
185         super.init(s0, t);
186         final TimeScale utc = dataContext.getTimeScales().getUTC();
187         this.currentYear = s0.getDate().getComponents(utc).getDate().getYear();
188         this.field = dataContext.getGeoMagneticFields().getField(model, currentYear);
189     }
190 
191     /** Compute the value of the detection function.
192      * <p>
193      * The returned value is the difference between the field intensity at spacecraft location,
194      * taking <code>atSeaLevel</code> switch into account, and the fixed threshold value.
195      * </p>
196      * @param s the current state information: date, kinematics, attitude
197      * @return difference between the field intensity at spacecraft location
198      *         and the fixed threshold value
199      */
200     public double g(final SpacecraftState s) {
201         final TimeScale utc = dataContext.getTimeScales().getUTC();
202         if (s.getDate().getComponents(utc).getDate().getYear() != currentYear) {
203             this.currentYear = s.getDate().getComponents(utc).getDate().getYear();
204             this.field = dataContext.getGeoMagneticFields().getField(model, currentYear);
205         }
206         final GeodeticPoint geoPoint = body.transform(s.getPosition(), s.getFrame(), s.getDate());
207         final double altitude = atSeaLevel ? 0. : geoPoint.getAltitude();
208         final double value = field.calculateField(geoPoint.getLatitude(), geoPoint.getLongitude(), altitude).getTotalIntensity();
209         return value - limit;
210     }
211 
212 }