GeographicLib
1.21
|
00001 /** 00002 * \file MagneticModel.hpp 00003 * \brief Header for GeographicLib::MagneticModel class 00004 * 00005 * Copyright (c) Charles Karney (2011) <charles@karney.com> and licensed under 00006 * the MIT/X11 License. For more information, see 00007 * http://geographiclib.sourceforge.net/ 00008 **********************************************************************/ 00009 00010 #if !defined(GEOGRAPHICLIB_MAGNETICMODEL_HPP) 00011 #define GEOGRAPHICLIB_MAGNETICMODEL_HPP \ 00012 "$Id: 7f8c59ee3cdbfce252d1172c1bb4d7db7cf5ef38 $" 00013 00014 #include <string> 00015 #include <sstream> 00016 #include <vector> 00017 #include <GeographicLib/Constants.hpp> 00018 #include <GeographicLib/Geocentric.hpp> 00019 #include <GeographicLib/SphericalHarmonic.hpp> 00020 00021 #if defined(_MSC_VER) 00022 // Squelch warnings about dll vs vector 00023 #pragma warning (push) 00024 #pragma warning (disable: 4251) 00025 #endif 00026 00027 namespace GeographicLib { 00028 00029 class MagneticCircle; 00030 00031 /** 00032 * \brief Model of the earth's magnetic field 00033 * 00034 * Evaluate the earth's magnetic field according to a model. At present only 00035 * internal magnetic fields are handled. These are due to the earth's code 00036 * and crust; these vary slowly (over many years). Excluded are the effects 00037 * of currents in the ionosphere and magnetosphere which have daily and 00038 * annual variations. 00039 * 00040 * See \ref magnetic for details of how to install the magnetic model and the 00041 * data format. 00042 * 00043 * See 00044 * - General information: 00045 * - http://geomag.org/models/index.html 00046 * - WMM2010: 00047 * - http://ngdc.noaa.gov/geomag/WMM/DoDWMM.shtml 00048 * - http://ngdc.noaa.gov/geomag/WMM/data/WMM2010/WMM2010COF.zip 00049 * - IGRF11: 00050 * - http://ngdc.noaa.gov/IAGA/vmod/igrf.html 00051 * - http://ngdc.noaa.gov/IAGA/vmod/igrf11coeffs.txt 00052 * - http://ngdc.noaa.gov/IAGA/vmod/geomag70_linux.tar.gz 00053 * - EMM2010: 00054 * - http://ngdc.noaa.gov/geomag/EMM/index.html 00055 * - http://ngdc.noaa.gov/geomag/EMM/data/geomag/EMM2010_Sph_Windows_Linux.zip 00056 * 00057 * Example of use: 00058 * \include example-MagneticModel.cpp 00059 * 00060 * <a href="MagneticField.1.html">MagneticField</a> is a command-line utility 00061 * providing access to the functionality of MagneticModel and MagneticCircle. 00062 **********************************************************************/ 00063 00064 class GEOGRAPHIC_EXPORT MagneticModel { 00065 private: 00066 typedef Math::real real; 00067 static const int idlength_ = 8; 00068 std::string _name, _dir, _description, _date, _filename, _id; 00069 real _t0, _dt0, _tmin, _tmax, _a, _hmin, _hmax; 00070 int _Nmodels; 00071 SphericalHarmonic::normalization _norm; 00072 Geocentric _earth; 00073 std::vector< std::vector<real> > _G; 00074 std::vector< std::vector<real> > _H; 00075 std::vector<SphericalHarmonic> _harm; 00076 void Field(real t, real lat, real lon, real h, bool diffp, 00077 real& Bx, real& By, real& Bz, 00078 real& Bxt, real& Byt, real& Bzt) const throw(); 00079 void ReadMetadata(const std::string& name); 00080 MagneticModel(const MagneticModel&); // copy constructor not allowed 00081 MagneticModel& operator=(const MagneticModel&); // nor copy assignment 00082 public: 00083 00084 /** \name Setting up the magnetic model 00085 **********************************************************************/ 00086 ///@{ 00087 /** 00088 * Construct a magnetic model. 00089 * 00090 * @param[in] name the name of the model. 00091 * @param[in] path (optional) directory for data file. 00092 * @param[in] earth (optional) Geocentric object for converting 00093 * coordinates; default Geocentric::WGS84. 00094 * 00095 * A filename is formed by appending ".wmm" (World Magnetic Model) to the 00096 * name. If \e path is specified (and is non-empty), then the file is 00097 * loaded from directory, \e path. Otherwise the path is given by the 00098 * DefaultMagneticPath(). This may throw an exception because the file 00099 * does not exist, is unreadable, or is corrupt. 00100 * 00101 * This file contains the metadata which specifies the properties of the 00102 * model. The coefficients for the spherical harmonic sums are obtained 00103 * from a file obtained by appending ".cof" to metadata file (so the 00104 * filename ends in ".wwm.cof"). 00105 * 00106 * The model is not tied to a particular ellipsoidal model of the earth. 00107 * The final earth argument to the constructor specify an ellipsoid to 00108 * allow geodetic coordinates to the transformed into the spherical 00109 * coordinates used in the spherical harmonic sum. 00110 **********************************************************************/ 00111 explicit MagneticModel(const std::string& name, 00112 const std::string& path = "", 00113 const Geocentric& earth = Geocentric::WGS84); 00114 ///@} 00115 00116 /** \name Compute the magnetic field 00117 **********************************************************************/ 00118 ///@{ 00119 /** 00120 * Evaluate the components of the geomagnetic field. 00121 * 00122 * @param[in] t the time (years). 00123 * @param[in] lat latitude of the point (degrees). 00124 * @param[in] lon longitude of the point (degrees). 00125 * @param[in] h the height of the point above the ellipsoid (meters). 00126 * @param[out] Bx the easterly component of the magnetic field (nanotesla). 00127 * @param[out] By the northerly component of the magnetic field (nanotesla). 00128 * @param[out] Bz the vertical (up) component of the magnetic field 00129 * (nanotesla). 00130 **********************************************************************/ 00131 void operator()(real t, real lat, real lon, real h, 00132 real& Bx, real& By, real& Bz) const throw() { 00133 real dummy; 00134 Field(t, lat, lon, h, false, Bx, By, Bz, dummy, dummy, dummy); 00135 } 00136 00137 /** 00138 * Evaluate the components of the geomagnetic field and their time 00139 * derivatives 00140 * 00141 * @param[in] t the time (years). 00142 * @param[in] lat latitude of the point (degrees). 00143 * @param[in] lon longitude of the point (degrees). 00144 * @param[in] h the height of the point above the ellipsoid (meters). 00145 * @param[out] Bx the easterly component of the magnetic field (nanotesla). 00146 * @param[out] By the northerly component of the magnetic field (nanotesla). 00147 * @param[out] Bz the vertical (up) component of the magnetic field 00148 * (nanotesla). 00149 * @param[out] Bxt the rate of change of \e Bx (nT/yr). 00150 * @param[out] Byt the rate of change of \e By (nT/yr). 00151 * @param[out] Bzt the rate of change of \e Bz (nT/yr). 00152 **********************************************************************/ 00153 void operator()(real t, real lat, real lon, real h, 00154 real& Bx, real& By, real& Bz, 00155 real& Bxt, real& Byt, real& Bzt) const throw() { 00156 Field(t, lat, lon, h, true, Bx, By, Bz, Bxt, Byt, Bzt); 00157 } 00158 00159 /** 00160 * Create a MagneticCircle object to allow the geomagnetic field at many 00161 * points with constant \e lat, \e h, and \e t and varying \e lon to be 00162 * computed efficiently. 00163 * 00164 * @param[in] t the time (years). 00165 * @param[in] lat latitude of the point (degrees). 00166 * @param[in] h the height of the point above the ellipsoid (meters). 00167 * @return a MagneticCircle object whose MagneticCircle::operator()(real 00168 * lon) member function computes the field at particular values of \e 00169 * lon. 00170 * 00171 * If the field at several points on a circle of latitude need to be 00172 * calculated then creating a MagneticCircle and using its member functions 00173 * will be substantially faster, especially for high-degree models. 00174 **********************************************************************/ 00175 MagneticCircle Circle(real t, real lat, real h) const; 00176 00177 /** 00178 * Compute various quantities dependent on the magnetic field. 00179 * 00180 * @param[in] Bx the \e x (easterly) component of the magnetic field (nT). 00181 * @param[in] By the \e y (northerly) component of the magnetic field (nT). 00182 * @param[in] Bz the \e z (vertical, up positive) component of the magnetic 00183 * field (nT). 00184 * @param[out] H the horizontal magnetic field (nT). 00185 * @param[out] F the total magnetic field (nT). 00186 * @param[out] D the declination of the field (degrees east of north). 00187 * @param[out] I the inclination of the field (degrees down from 00188 * horizontal). 00189 **********************************************************************/ 00190 static void FieldComponents(real Bx, real By, real Bz, 00191 real& H, real& F, real& D, real& I) throw() { 00192 real Ht, Ft, Dt, It; 00193 FieldComponents(Bx, By, Bz, real(0), real(1), real(0), 00194 H, F, D, I, Ht, Ft, Dt, It); 00195 } 00196 00197 /** 00198 * Compute various quantities dependent on the magnetic field and its rate 00199 * of change. 00200 * 00201 * @param[in] Bx the \e x (easterly) component of the magnetic field (nT). 00202 * @param[in] By the \e y (northerly) component of the magnetic field (nT). 00203 * @param[in] Bz the \e z (vertical, up positive) component of the magnetic 00204 * field (nT). 00205 * @param[in] Bxt the rate of change of \e Bx (nT/yr). 00206 * @param[in] Byt the rate of change of \e By (nT/yr). 00207 * @param[in] Bzt the rate of change of \e Bz (nT/yr). 00208 * @param[out] H the horizontal magnetic field (nT). 00209 * @param[out] F the total magnetic field (nT). 00210 * @param[out] D the declination of the field (degrees east of north). 00211 * @param[out] I the inclination of the field (degrees down from 00212 * horizontal). 00213 * @param[out] Ht the rate of change of \e H (nT/yr). 00214 * @param[out] Ft the rate of change of \e F (nT/yr). 00215 * @param[out] Dt the rate of change of \e D (degrees/yr). 00216 * @param[out] It the rate of change of \e I (degrees/yr). 00217 **********************************************************************/ 00218 static void FieldComponents(real Bx, real By, real Bz, 00219 real Bxt, real Byt, real Bzt, 00220 real& H, real& F, real& D, real& I, 00221 real& Ht, real& Ft, real& Dt, real& It) throw(); 00222 ///@} 00223 00224 /** \name Inspector functions 00225 **********************************************************************/ 00226 ///@{ 00227 /** 00228 * @return the description of the magnetic model, if available, from the 00229 * Description file in the data file; if absent, return "NONE". 00230 **********************************************************************/ 00231 const std::string& Description() const throw() { return _description; } 00232 00233 /** 00234 * @return date of the model, if available, from the ReleaseDate field in 00235 * the data file; if absent, return "UNKNOWN". 00236 **********************************************************************/ 00237 const std::string& DateTime() const throw() { return _date; } 00238 00239 /** 00240 * @return full file name used to load the magnetic model. 00241 **********************************************************************/ 00242 const std::string& MagneticFile() const throw() { return _filename; } 00243 00244 /** 00245 * @return "name" used to load the magnetic model (from the first argument 00246 * of the constructor, but this may be overridden by the model file). 00247 **********************************************************************/ 00248 const std::string& MagneticModelName() const throw() { return _name; } 00249 00250 /** 00251 * @return directory used to load the magnetic model. 00252 **********************************************************************/ 00253 const std::string& MagneticModelDirectory() const throw() { return _dir; } 00254 00255 /** 00256 * @return the minimum height above the ellipsoid (in meters) for which 00257 * this MagneticModel should be used. 00258 * 00259 * Because the model will typically provide useful results 00260 * slightly outside the range of allowed heights, no check of \e t 00261 * argument is made by MagneticModel::operator()() or 00262 * MagneticModel::Circle. 00263 **********************************************************************/ 00264 Math::real MinHeight() const throw() { return _hmin; } 00265 00266 /** 00267 * @return the maximum height above the ellipsoid (in meters) for which 00268 * this MagneticModel should be used. 00269 * 00270 * Because the model will typically provide useful results 00271 * slightly outside the range of allowed heights, no check of \e t 00272 * argument is made by MagneticModel::operator()() or 00273 * MagneticModel::Circle. 00274 **********************************************************************/ 00275 Math::real MaxHeight() const throw() { return _hmax; } 00276 00277 /** 00278 * @return the minimum time (in years) for which this MagneticModel should 00279 * be used. 00280 * 00281 * Because the model will typically provide useful results 00282 * slightly outside the range of allowed times, no check of \e t 00283 * argument is made by MagneticModel::operator()() or 00284 * MagneticModel::Circle. 00285 **********************************************************************/ 00286 Math::real MinTime() const throw() { return _tmin; } 00287 00288 /** 00289 * @return the maximum time (in years) for which this MagneticModel should 00290 * be used. 00291 * 00292 * Because the model will typically provide useful results 00293 * slightly outside the range of allowed times, no check of \e t 00294 * argument is made by MagneticModel::operator()() or 00295 * MagneticModel::Circle. 00296 **********************************************************************/ 00297 Math::real MaxTime() const throw() { return _tmax; } 00298 00299 /** 00300 * @return \e a the equatorial radius of the ellipsoid (meters). This is 00301 * the value of \e a inherited from the Geocentric object used in the 00302 * constructor. 00303 **********************************************************************/ 00304 Math::real MajorRadius() const throw() { return _earth.MajorRadius(); } 00305 00306 /** 00307 * @return \e f the flattening of the ellipsoid. This is the value 00308 * inherited from the Geocentric object used in the constructor. 00309 **********************************************************************/ 00310 Math::real Flattening() const throw() { return _earth.Flattening(); } 00311 ///@} 00312 00313 /** 00314 * @return the default path for magnetic model data files. 00315 * 00316 * This is the value of the environment variable MAGNETIC_PATH, if set; 00317 * otherwise, it is $GEOGRAPHICLIB_DATA/magnetic if the environment 00318 * variable GEOGRAPHICLIB_DATA is set; otherwise, it is a compile-time 00319 * default (/usr/local/share/GeographicLib/magnetic on non-Windows systems 00320 * and C:/Documents and Settings/All Users/Application 00321 * Data/GeographicLib/magnetic on Windows systems). 00322 **********************************************************************/ 00323 static std::string DefaultMagneticPath(); 00324 00325 /** 00326 * @return the default name for the magnetic model. 00327 * 00328 * This is the value of the environment variable MAGNETIC_NAME, if set, 00329 * otherwise, it is "wmm2010". The MagneticModel class does not use this 00330 * function; it is just provided as a convenience for a calling program 00331 * when constructing a MagneticModel object. 00332 **********************************************************************/ 00333 static std::string DefaultMagneticName(); 00334 }; 00335 00336 } // namespace GeographicLib 00337 00338 #if defined(_MSC_VER) 00339 #pragma warning (pop) 00340 #endif 00341 00342 #endif // GEOGRAPHICLIB_MAGNETICMODEL_HPP