GeographicLib  1.50.1
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-2019) <charles@karney.com> 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 
23 namespace 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  * - IGRF11:
55  * - https://ngdc.noaa.gov/IAGA/vmod/igrf.html
56  * - https://ngdc.noaa.gov/IAGA/vmod/igrf11coeffs.txt
57  * - https://ngdc.noaa.gov/IAGA/vmod/geomag70_linux.tar.gz
58  * - EMM2010:
59  * - https://ngdc.noaa.gov/geomag/EMM/index.html
60  * - https://ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2010_Sph_Windows_Linux.zip
61  * - EMM2015:
62  * - https://ngdc.noaa.gov/geomag/EMM/index.html
63  * - https://www.ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2015_Sph_Linux.zip
64  * - EMM2017:
65  * - https://ngdc.noaa.gov/geomag/EMM/index.html
66  * - https://www.ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2017_Sph_Linux.zip
67  *
68  * Example of use:
69  * \include example-MagneticModel.cpp
70  *
71  * <a href="MagneticField.1.html">MagneticField</a> is a command-line utility
72  * providing access to the functionality of MagneticModel and MagneticCircle.
73  **********************************************************************/
74 
76  private:
77  typedef Math::real real;
78  static const int idlength_ = 8;
79  std::string _name, _dir, _description, _date, _filename, _id;
80  real _t0, _dt0, _tmin, _tmax, _a, _hmin, _hmax;
81  int _Nmodels, _Nconstants, _nmx, _mmx;
83  Geocentric _earth;
84  std::vector< std::vector<real> > _G;
85  std::vector< std::vector<real> > _H;
86  std::vector<SphericalHarmonic> _harm;
87  void Field(real t, real lat, real lon, real h, bool diffp,
88  real& Bx, real& By, real& Bz,
89  real& Bxt, real& Byt, real& Bzt) const;
90  void ReadMetadata(const std::string& name);
91  MagneticModel(const MagneticModel&); // copy constructor not allowed
92  MagneticModel& operator=(const MagneticModel&); // nor copy assignment
93  public:
94 
95  /** \name Setting up the magnetic model
96  **********************************************************************/
97  ///@{
98  /**
99  * Construct a magnetic model.
100  *
101  * @param[in] name the name of the model.
102  * @param[in] path (optional) directory for data file.
103  * @param[in] earth (optional) Geocentric object for converting
104  * coordinates; default Geocentric::WGS84().
105  * @param[in] Nmax (optional) if non-negative, truncate the degree of the
106  * model this value.
107  * @param[in] Mmax (optional) if non-negative, truncate the order of the
108  * model this value.
109  * @exception GeographicErr if the data file cannot be found, is
110  * unreadable, or is corrupt, or if \e Mmax > \e Nmax.
111  * @exception std::bad_alloc if the memory necessary for storing the model
112  * can't be allocated.
113  *
114  * A filename is formed by appending ".wmm" (World Magnetic Model) to the
115  * name. If \e path is specified (and is non-empty), then the file is
116  * loaded from directory, \e path. Otherwise the path is given by the
117  * DefaultMagneticPath().
118  *
119  * This file contains the metadata which specifies the properties of the
120  * model. The coefficients for the spherical harmonic sums are obtained
121  * from a file obtained by appending ".cof" to metadata file (so the
122  * filename ends in ".wwm.cof").
123  *
124  * The model is not tied to a particular ellipsoidal model of the earth.
125  * The final earth argument to the constructor specifies an ellipsoid to
126  * allow geodetic coordinates to the transformed into the spherical
127  * coordinates used in the spherical harmonic sum.
128  *
129  * If \e Nmax &ge; 0 and \e Mmax < 0, then \e Mmax is set to \e Nmax.
130  * After the model is loaded, the maximum degree and order of the model can
131  * be found by the Degree() and Order() methods.
132  **********************************************************************/
133  explicit MagneticModel(const std::string& name,
134  const std::string& path = "",
135  const Geocentric& earth = Geocentric::WGS84(),
136  int Nmax = -1, int Mmax = -1);
137  ///@}
138 
139  /** \name Compute the magnetic field
140  **********************************************************************/
141  ///@{
142  /**
143  * Evaluate the components of the geomagnetic field.
144  *
145  * @param[in] t the time (years).
146  * @param[in] lat latitude of the point (degrees).
147  * @param[in] lon longitude of the point (degrees).
148  * @param[in] h the height of the point above the ellipsoid (meters).
149  * @param[out] Bx the easterly component of the magnetic field (nanotesla).
150  * @param[out] By the northerly component of the magnetic field
151  * (nanotesla).
152  * @param[out] Bz the vertical (up) component of the magnetic field
153  * (nanotesla).
154  **********************************************************************/
155  void operator()(real t, real lat, real lon, real h,
156  real& Bx, real& By, real& Bz) const {
157  real dummy;
158  Field(t, lat, lon, h, false, Bx, By, Bz, dummy, dummy, dummy);
159  }
160 
161  /**
162  * Evaluate the components of the geomagnetic field and their time
163  * derivatives
164  *
165  * @param[in] t the time (years).
166  * @param[in] lat latitude of the point (degrees).
167  * @param[in] lon longitude of the point (degrees).
168  * @param[in] h the height of the point above the ellipsoid (meters).
169  * @param[out] Bx the easterly component of the magnetic field (nanotesla).
170  * @param[out] By the northerly component of the magnetic field
171  * (nanotesla).
172  * @param[out] Bz the vertical (up) component of the magnetic field
173  * (nanotesla).
174  * @param[out] Bxt the rate of change of \e Bx (nT/yr).
175  * @param[out] Byt the rate of change of \e By (nT/yr).
176  * @param[out] Bzt the rate of change of \e Bz (nT/yr).
177  **********************************************************************/
178  void operator()(real t, real lat, real lon, real h,
179  real& Bx, real& By, real& Bz,
180  real& Bxt, real& Byt, real& Bzt) const {
181  Field(t, lat, lon, h, true, Bx, By, Bz, Bxt, Byt, Bzt);
182  }
183 
184  /**
185  * Create a MagneticCircle object to allow the geomagnetic field at many
186  * points with constant \e lat, \e h, and \e t and varying \e lon to be
187  * computed efficiently.
188  *
189  * @param[in] t the time (years).
190  * @param[in] lat latitude of the point (degrees).
191  * @param[in] h the height of the point above the ellipsoid (meters).
192  * @exception std::bad_alloc if the memory necessary for creating a
193  * MagneticCircle can't be allocated.
194  * @return a MagneticCircle object whose MagneticCircle::operator()(real
195  * lon) member function computes the field at particular values of \e
196  * lon.
197  *
198  * If the field at several points on a circle of latitude need to be
199  * calculated then creating a MagneticCircle and using its member functions
200  * will be substantially faster, especially for high-degree models.
201  **********************************************************************/
202  MagneticCircle Circle(real t, real lat, real h) const;
203 
204  /**
205  * Compute various quantities dependent on the magnetic field.
206  *
207  * @param[in] Bx the \e x (easterly) component of the magnetic field (nT).
208  * @param[in] By the \e y (northerly) component of the magnetic field (nT).
209  * @param[in] Bz the \e z (vertical, up positive) component of the magnetic
210  * field (nT).
211  * @param[out] H the horizontal magnetic field (nT).
212  * @param[out] F the total magnetic field (nT).
213  * @param[out] D the declination of the field (degrees east of north).
214  * @param[out] I the inclination of the field (degrees down from
215  * horizontal).
216  **********************************************************************/
217  static void FieldComponents(real Bx, real By, real Bz,
218  real& H, real& F, real& D, real& I) {
219  real Ht, Ft, Dt, It;
220  FieldComponents(Bx, By, Bz, real(0), real(1), real(0),
221  H, F, D, I, Ht, Ft, Dt, It);
222  }
223 
224  /**
225  * Compute various quantities dependent on the magnetic field and its rate
226  * of change.
227  *
228  * @param[in] Bx the \e x (easterly) component of the magnetic field (nT).
229  * @param[in] By the \e y (northerly) component of the magnetic field (nT).
230  * @param[in] Bz the \e z (vertical, up positive) component of the magnetic
231  * field (nT).
232  * @param[in] Bxt the rate of change of \e Bx (nT/yr).
233  * @param[in] Byt the rate of change of \e By (nT/yr).
234  * @param[in] Bzt the rate of change of \e Bz (nT/yr).
235  * @param[out] H the horizontal magnetic field (nT).
236  * @param[out] F the total magnetic field (nT).
237  * @param[out] D the declination of the field (degrees east of north).
238  * @param[out] I the inclination of the field (degrees down from
239  * horizontal).
240  * @param[out] Ht the rate of change of \e H (nT/yr).
241  * @param[out] Ft the rate of change of \e F (nT/yr).
242  * @param[out] Dt the rate of change of \e D (degrees/yr).
243  * @param[out] It the rate of change of \e I (degrees/yr).
244  **********************************************************************/
245  static void FieldComponents(real Bx, real By, real Bz,
246  real Bxt, real Byt, real Bzt,
247  real& H, real& F, real& D, real& I,
248  real& Ht, real& Ft, real& Dt, real& It);
249  ///@}
250 
251  /** \name Inspector functions
252  **********************************************************************/
253  ///@{
254  /**
255  * @return the description of the magnetic model, if available, from the
256  * Description file in the data file; if absent, return "NONE".
257  **********************************************************************/
258  const std::string& Description() const { return _description; }
259 
260  /**
261  * @return date of the model, if available, from the ReleaseDate field in
262  * the data file; if absent, return "UNKNOWN".
263  **********************************************************************/
264  const std::string& DateTime() const { return _date; }
265 
266  /**
267  * @return full file name used to load the magnetic model.
268  **********************************************************************/
269  const std::string& MagneticFile() const { return _filename; }
270 
271  /**
272  * @return "name" used to load the magnetic model (from the first argument
273  * of the constructor, but this may be overridden by the model file).
274  **********************************************************************/
275  const std::string& MagneticModelName() const { return _name; }
276 
277  /**
278  * @return directory used to load the magnetic model.
279  **********************************************************************/
280  const std::string& MagneticModelDirectory() const { return _dir; }
281 
282  /**
283  * @return the minimum height above the ellipsoid (in meters) for which
284  * this MagneticModel should be used.
285  *
286  * Because the model will typically provide useful results
287  * slightly outside the range of allowed heights, no check of \e t
288  * argument is made by MagneticModel::operator()() or
289  * MagneticModel::Circle.
290  **********************************************************************/
291  Math::real MinHeight() const { return _hmin; }
292 
293  /**
294  * @return the maximum height above the ellipsoid (in meters) for which
295  * this MagneticModel should be used.
296  *
297  * Because the model will typically provide useful results
298  * slightly outside the range of allowed heights, no check of \e t
299  * argument is made by MagneticModel::operator()() or
300  * MagneticModel::Circle.
301  **********************************************************************/
302  Math::real MaxHeight() const { return _hmax; }
303 
304  /**
305  * @return the minimum time (in years) for which this MagneticModel should
306  * be used.
307  *
308  * Because the model will typically provide useful results
309  * slightly outside the range of allowed times, no check of \e t
310  * argument is made by MagneticModel::operator()() or
311  * MagneticModel::Circle.
312  **********************************************************************/
313  Math::real MinTime() const { return _tmin; }
314 
315  /**
316  * @return the maximum time (in years) for which this MagneticModel should
317  * be used.
318  *
319  * Because the model will typically provide useful results
320  * slightly outside the range of allowed times, no check of \e t
321  * argument is made by MagneticModel::operator()() or
322  * MagneticModel::Circle.
323  **********************************************************************/
324  Math::real MaxTime() const { return _tmax; }
325 
326  /**
327  * @return \e a the equatorial radius of the ellipsoid (meters). This is
328  * the value of \e a inherited from the Geocentric object used in the
329  * constructor.
330  **********************************************************************/
331  Math::real EquatorialRadius() const { return _earth.EquatorialRadius(); }
332 
333  /**
334  * @return \e f the flattening of the ellipsoid. This is the value
335  * inherited from the Geocentric object used in the constructor.
336  **********************************************************************/
337  Math::real Flattening() const { return _earth.Flattening(); }
338 
339  /**
340  * @return \e Nmax the maximum degree of the components of the model.
341  **********************************************************************/
342  int Degree() const { return _nmx; }
343 
344  /**
345  * @return \e Mmax the maximum order of the components of the model.
346  **********************************************************************/
347  int Order() const { return _mmx; }
348 
349  /**
350  * \deprecated An old name for EquatorialRadius().
351  **********************************************************************/
352  // GEOGRAPHICLIB_DEPRECATED("Use EquatorialRadius()")
353  Math::real MajorRadius() const { return EquatorialRadius(); }
354  ///@}
355 
356  /**
357  * @return the default path for magnetic model data files.
358  *
359  * This is the value of the environment variable
360  * GEOGRAPHICLIB_MAGNETIC_PATH, if set; otherwise, it is
361  * $GEOGRAPHICLIB_DATA/magnetic if the environment variable
362  * GEOGRAPHICLIB_DATA is set; otherwise, it is a compile-time default
363  * (/usr/local/share/GeographicLib/magnetic on non-Windows systems and
364  * C:/ProgramData/GeographicLib/magnetic on Windows systems).
365  **********************************************************************/
366  static std::string DefaultMagneticPath();
367 
368  /**
369  * @return the default name for the magnetic model.
370  *
371  * This is the value of the environment variable
372  * GEOGRAPHICLIB_MAGNETIC_NAME, if set; otherwise, it is "wmm2020". The
373  * MagneticModel class does not use this function; it is just provided as a
374  * convenience for a calling program when constructing a MagneticModel
375  * object.
376  **********************************************************************/
377  static std::string DefaultMagneticName();
378  };
379 
380 } // namespace GeographicLib
381 
382 #if defined(_MSC_VER)
383 # pragma warning (pop)
384 #endif
385 
386 #endif // GEOGRAPHICLIB_MAGNETICMODEL_HPP
GeographicLib::MagneticModel::DateTime
const std::string & DateTime() const
Definition: MagneticModel.hpp:264
real
GeographicLib::Math::real real
Definition: GeodSolve.cpp:31
GeographicLib::MagneticModel::MinTime
Math::real MinTime() const
Definition: MagneticModel.hpp:313
GeographicLib
Namespace for GeographicLib.
Definition: Accumulator.cpp:12
GeographicLib::MagneticModel::operator()
void operator()(real t, real lat, real lon, real h, real &Bx, real &By, real &Bz, real &Bxt, real &Byt, real &Bzt) const
Definition: MagneticModel.hpp:178
GeographicLib::MagneticModel
Model of the earth's magnetic field.
Definition: MagneticModel.hpp:75
GeographicLib::MagneticModel::MaxTime
Math::real MaxTime() const
Definition: MagneticModel.hpp:324
GeographicLib::MagneticModel::MagneticModelDirectory
const std::string & MagneticModelDirectory() const
Definition: MagneticModel.hpp:280
GEOGRAPHICLIB_EXPORT
#define GEOGRAPHICLIB_EXPORT
Definition: Constants.hpp:92
GeographicLib::MagneticModel::Description
const std::string & Description() const
Definition: MagneticModel.hpp:258
GeographicLib::SphericalHarmonic::normalization
normalization
Definition: SphericalHarmonic.hpp:74
GeographicLib::Geocentric
Geocentric coordinates
Definition: Geocentric.hpp:67
GeographicLib::Math::real
double real
Definition: Math.hpp:121
Geocentric.hpp
Header for GeographicLib::Geocentric class.
GeographicLib::MagneticModel::MagneticModelName
const std::string & MagneticModelName() const
Definition: MagneticModel.hpp:275
GeographicLib::MagneticModel::Flattening
Math::real Flattening() const
Definition: MagneticModel.hpp:337
GeographicLib::MagneticModel::MaxHeight
Math::real MaxHeight() const
Definition: MagneticModel.hpp:302
GeographicLib::Geocentric::EquatorialRadius
Math::real EquatorialRadius() const
Definition: Geocentric.hpp:248
GeographicLib::MagneticModel::FieldComponents
static void FieldComponents(real Bx, real By, real Bz, real &H, real &F, real &D, real &I)
Definition: MagneticModel.hpp:217
Constants.hpp
Header for GeographicLib::Constants class.
GeographicLib::MagneticModel::operator()
void operator()(real t, real lat, real lon, real h, real &Bx, real &By, real &Bz) const
Definition: MagneticModel.hpp:155
GeographicLib::MagneticModel::MagneticFile
const std::string & MagneticFile() const
Definition: MagneticModel.hpp:269
GeographicLib::MagneticModel::Order
int Order() const
Definition: MagneticModel.hpp:347
GeographicLib::Geocentric::Flattening
Math::real Flattening() const
Definition: Geocentric.hpp:255
GeographicLib::MagneticModel::EquatorialRadius
Math::real EquatorialRadius() const
Definition: MagneticModel.hpp:331
GeographicLib::MagneticModel::Degree
int Degree() const
Definition: MagneticModel.hpp:342
GeographicLib::MagneticModel::MinHeight
Math::real MinHeight() const
Definition: MagneticModel.hpp:291
GeographicLib::MagneticModel::MajorRadius
Math::real MajorRadius() const
Definition: MagneticModel.hpp:353
GeographicLib::Geocentric::WGS84
static const Geocentric & WGS84()
Definition: Geocentric.cpp:31
GeographicLib::MagneticCircle
Geomagnetic field on a circle of latitude.
Definition: MagneticCircle.hpp:37
SphericalHarmonic.hpp
Header for GeographicLib::SphericalHarmonic class.