GeographicLib 2.5
MagneticModel.hpp
Go to the documentation of this file.
1/**
2 * \file MagneticModel.hpp
3 * \brief Header for GeographicLib::MagneticModel class
4 *
5 * Copyright (c) Charles Karney (2011-2022) <karney@alum.mit.edu> and licensed
6 * under the MIT/X11 License. For more information, see
7 * https://geographiclib.sourceforge.io/
8 **********************************************************************/
9
10#if !defined(GEOGRAPHICLIB_MAGNETICMODEL_HPP)
11#define GEOGRAPHICLIB_MAGNETICMODEL_HPP 1
12
16
17#if defined(_MSC_VER)
18// Squelch warnings about dll vs vector
19# pragma warning (push)
20# pragma warning (disable: 4251)
21#endif
22
23namespace GeographicLib {
24
25 class MagneticCircle;
26
27 /**
28 * \brief Model of the earth's magnetic field
29 *
30 * Evaluate the earth's magnetic field according to a model. At present only
31 * internal magnetic fields are handled. These are due to the earth's code
32 * and crust; these vary slowly (over many years). Excluded are the effects
33 * of currents in the ionosphere and magnetosphere which have daily and
34 * annual variations.
35 *
36 * See \ref magnetic for details of how to install the magnetic models and
37 * the data format.
38 *
39 * See
40 * - General information:
41 * - http://geomag.org/models/index.html
42 * - WMM2010:
43 * - https://ngdc.noaa.gov/geomag/WMM/DoDWMM.shtml
44 * - https://ngdc.noaa.gov/geomag/WMM/data/WMM2010/WMM2010COF.zip
45 * - WMM2015 (deprecated):
46 * - https://ngdc.noaa.gov/geomag/WMM/DoDWMM.shtml
47 * - https://ngdc.noaa.gov/geomag/WMM/data/WMM2015/WMM2015COF.zip
48 * - WMM2015V2:
49 * - https://ngdc.noaa.gov/geomag/WMM/DoDWMM.shtml
50 * - https://ngdc.noaa.gov/geomag/WMM/data/WMM2015/WMM2015v2COF.zip
51 * - WMM2020:
52 * - https://ngdc.noaa.gov/geomag/WMM/DoDWMM.shtml
53 * - https://ngdc.noaa.gov/geomag/WMM/data/WMM2020/WMM2020COF.zip
54 * - WMM2025:
55 * - https://www.ncei.noaa.gov/products/world-magnetic-model
56 * - WMMHR2025:
57 * - https://www.ncei.noaa.gov/products/world-magnetic-model-high-resolution
58 * - IGRF11:
59 * - https://ngdc.noaa.gov/IAGA/vmod/igrf.html
60 * - https://ngdc.noaa.gov/IAGA/vmod/igrf11coeffs.txt
61 * - https://ngdc.noaa.gov/IAGA/vmod/geomag70_linux.tar.gz
62 * - EMM2010:
63 * - https://ngdc.noaa.gov/geomag/EMM/index.html
64 * - https://ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2010_Sph_Windows_Linux.zip
65 * - EMM2015:
66 * - https://ngdc.noaa.gov/geomag/EMM/index.html
67 * - https://www.ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2015_Sph_Linux.zip
68 * - EMM2017:
69 * - https://ngdc.noaa.gov/geomag/EMM/index.html
70 * - https://www.ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2017_Sph_Linux.zip
71 *
72 * Example of use:
73 * \include example-MagneticModel.cpp
74 *
75 * <a href="MagneticField.1.html">MagneticField</a> is a command-line utility
76 * providing access to the functionality of MagneticModel and MagneticCircle.
77 **********************************************************************/
78
80 private:
81 typedef Math::real real;
82 static const int idlength_ = 8;
83 std::string _name, _dir, _description, _date, _filename, _id;
84 real _t0, _dt0, _tmin, _tmax, _a, _hmin, _hmax;
85 int _nNmodels, _nNconstants, _nmx, _mmx;
87 Geocentric _earth;
88 std::vector< std::vector<real> > _gG;
89 std::vector< std::vector<real> > _hH;
90 std::vector<SphericalHarmonic> _harm;
91 void Field(real t, real lat, real lon, real h, bool diffp,
92 real& Bx, real& By, real& Bz,
93 real& Bxt, real& Byt, real& Bzt) const;
94 void ReadMetadata(const std::string& name);
95 // copy constructor not allowed
96 MagneticModel(const MagneticModel&) = delete;
97 // nor copy assignment
98 MagneticModel& operator=(const MagneticModel&) = delete;
99 public:
100
101 /**
102 * Move constructs a magnetic model.
103 **********************************************************************/
105
106 /**
107 * Move assigns a magnetic model.
108 **********************************************************************/
110
111 /** \name Setting up the magnetic model
112 **********************************************************************/
113 ///@{
114 /**
115 * Construct a magnetic model.
116 *
117 * @param[in] name the name of the model.
118 * @param[in] path (optional) directory for data file.
119 * @param[in] earth (optional) Geocentric object for converting
120 * coordinates; default Geocentric::WGS84().
121 * @param[in] Nmax (optional) if non-negative, truncate the degree of the
122 * model this value.
123 * @param[in] Mmax (optional) if non-negative, truncate the order of the
124 * model this value.
125 * @exception GeographicErr if the data file cannot be found, is
126 * unreadable, or is corrupt, or if \e Mmax > \e Nmax.
127 * @exception std::bad_alloc if the memory necessary for storing the model
128 * can't be allocated.
129 *
130 * A filename is formed by appending ".wmm" (World Magnetic Model) to the
131 * name. If \e path is specified (and is non-empty), then the file is
132 * loaded from directory, \e path. Otherwise the path is given by the
133 * DefaultMagneticPath().
134 *
135 * This file contains the metadata which specifies the properties of the
136 * model. The coefficients for the spherical harmonic sums are obtained
137 * from a file obtained by appending ".cof" to metadata file (so the
138 * filename ends in ".wwm.cof").
139 *
140 * The model is not tied to a particular ellipsoidal model of the earth.
141 * The final earth argument to the constructor specifies an ellipsoid to
142 * allow geodetic coordinates to the transformed into the spherical
143 * coordinates used in the spherical harmonic sum.
144 *
145 * If \e Nmax &ge; 0 and \e Mmax < 0, then \e Mmax is set to \e Nmax.
146 * After the model is loaded, the maximum degree and order of the model can
147 * be found by the Degree() and Order() methods.
148 **********************************************************************/
149 explicit MagneticModel(const std::string& name,
150 const std::string& path = "",
151 const Geocentric& earth = Geocentric::WGS84(),
152 int Nmax = -1, int Mmax = -1);
153 ///@}
154
155 /** \name Compute the magnetic field
156 **********************************************************************/
157 ///@{
158 /**
159 * Evaluate the components of the geomagnetic field.
160 *
161 * @param[in] t the time (fractional years).
162 * @param[in] lat latitude of the point (degrees).
163 * @param[in] lon longitude of the point (degrees).
164 * @param[in] h the height of the point above the ellipsoid (meters).
165 * @param[out] Bx the easterly component of the magnetic field (nanotesla).
166 * @param[out] By the northerly component of the magnetic field
167 * (nanotesla).
168 * @param[out] Bz the vertical (up) component of the magnetic field
169 * (nanotesla).
170 *
171 * Use Utility::fractionalyear to convert a date of the form yyyy-mm or
172 * yyyy-mm-dd into a fractional year.
173 **********************************************************************/
174 void operator()(real t, real lat, real lon, real h,
175 real& Bx, real& By, real& Bz) const {
176 real dummy;
177 Field(t, lat, lon, h, false, Bx, By, Bz, dummy, dummy, dummy);
178 }
179
180 /**
181 * Evaluate the components of the geomagnetic field and their time
182 * derivatives
183 *
184 * @param[in] t the time (fractional years).
185 * @param[in] lat latitude of the point (degrees).
186 * @param[in] lon longitude of the point (degrees).
187 * @param[in] h the height of the point above the ellipsoid (meters).
188 * @param[out] Bx the easterly component of the magnetic field (nanotesla).
189 * @param[out] By the northerly component of the magnetic field
190 * (nanotesla).
191 * @param[out] Bz the vertical (up) component of the magnetic field
192 * (nanotesla).
193 * @param[out] Bxt the rate of change of \e Bx (nT/yr).
194 * @param[out] Byt the rate of change of \e By (nT/yr).
195 * @param[out] Bzt the rate of change of \e Bz (nT/yr).
196 *
197 * Use Utility::fractionalyear to convert a date of the form yyyy-mm or
198 * yyyy-mm-dd into a fractional year.
199 **********************************************************************/
200 void operator()(real t, real lat, real lon, real h,
201 real& Bx, real& By, real& Bz,
202 real& Bxt, real& Byt, real& Bzt) const {
203 Field(t, lat, lon, h, true, Bx, By, Bz, Bxt, Byt, Bzt);
204 }
205
206 /**
207 * Create a MagneticCircle object to allow the geomagnetic field at many
208 * points with constant \e lat, \e h, and \e t and varying \e lon to be
209 * computed efficiently.
210 *
211 * @param[in] t the time (fractional years).
212 * @param[in] lat latitude of the point (degrees).
213 * @param[in] h the height of the point above the ellipsoid (meters).
214 * @exception std::bad_alloc if the memory necessary for creating a
215 * MagneticCircle can't be allocated.
216 * @return a MagneticCircle object whose MagneticCircle::operator()(real
217 * lon) member function computes the field at particular values of \e
218 * lon.
219 *
220 * If the field at several points on a circle of latitude need to be
221 * calculated then creating a MagneticCircle and using its member functions
222 * will be substantially faster, especially for high-degree models.
223 *
224 * Use Utility::fractionalyear to convert a date of the form yyyy-mm or
225 * yyyy-mm-dd into a fractional year.
226 **********************************************************************/
227 MagneticCircle Circle(real t, real lat, real h) const;
228
229 /**
230 * Compute the magnetic field in geocentric coordinate.
231 *
232 * @param[in] t the time (fractional years).
233 * @param[in] X geocentric coordinate (meters).
234 * @param[in] Y geocentric coordinate (meters).
235 * @param[in] Z geocentric coordinate (meters).
236 * @param[out] BX the \e X component of the magnetic field (nT).
237 * @param[out] BY the \e Y component of the magnetic field (nT).
238 * @param[out] BZ the \e Z component of the magnetic field (nT).
239 * @param[out] BXt the rate of change of \e BX (nT/yr).
240 * @param[out] BYt the rate of change of \e BY (nT/yr).
241 * @param[out] BZt the rate of change of \e BZ (nT/yr).
242 *
243 * Use Utility::fractionalyear to convert a date of the form yyyy-mm or
244 * yyyy-mm-dd into a fractional year.
245 **********************************************************************/
246 void FieldGeocentric(real t, real X, real Y, real Z,
247 real& BX, real& BY, real& BZ,
248 real& BXt, real& BYt, real& BZt) const;
249
250 /**
251 * Compute various quantities dependent on the magnetic field.
252 *
253 * @param[in] Bx the \e x (easterly) component of the magnetic field (nT).
254 * @param[in] By the \e y (northerly) component of the magnetic field (nT).
255 * @param[in] Bz the \e z (vertical, up positive) component of the magnetic
256 * field (nT).
257 * @param[out] H the horizontal magnetic field (nT).
258 * @param[out] F the total magnetic field (nT).
259 * @param[out] D the declination of the field (degrees east of north).
260 * @param[out] I the inclination of the field (degrees down from
261 * horizontal).
262 **********************************************************************/
263 static void FieldComponents(real Bx, real By, real Bz,
264 real& H, real& F, real& D, real& I) {
265 real Ht, Ft, Dt, It;
266 FieldComponents(Bx, By, Bz, real(0), real(1), real(0),
267 H, F, D, I, Ht, Ft, Dt, It);
268 }
269
270 /**
271 * Compute various quantities dependent on the magnetic field and its rate
272 * of change.
273 *
274 * @param[in] Bx the \e x (easterly) component of the magnetic field (nT).
275 * @param[in] By the \e y (northerly) component of the magnetic field (nT).
276 * @param[in] Bz the \e z (vertical, up positive) component of the magnetic
277 * field (nT).
278 * @param[in] Bxt the rate of change of \e Bx (nT/yr).
279 * @param[in] Byt the rate of change of \e By (nT/yr).
280 * @param[in] Bzt the rate of change of \e Bz (nT/yr).
281 * @param[out] H the horizontal magnetic field (nT).
282 * @param[out] F the total magnetic field (nT).
283 * @param[out] D the declination of the field (degrees east of north).
284 * @param[out] I the inclination of the field (degrees down from
285 * horizontal).
286 * @param[out] Ht the rate of change of \e H (nT/yr).
287 * @param[out] Ft the rate of change of \e F (nT/yr).
288 * @param[out] Dt the rate of change of \e D (degrees/yr).
289 * @param[out] It the rate of change of \e I (degrees/yr).
290 **********************************************************************/
291 static void FieldComponents(real Bx, real By, real Bz,
292 real Bxt, real Byt, real Bzt,
293 real& H, real& F, real& D, real& I,
294 real& Ht, real& Ft, real& Dt, real& It);
295 ///@}
296
297 /** \name Inspector functions
298 **********************************************************************/
299 ///@{
300 /**
301 * @return the description of the magnetic model, if available, from the
302 * Description file in the data file; if absent, return "NONE".
303 **********************************************************************/
304 const std::string& Description() const { return _description; }
305
306 /**
307 * @return date of the model, if available, from the ReleaseDate field in
308 * the data file; if absent, return "UNKNOWN".
309 **********************************************************************/
310 const std::string& DateTime() const { return _date; }
311
312 /**
313 * @return full file name used to load the magnetic model.
314 **********************************************************************/
315 const std::string& MagneticFile() const { return _filename; }
316
317 /**
318 * @return "name" used to load the magnetic model (from the first argument
319 * of the constructor, but this may be overridden by the model file).
320 **********************************************************************/
321 const std::string& MagneticModelName() const { return _name; }
322
323 /**
324 * @return directory used to load the magnetic model.
325 **********************************************************************/
326 const std::string& MagneticModelDirectory() const { return _dir; }
327
328 /**
329 * @return the minimum height above the ellipsoid (in meters) for which
330 * this MagneticModel should be used.
331 *
332 * Because the model will typically provide useful results
333 * slightly outside the range of allowed heights, no check of \e t
334 * argument is made by MagneticModel::operator()() or
335 * MagneticModel::Circle.
336 **********************************************************************/
337 Math::real MinHeight() const { return _hmin; }
338
339 /**
340 * @return the maximum height above the ellipsoid (in meters) for which
341 * this MagneticModel should be used.
342 *
343 * Because the model will typically provide useful results
344 * slightly outside the range of allowed heights, no check of \e t
345 * argument is made by MagneticModel::operator()() or
346 * MagneticModel::Circle.
347 **********************************************************************/
348 Math::real MaxHeight() const { return _hmax; }
349
350 /**
351 * @return the minimum time (in years) for which this MagneticModel should
352 * be used.
353 *
354 * Because the model will typically provide useful results
355 * slightly outside the range of allowed times, no check of \e t
356 * argument is made by MagneticModel::operator()() or
357 * MagneticModel::Circle.
358 **********************************************************************/
359 Math::real MinTime() const { return _tmin; }
360
361 /**
362 * @return the maximum time (in years) for which this MagneticModel should
363 * be used.
364 *
365 * Because the model will typically provide useful results
366 * slightly outside the range of allowed times, no check of \e t
367 * argument is made by MagneticModel::operator()() or
368 * MagneticModel::Circle.
369 **********************************************************************/
370 Math::real MaxTime() const { return _tmax; }
371
372 /**
373 * @return \e a the equatorial radius of the ellipsoid (meters). This is
374 * the value of \e a inherited from the Geocentric object used in the
375 * constructor.
376 **********************************************************************/
377 Math::real EquatorialRadius() const { return _earth.EquatorialRadius(); }
378
379 /**
380 * @return \e f the flattening of the ellipsoid. This is the value
381 * inherited from the Geocentric object used in the constructor.
382 **********************************************************************/
383 Math::real Flattening() const { return _earth.Flattening(); }
384
385 /**
386 * @return \e Nmax the maximum degree of the components of the model.
387 **********************************************************************/
388 int Degree() const { return _nmx; }
389
390 /**
391 * @return \e Mmax the maximum order of the components of the model.
392 **********************************************************************/
393 int Order() const { return _mmx; }
394 ///@}
395
396 /**
397 * @return the default path for magnetic model data files.
398 *
399 * This is the value of the environment variable
400 * GEOGRAPHICLIB_MAGNETIC_PATH, if set; otherwise, it is
401 * $GEOGRAPHICLIB_DATA/magnetic if the environment variable
402 * GEOGRAPHICLIB_DATA is set; otherwise, it is a compile-time default
403 * (/usr/local/share/GeographicLib/magnetic on non-Windows systems and
404 * C:/ProgramData/GeographicLib/magnetic on Windows systems).
405 **********************************************************************/
406 static std::string DefaultMagneticPath();
407
408 /**
409 * @return the default name for the magnetic model.
410 *
411 * This is the value of the environment variable
412 * GEOGRAPHICLIB_MAGNETIC_NAME, if set; otherwise, it is "wmm2025". The
413 * MagneticModel class does not use this function; it is just provided as a
414 * convenience for a calling program when constructing a MagneticModel
415 * object.
416 **********************************************************************/
417 static std::string DefaultMagneticName();
418 };
419
420} // namespace GeographicLib
421
422#if defined(_MSC_VER)
423# pragma warning (pop)
424#endif
425
426#endif // GEOGRAPHICLIB_MAGNETICMODEL_HPP
Header for GeographicLib::Constants class.
#define GEOGRAPHICLIB_EXPORT
Definition Constants.hpp:67
Header for GeographicLib::Geocentric class.
GeographicLib::Math::real real
Definition GeodSolve.cpp:28
Header for GeographicLib::SphericalHarmonic class.
Geocentric coordinates
Math::real Flattening() const
Math::real EquatorialRadius() const
Geomagnetic field on a circle of latitude.
Model of the earth's magnetic field.
Math::real EquatorialRadius() const
void operator()(real t, real lat, real lon, real h, real &Bx, real &By, real &Bz) const
const std::string & Description() const
const std::string & MagneticFile() const
const std::string & DateTime() const
const std::string & MagneticModelName() const
MagneticModel(MagneticModel &&)=default
const std::string & MagneticModelDirectory() const
MagneticModel & operator=(MagneticModel &&)=default
void operator()(real t, real lat, real lon, real h, real &Bx, real &By, real &Bz, real &Bxt, real &Byt, real &Bzt) const
static void FieldComponents(real Bx, real By, real Bz, real &H, real &F, real &D, real &I)
Namespace for GeographicLib.