Hotfix: IGRF-13

eoxmagmod/eoxmagmod/__init__.py

@@ -140,7 +140,7 @@
 except ImportError:
     pass
 
-__version__ = '0.9.6'
+__version__ = '0.9.7'
 __author__ = 'Martin Paces ([email protected])'
 __copyright__ = 'Copyright (C) 2014 EOX IT Services GmbH'
 __licence__ = 'EOX licence (MIT style)'

eoxmagmod/eoxmagmod/data/__init__.py

@@ -44,6 +44,7 @@
 CHAOS_CORE_LATEST = CHAOS7_CORE_LATEST
 IGRF11 = join(_DIRNAME, 'igrf11coeffs.txt')
 IGRF12 = join(_DIRNAME, 'IGRF12.shc')
+IGRF13 = join(_DIRNAME, 'IGRF13.shc')
 SIFM = join(_DIRNAME, 'SIFM.shc')
 LCS1 = join(_DIRNAME, 'LCS-1.shc')
 MF7 = join(_DIRNAME, 'MF7.shc')

eoxmagmod/eoxmagmod/magnetic_model/loader_igrf.py

@@ -35,8 +35,11 @@
 def load_model_igrf(path):
     """ Load model from an IGRF coefficient file.
 
-    Note: IGRF file format was discontinued after IGRF11. Starting by IGRF12
-    the IRFF models uses SHC format which is loaded by the SHC file loader.
+    This loader can be used load to load the models in the original IGRF
+    coefficient file format.
+
+    Starting by IGRF12, the eoxmagmod package contains IGRF models
+    converted to the common SHC format which are loaded by the SHC file loader.
     """
     return SphericalHarmomicGeomagneticModel(load_coeff_igrf(path))
 

eoxmagmod/eoxmagmod/magnetic_model/tests/coefficient_loaders.py

@@ -30,13 +30,11 @@
 from unittest import TestCase, main
 from numpy import inf
 from numpy.testing import assert_allclose
-from eoxmagmod.time_util import (
-    decimal_year_to_mjd2000, decimal_year_to_mjd2000_simple,
-)
+from eoxmagmod.time_util import decimal_year_to_mjd2000
 from eoxmagmod.data import (
     CHAOS6_CORE_LATEST, CHAOS6_STATIC,
     CHAOS7_CORE_LATEST, CHAOS7_STATIC,
-    IGRF11, IGRF12, SIFM, WMM_2015,
+    IGRF11, IGRF12, IGRF13, SIFM, WMM_2015,
     EMM_2010_STATIC, EMM_2010_SECVAR,
     LCS1, MF7,
 )
@@ -186,6 +184,17 @@ class TestCoeffSIFM(TestCase, ShcTestMixIn):
     validity = decimal_year_to_mjd2000((2013.4976, 2015.4962))
 
 
+class TestCoeffIGRF13(TestCase, ShcTestMixIn):
+    class_ = SparseSHCoefficientsTimeDependentDecimalYear
+    path = IGRF13
+    degree = 13
+    min_degree = 1
+    kwargs = {
+        "interpolate_in_decimal_years": True,
+    }
+    validity = decimal_year_to_mjd2000((1900.0, 2025.0))
+
+
 class TestCoeffIGRF12(TestCase, ShcTestMixIn):
     class_ = SparseSHCoefficientsTimeDependentDecimalYear
     path = IGRF12

eoxmagmod/eoxmagmod/magnetic_model/tests/model_loaders.py

@@ -34,9 +34,7 @@
 from numpy.random import uniform
 from numpy.testing import assert_allclose
 from eoxmagmod.magnetic_time import mjd2000_to_magnetic_universal_time
-from eoxmagmod.time_util import (
-    decimal_year_to_mjd2000, decimal_year_to_mjd2000_simple,
-)
+from eoxmagmod.time_util import decimal_year_to_mjd2000
 from eoxmagmod.magnetic_model.loader_shc import (
     load_model_shc, load_model_shc_combined,
 )
@@ -59,7 +57,7 @@
     EMM_2010_STATIC, EMM_2010_SECVAR, WMM_2015,
     CHAOS6_CORE_LATEST, CHAOS6_STATIC,
     CHAOS7_CORE_LATEST, CHAOS7_STATIC,
