Add IGRF-14 (#17)

* Change tests to distinguish DGRF and IGRF generations * Add IGRF14 coefficients * Make IGRF-14 the default * Include IGRF14.shc in package * Bump version to 2.0.0-alpha
IAGA-VMOD · Nov 25, 2024 · c06878c · c06878c
1 parent 4dc503d
commit c06878c
Show file tree

Hide file tree

Showing 14 changed files with 174,160 additions and 19 deletions.
diff --git a/setup.py b/setup.py
@@ -17,7 +17,7 @@
 
 setuptools.setup(
     name="ppigrf",
-    version="1.1.0",
+    version="2.0.0-alpha",
     author="Karl Laundal",
     author_email="readme@file.md",
     description="Pure Python IGRF",
@@ -34,7 +34,7 @@
         'pandas>=1.3.5'
     ],
     package_dir={"": "src"},
-    package_data={'':['IGRF13.shc']},
+    package_data={'':['IGRF13.shc', 'IGRF14.shc']},
     packages=setuptools.find_packages(where="src"),
     python_requires=">=3.6",
 )
diff --git a/src/ppigrf/IGRF14.shc b/src/ppigrf/IGRF14.shc
diff --git a/src/ppigrf/ppigrf.py b/src/ppigrf/ppigrf.py
@@ -77,7 +77,9 @@
 import os
 
 basepath = os.path.dirname(__file__)
-shc_fn = basepath + '/IGRF13.shc' # Default shc file
+shc_fn = basepath + '/IGRF14.shc' # Default shc file
+shc_fn_igrf13 = basepath + '/IGRF13.shc'
+shc_fn_igrf14 = basepath + '/IGRF14.shc'
 
 # Geomagnetic reference radius:
 RE = 6371.2 # km

diff --git a/test/data/README.md b/test/data/README.md
@@ -15,3 +15,18 @@ following the `YYYY-MM-DD` format.
 
 The `b_z.csv` contains the **downward** components of the magnetic vector on
 each location.
+
+Note that each IGRF release includes the DGRF (definitive field) for the previous epoch. IGRF-14 includes the DGRF for 2020.0, overriding the IGRF for 2020.0 issued within IGRF-13. This means that moving from IGRF-13 to IGRF-14 changes the results for calculations within the 2015-2025 period. Calculations within the definitive period (i.e. for IGRF-14 and lower: up to 2015.0) should remain the same.
+
+```
+Example requests
+
+dgrf-2010-01-01/b_n.csv:
+https://www.ngdc.noaa.gov/geomag-web/calculators/calculateIgrfgrid?browserRequest=true&key=gFE5W&lat1=80&lat1Hemisphere=S&lat2=80&lat2Hemisphere=N&latStepSize=1.0&lon1=179&lon1Hemisphere=W&lon2=180&lon2Hemisphere=E&lonStepSize=1.0&coordinateSystem=D&elevation=5&elevationUnits=K&magneticComponent=x&model=IGRF&startYear=2010&startMonth=01&startDay=1&endYear=2010&endMonth=1&endDay=1&dateStepSize=1.0&resultFormat=csv
+
+dgrf-2010-01-01/b_e.csv:
+https://www.ngdc.noaa.gov/geomag-web/calculators/calculateIgrfgrid?browserRequest=true&key=gFE5W&lat1=80&lat1Hemisphere=S&lat2=80&lat2Hemisphere=N&latStepSize=1.0&lon1=179&lon1Hemisphere=W&lon2=180&lon2Hemisphere=E&lonStepSize=1.0&coordinateSystem=D&elevation=5&elevationUnits=K&magneticComponent=y&model=IGRF&startYear=2010&startMonth=01&startDay=1&endYear=2010&endMonth=1&endDay=1&dateStepSize=1.0&resultFormat=csv
+
+dgrf-2010-01-01/b_z.csv:
+https://www.ngdc.noaa.gov/geomag-web/calculators/calculateIgrfgrid?browserRequest=true&key=gFE5W&lat1=80&lat1Hemisphere=S&lat2=80&lat2Hemisphere=N&latStepSize=1.0&lon1=179&lon1Hemisphere=W&lon2=180&lon2Hemisphere=E&lonStepSize=1.0&coordinateSystem=D&elevation=5&elevationUnits=K&magneticComponent=z&model=IGRF&startYear=2010&startMonth=01&startDay=1&endYear=2010&endMonth=1&endDay=1&dateStepSize=1.0&resultFormat=csv
+```
diff --git a/test/data/dgrf-2010-01-01/b_e.csv b/test/data/dgrf-2010-01-01/b_e.csv
diff --git a/test/data/dgrf-2010-01-01/b_n.csv b/test/data/dgrf-2010-01-01/b_n.csv
diff --git a/test/data/dgrf-2010-01-01/b_z.csv b/test/data/dgrf-2010-01-01/b_z.csv
diff --git a/test/data/2020-01-01/b_e.csv → test/data/igrf13-2020-01-01/b_e.csv b/test/data/2020-01-01/b_e.csv → test/data/igrf13-2020-01-01/b_e.csv
diff --git a/test/data/2020-01-01/b_n.csv → test/data/igrf13-2020-01-01/b_n.csv b/test/data/2020-01-01/b_n.csv → test/data/igrf13-2020-01-01/b_n.csv
diff --git a/test/data/2020-01-01/b_z.csv → test/data/igrf13-2020-01-01/b_z.csv b/test/data/2020-01-01/b_z.csv → test/data/igrf13-2020-01-01/b_z.csv
diff --git a/test/data/2022-10-05/b_e.csv → test/data/igrf13-2022-10-05/b_e.csv b/test/data/2022-10-05/b_e.csv → test/data/igrf13-2022-10-05/b_e.csv
diff --git a/test/data/2022-10-05/b_n.csv → test/data/igrf13-2022-10-05/b_n.csv b/test/data/2022-10-05/b_n.csv → test/data/igrf13-2022-10-05/b_n.csv
diff --git a/test/data/2022-10-05/b_z.csv → test/data/igrf13-2022-10-05/b_z.csv b/test/data/2022-10-05/b_z.csv → test/data/igrf13-2022-10-05/b_z.csv
diff --git a/test/test_igrf.py b/test/test_igrf.py
@@ -10,25 +10,20 @@
 from datetime import datetime
 
 from ppigrf import igrf
