‹‹‹ Blog Overview

International Geomagnetic Reference Field vs. Python
On the proliferation of implementations - plus a new one

There are literally a ton of Python-implementations of the International Geomagnetic Reference Field (IGRF) out there. I just added a new one to the mix. Why?!?

  • IGRF
  • Python
  • geomagnetism
  • geophysics
  • magnetic field

Where to start?

This discussion on Github beautifully summarizes the problem: "Too many IGRFs"

In all seriousness, it is a lot of implementations, at least 10 as of writing this text. This does not need to be a bad thing - competition is usually good and I like when science gets reproduced. However, what I found when I entered this arena was that few of the many implementations had been properly tested or verified. Besides, almost everyone of them had its smaller or bigger issues, anything from structure, design, speed, packaging, documentation etc. pyIGRF by zzyztyy was by far not the best, but it was (a) the one closest to IAGA's original Fortran implementation and, given the previous constraint, (b) the most "pure Python" one. So I based my work on it. I needed something I could trust, fast and verified.

What are the options?

There are other, independent implementations named pyIGRF, not to be confused with the previously mentioned one:

The official IGRF implementations by IAGA are available via the NOAA website:

A good, modern and pure Python implementation:

Another Python implementation is part of the navtools package:

Another Python wrapper around a cleaned-up version of IAGA's Fortran code is maintained as part of the space physics project:

A more or less undocumented, pure Python implementation can be found here:

All of the above models rely on the IGRF-13 coefficients1 as provided by NOAA.

The IGRF is not to be confused with the World Magnetic Model (WMM). Recent implementations of the WMM accessible via Python are maintained as part of the space physics project:

A new implementation

Just to make a small difference, I named my implementation pyCRGI, using the French translation of IGRF's name: Champ de Référence Géomagnétique International.

pyCRGI is a massively cleaned-up and modernized fork of pyIGRF by zzyztyy. Be aware that there are a number of small function and module name differences to the original pyIGRF package. The fork's main goals are verified results, tests, speed and ease of maintainability. This is work in progress. The package offers different implementations, pure Python, Python JIT-compiled via numba and array-oriented via numpy and numba.

The current project status can be summarized as follows:

  • [x] package structure cleanup, allowing proper testing and packaging
  • [x] type annotations
  • [x] doc strings completed and prepared for Sphinx autodoc
  • [x] debug mode via environment variable PYCGIR_DEBUG=1
  • [x] pure Python 3 implementation without dependency to numpy
  • [x] JIT-compiled implementation depending on numba and numpy, installation target jited
  • [x] array implementation depending on numba and numpy, installation target array
  • [x] unit tests via test makefile target
  • [x] verification of field's values against original Fortran implementation
  • [ ] verification of field's annual/seasonal variation -> huge differences to original Fortran implementation
  • [x] benchmark
  • [ ] clean up and add missing relevant coordinate system conversions

How does it compare?

Interactive plot 1: Benchmark of pyCRGI for different quantities of data, different years and different itypes (1==geodetic, 2==geocentric) carried out on an AMD Epyc 7443p

How to Install?

Use pip to install the latest development version from Github:

pip install git+https://github.com/pleiszenburg/pyCRGI.git@master
Source code 1: Installing the latest development version directly from Github

How to Use it?

First import the package, either as pure Python 3 or JIT-compiled via numba and numpy:

from pyCRGI.pure import get_value, get_variation # pure Python 3
# or
from pyCRGI.jited import get_value, get_variation # JIT-compiled via `numba` and `numpy`
# or
from pyCRGI.array import get_value, get_variation # array implementation via `numba` and `numpy`
Source code 2: Importing different implementations for different use-cases

You can calculate the magnetic field's intensity:

get_value(lat, lon, alt, year)
Source code 3: Calculating the field's intensity

You can calculate the annual variation of the magnetic field's intensity:

get_variation(lat, lon, alt, year)
Source code 4: Calculating the variation of the field's intensity

The return value is a tuple (or array) of seven floating point numbers representing the local magnetic field:

  • D: declination (+ve east) [degree]
  • I: inclination (+ve down) [degree]
  • H: horizontal intensity [nT]
  • X: north component [nT]
  • Y: east component [nT]
  • Z: vertical component (+ve down) [nT]
  • F: total intensity [nT]

If you want to use the IGRF-13 model in a more flexible manner, you can use the functions geodetic2geocentric and get_syn. They are somewhat closer to the original Fortran implementation.

Another function, get_coeffs, can be used to get g[m][n] or h[m][n] corresponding to the IGRF's formula.

‹‹‹ Blog Overview