-    IGRF11, IGRF12, SIFM, LCS1, MF7,
+    IGRF11, IGRF12, IGRF13, SIFM, LCS1, MF7,
 )
 from eoxmagmod.magnetic_model.tests.data import (
     SWARM_MMA_SHA_2C_TEST_DATA,
@@ -475,6 +473,18 @@ def load(self):
         return load_model_shc(IGRF12, interpolate_in_decimal_years=True)
 
 
+class TestIGRF13(TestCase, SHModelTestMixIn):
+    reference_values = (
+        -15522.5, (30.0, 40.0, 8000.0),
+        (15259.57386772841, -159.00767967612023, -14015.952721753336)
+    )
+    degree = 13
+    min_degree = 1
+    validity = decimal_year_to_mjd2000((1900.0, 2025.0))
+
+    def load(self):
+        return load_model_shc(IGRF13, interpolate_in_decimal_years=True)
+
 class TestSIFM(TestCase, SHModelTestMixIn):
     reference_values = (
         5295.36, (30.0, 40.0, 8000.0),

eoxmagmod/eoxmagmod/magnetic_model/tests/parser_shc.py

@@ -34,7 +34,7 @@
 from eoxmagmod.data import (
     CHAOS6_CORE_LATEST, CHAOS6_STATIC,
     CHAOS7_CORE_LATEST, CHAOS7_STATIC,
-    IGRF12, SIFM, LCS1, MF7,
+    IGRF12, IGRF13, SIFM, LCS1, MF7,
 )
 
 
@@ -94,6 +94,16 @@ def test_parse_shc_file_igrf12(self):
             "nstep": 1,
         })
 
+    def test_parse_shc_file_igrf13(self):
+        data = self.parse(IGRF13)
+        self._assert_valid(data, {
+            "degree_min": 1,
+            "degree_max": 13,
+            "spline_order": 2,
+            "ntime": 26,
+            "nstep": 1,
+        })
+
     def test_parse_shc_header_igrf12(self):
         self._test_header(IGRF12)
 

eoxmagmod/eoxmagmod/tests/data/update_quasi_dipole_test_data.py

@@ -0,0 +1,96 @@
+#-------------------------------------------------------------------------------
+#
+#  Quasi-Dipole coordinates - test
+#
+# Author: Martin Paces <[email protected]>
+#
+#-------------------------------------------------------------------------------
+# Copyright (C) 2018 EOX IT Services GmbH
+#
+# Permission is hereby granted, free of charge, to any person obtaining a copy
+# of this software and associated documentation files (the "Software"), to deal
+# in the Software without restriction, including without limitation the rights
+# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+# copies of the Software, and to permit persons to whom the Software is
+# furnished to do so, subject to the following conditions:
+#
+# The above copyright notice and this permission notice shall be included in
+# all copies of this Software or works derived from this Software.
+#
+# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
+# THE SOFTWARE.
+#-------------------------------------------------------------------------------
+#pylint: disable=missing-docstring
+
+import sys
+try:
+    # Python 2
+    from itertools import izip as zip
+except ImportError:
+    pass
+from numpy import array
+from eoxmagmod import mjd2000_to_decimal_year
+from eoxmagmod.quasi_dipole_coordinates import (
+    eval_mlt, eval_subsol, eval_qdlatlon_with_base_vectors,
+)
+
+EARTH_RADIUS = 6371.2 # km
+START_TIME = 0.0    # MJD2000 / 2000-01-01T00:00:00Z
+END_TIME = 7305.0   # MJD2000 / 2020-01-01T00:00:00Z
+
+
+def load_coordinates(file_in):
+
+    _ = next(file_in)
+
+    times, lats, lons, rads = [], [], [], []
+
+    for line in file_in:
+        time, _, lat, lon, rad = [float(v) for v in line.split()[:5]]
+        times.append(time)
+        lats.append(lat)
+        lons.append(lon)
+        rads.append(rad)
+
+    return array(times), array(lats), array(lons), array(rads)
+
+
+def update_test_data(file_out, time_mjd2000, lat, lon, rad):
+    """ Generate test dataset. """
+    time_decimal_year = mjd2000_to_decimal_year(time_mjd2000)
+
+    qdlat, qdlon, f11, f12, f21, f22, f__ = eval_qdlatlon_with_base_vectors(
+        lat, lon, rad, time_decimal_year
+    )
+    mlt = eval_mlt(qdlon, time_mjd2000)
+    sol_lat, sol_lon = eval_subsol(time_mjd2000)
+
+    header = [
+        "MJD2000", "DecimalYear", "Latitude", "Longitude", "Radius",
+        "QDLatitude", "QDLongitude", "F11", "F12", "F21", "F22", "F",
+        "MagneticLocalTime", "SubsolarLatitude", "SubsolarLongitude",
+    ]
+    records = zip(
+        time_mjd2000, time_decimal_year, lat, lon, rad,
+        qdlat, qdlon, f11, f12, f21, f22, f__,
+        mlt, sol_lat, sol_lon
+    )
+
+    file_out.write("\t".join(header) + "\n")
+    for record in records:
+        file_out.write("\t".join("%.14e" % value for value in record) + "\n")
+
+
+def main(filename):
+    with open(filename) as file_in:
+        time, lat, lon, rad = load_coordinates(file_in)
+    update_test_data(sys.stdout, time, lat, lon, rad)
+
+
+if __name__ == "__main__":
+    main(sys.argv[1])