-from ppigrf.ppigrf import yearfrac_to_datetime
+from ppigrf.ppigrf import yearfrac_to_datetime, shc_fn, shc_fn_igrf13, shc_fn_igrf14
 
 # Define paths to test directory and test data directory
 TEST_DIR = Path(os.path.dirname(__file__))
 TEST_DATA_DIR = TEST_DIR / "data"
 
 
-def load_precomputed_igrf(date):
+def load_precomputed_igrf(directory):
     """
-    Loads the precomputed IGRF files from a given date
-
-    Available dates:
-        * 2020-01-01
-        * 2022-10-05
+    Loads the precomputed IGRF files from a given directory
 
     Parameters
     ----------
-    date : :class:`datetime.datetime` object
-        Date of the precomputed IGRF field files that will be loaded.
+    directory : :class:`pathlib.Path`
 
     Returns
     -------
@@ -37,13 +32,12 @@ def load_precomputed_igrf(date):
         date.
     """
     # Read the csv files
-    date_dir = TEST_DATA_DIR / date.strftime("%Y-%m-%d")
     first_columns = ["date", "latitude", "longitude", "altitude_km"]
     components = ("b_e", "b_n", "b_z")
     dataframes = []
     for component in components:
         columns = first_columns + [component, component + "_sv"]
-        fname = date_dir / f"{component}.csv"
+        fname = directory / f"{component}.csv"
         df = pd.read_csv(fname, skiprows=13, names=columns)
         dataframes.append(df)
     # Merge the dataframes
@@ -64,11 +58,21 @@ class TestIGRFKnownValues:
     rtol = 1e-2  # 1% of error
 
     @pytest.mark.parametrize(
-        "date, atol",
-        [[datetime(2020, 1, 1), 1], [datetime(2022, 10, 5), 4]],
-        ids=["2020-01-01", "2022-10-05"],
+        "date, subdirectory, igrf_version, atol",
+        [
+            [datetime(2010, 1, 1), "dgrf-2010-01-01", shc_fn_igrf13, 1],
+            [datetime(2010, 1, 1), "dgrf-2010-01-01", shc_fn_igrf14, 1],
+            [datetime(2020, 1, 1), "igrf13-2020-01-01", shc_fn_igrf13, 1],
+            [datetime(2022, 10, 5), "igrf13-2022-10-05", shc_fn_igrf13, 4],
+        ],
+        ids=[
+            "IGRF-13: dgrf-2010-01-01",
+            "IGRF-14: dgrf-2010-01-01",
+            "IGRF-13: igrf13-2020-01-01",
+            "IGRF-13: igrf13-2022-10-05",
+        ],
     )
-    def test_igrf(self, date, atol):
+    def test_igrf(self, date, subdirectory, igrf_version, atol):
         """
         Test IGRF against the precomputed values
 
@@ -79,7 +83,7 @@ def test_igrf(self, date, atol):
         differences due to different types of dates interpolations.
         """
         # Get precomputed IGRF field
-        precomputed_igrf = load_precomputed_igrf(date)
+        precomputed_igrf = load_precomputed_igrf(TEST_DATA_DIR/subdirectory)
         # Overwrite the date with the one in the data file
         # date = precomputed_igrf.date.values[0]
         # Compute igrf using ppigrf
@@ -88,6 +92,7 @@ def test_igrf(self, date, atol):
             precomputed_igrf.latitude,
             precomputed_igrf.altitude_km,
             date,
+            coeff_fn=igrf_version,
         )
         # Ravel the arrays
         b_e, b_n, b_u = tuple(np.ravel(component) for component in (b_e, b_n, b_u))