qdipole/Makefile.in

@@ -25,11 +25,14 @@ INC=
 HDR=cqdipole.h
 OBJ=cqdipole.o qdipole.o apex.o apexsh.o makeapexsh.o eval_subsol.o eval_qdlatlon.o eval_mlt.o
 DOC=README
-DATA=apexsh_1980-2020.txt apexsh_1995-2015.txt \
+DATA=apexsh_1980-2025.txt apexsh_1980-2020.txt apexsh_1995-2015.txt \
+ test_result.apexsh_1980-2025.txt \
  test_result.apexsh_1980-2020.txt \
  test_result.apexsh_1995-2015.txt \
+ test_result.qdlatlon.apexsh_1980-2025.txt \
  test_result.qdlatlon.apexsh_1980-2020.txt \
  test_result.qdlatlon.apexsh_1995-2015.txt \
+ test_result.mlt.apexsh_1980-2025.txt \
  test_result.mlt.apexsh_1980-2020.txt \
  test_result.mlt.apexsh_1995-2015.txt \
  test_result.subsol.txt

qdipole/README

@@ -12,6 +12,9 @@ EOX IT Services GmbH, Vienna, Austria
 CHANGE LOG
 ----------
 
+2020-01-17 - 0.3.0 Update
+ - update of the IGRF coefficients
+
 2016-06-09 - 0.2.0 Update
  - review and re-factoring of the C interface
  - separating MLT calculation from the DQ latitude/longitude evaluation
@@ -30,7 +33,7 @@ INSTALLATION
 
 Step 1 - configure:
 
-  $ ./configure --prefix=/usr/ --libdir='${exec_prefix}/lib64'
+  $ ./configure --prefix=/usr/
 
 For more detail see
 
@@ -52,8 +55,11 @@ TESTING
 
   $ ./qdipole_test apexsh_1995-2015.txt | diff - test_result.apexsh_1995-2015.txt
   $ ./qdipole_test apexsh_1980-2020.txt | diff - test_result.apexsh_1980-2020.txt
+  $ ./qdipole_test apexsh_1980-2025.txt | diff - test_result.apexsh_1980-2025.txt
   $ ./qdlatlon_test apexsh_1995-2015.txt | diff - test_result.qdlatlon.apexsh_1995-2015.txt
   $ ./qdlatlon_test apexsh_1980-2020.txt | diff - test_result.qdlatlon.apexsh_1980-2020.txt
+  $ ./qdlatlon_test apexsh_1980-2025.txt | diff - test_result.qdlatlon.apexsh_1980-2025.txt
   $ ./mlt_test apexsh_1995-2015.txt | diff - test_result.mlt.apexsh_1995-2015.txt
   $ ./mlt_test apexsh_1980-2020.txt | diff - test_result.mlt.apexsh_1980-2020.txt
+  $ ./mlt_test apexsh_1980-2025.txt | diff - test_result.mlt.apexsh_1980-2025.txt
   $ ./subsol_test | diff - test_result.subsol.txt