Merge pull request #8 from ESA-VirES/swarm_shc_models

eoxmagmod 0.5.1 update

eoxmagmod/eoxmagmod/__init__.py

@@ -96,6 +96,10 @@
 from .magnetic_model.loader_mma import (
     load_model_swarm_mma_2c_internal,
     load_model_swarm_mma_2c_external,
+    load_model_swarm_mma_2f_geo_internal,
+    load_model_swarm_mma_2f_geo_external,
+    load_model_swarm_mma_2f_sm_internal,
+    load_model_swarm_mma_2f_sm_external,
 )
 from .magnetic_model.loader_mio import (
     load_model_swarm_mio_internal,
@@ -164,12 +168,16 @@
     'load_model_emm',
     'load_model_swarm_mma_2c_internal',
     'load_model_swarm_mma_2c_external',
+    'load_model_swarm_mma_2f_geo_internal',
+    'load_model_swarm_mma_2f_geo_external',
+    'load_model_swarm_mma_2f_sm_internal',
+    'load_model_swarm_mma_2f_sm_external',
     'load_model_swarm_mio_internal',
     'load_model_swarm_mio_external',
     'trace_field_line',
 ]
 
-__version__ = '0.5.0'
+__version__ = '0.5.1'
 __author__ = 'Martin Paces ([email protected])'
 __copyright__ = 'Copyright (C) 2014 EOX IT Services GmbH'
 __licence__ = 'EOX licence (MIT style)'

eoxmagmod/eoxmagmod/magnetic_model/loader_mio.py

@@ -27,7 +27,10 @@
 #-------------------------------------------------------------------------------
 
 from numpy import arange
-from .model_mio import DipoleMIOGeomagneticModel, MIO_EARTH_RADIUS
+from .model_mio import (
+    DipoleMIOPrimaryGeomagneticModel, DipoleMIOGeomagneticModel,
+    MIO_EARTH_RADIUS,
+)
 from .coefficients_mio import SparseSHCoefficientsMIO
 from .parser_mio import parse_swarm_mio_file
 
@@ -39,11 +42,31 @@ def load_model_swarm_mio_internal(path):
     return _create_mio_model(coefficients, params)
 
 
-def load_model_swarm_mio_external(path, above_ionosphere=True):
+def load_model_swarm_mio_external(path, above_ionosphere=None):
     """ Load external (primary field) model from a Swarm MIO_SHA_2* product.
     """
-    coefficients, params = load_coeff_swarm_mio_external(path, above_ionosphere)
-    return _create_mio_model(coefficients, params)
+    with open(path, "rb") as file_in:
+        params = parse_swarm_mio_file(file_in)
+
+    if above_ionosphere is None:
+        return _create_composed_mio_model(
+            _get_coeff_swarm_mio_external(params, False),
+            _get_coeff_swarm_mio_external(params, True),
+            params
+        )
+    else:
+        return _create_mio_model(
+            _get_coeff_swarm_mio_external(params, above_ionosphere), params
+        )
+
+
+def _create_composed_mio_model(coefficients_below_ionosphere,
+                               coefficients_above_ionosphere, params):
+    return DipoleMIOPrimaryGeomagneticModel(
+        _create_mio_model(coefficients_below_ionosphere, params),
+        _create_mio_model(coefficients_above_ionosphere, params),
+        height=params["height"],
+    )
 
 
 def _create_mio_model(coefficients, params):
@@ -74,6 +97,11 @@ def load_coeff_swarm_mio_external(path, above_ionosphere=True):
     with open(path, "rb") as file_in:
         data = parse_swarm_mio_file(file_in)
 
+    return _get_coeff_swarm_mio_external(data, above_ionosphere), data
+
+
+def _get_coeff_swarm_mio_external(data, above_ionosphere):
+    """ Create coefficient object for the given source data. """
     indices = data["nm"]
     coefficients = data["qs"]
     if above_ionosphere:
@@ -88,7 +116,7 @@ def load_coeff_swarm_mio_external(path, above_ionosphere=True):
         indices, coefficients,
         ps_extent=(data["pmin"], data["pmax"], data["smin"], data["smax"]),
         is_internal=is_internal,
-    ), data
+    )
 
 
 def convert_external_mio_coeff(degree, indices, coefficients, height):

eoxmagmod/eoxmagmod/magnetic_model/loader_mma.py

@@ -28,12 +28,17 @@
 # pylint: disable=invalid-name
 
 from spacepy import pycdf
-from .model import DipoleSphericalHarmomicGeomagneticModel
+from .model import (
+    SphericalHarmomicGeomagneticModel,
+    DipoleSphericalHarmomicGeomagneticModel,
+)
 from .coefficients import (
     SparseSHCoefficientsTimeDependent, CombinedSHCoefficients,
 )
 from .parser_mma import (
     read_swarm_mma_2c_internal, read_swarm_mma_2c_external,
+    read_swarm_mma_2f_geo_internal, read_swarm_mma_2f_geo_external,
+    read_swarm_mma_2f_sm_internal, read_swarm_mma_2f_sm_external,
 )
 
 # North geomagnetic coordinates used by the MMA products (IGRF-11, 2010.0)
@@ -61,27 +66,115 @@ def load_model_swarm_mma_2c_external(path, lat_ngp=MMA2C_NGP_LATITUDE,
     )
 
 
+def load_model_swarm_mma_2f_geo_internal(path):
+    """ Load geographic frame internal (secondary field) model from a Swarm
+    MMA_SHA_2F product.
+    """
+    return SphericalHarmomicGeomagneticModel(
+        load_coeff_swarm_mma_2f_geo_internal(path)
+    )
+
+
+def load_model_swarm_mma_2f_geo_external(path):
+    """ Load geographic frame internal (primary field) model from a Swarm
+    MMA_SHA_2F product.
+    """
+    return SphericalHarmomicGeomagneticModel(
+        load_coeff_swarm_mma_2f_geo_external(path)
+    )
+
+
+def load_model_swarm_mma_2f_sm_internal(path):
+    """ Load solar magnetic frame internal (secondary field) model from a Swarm
+    MMA_SHA_2F product.
+    """
+    # FIXME: solar-magnetic frame model
+    return SphericalHarmomicGeomagneticModel(
+        load_coeff_swarm_mma_2f_sm_internal(path),
+    )
+
+
+def load_model_swarm_mma_2f_sm_external(path):
+    """ Load solar magnetic frame internal (primary field) model from a Swarm
+    MMA_SHA_2F product.
+    """
+    # FIXME: solar-magnetic frame model
+    return SphericalHarmomicGeomagneticModel(
+        load_coeff_swarm_mma_2f_sm_external(path)
+    )
+
+
 def load_coeff_swarm_mma_2c_internal(path):
     """ Load internal model coefficients from a Swarm MMA_SHA_2C product file.
     """
-    with pycdf.CDF(path) as cdf:
-        data = read_swarm_mma_2c_internal(cdf)
-
-    return CombinedSHCoefficients(*[
-        SparseSHCoefficientsTimeDependent(
-            item["nm"], item["gh"], item["t"], is_internal=True
-        ) for item in data
-    ])
+    return _load_coeff_swarm_mma_2c(
+        path, read_swarm_mma_2c_internal, "gh", is_internal=True
+    )
 
 
 def load_coeff_swarm_mma_2c_external(path):
+    """ Load external model coefficients from a Swarm MMA_SHA_2C product file.
+    """
+    return _load_coeff_swarm_mma_2c(
+        path, read_swarm_mma_2c_external, "qs", is_internal=False
+    )
+
+
+def load_coeff_swarm_mma_2f_geo_internal(path):
+    """ Load internal geographic frame model coefficients from a Swarm
+    MMA_SHA_2F product file.
+    """
+    return _load_coeff_swarm_mma_2f(
+        path, read_swarm_mma_2f_geo_internal, "gh", is_internal=False
+    )
+
+
+def load_coeff_swarm_mma_2f_geo_external(path):
+    """ Load external geographic frame model coefficients from a Swarm
+    MMA_SHA_2F product file.
+    """
+    return _load_coeff_swarm_mma_2f(
+        path, read_swarm_mma_2f_geo_external, "qs", is_internal=False
+    )
+
+
+def load_coeff_swarm_mma_2f_sm_internal(path):
+    """ Load internal solar magnetic frame model coefficients from a Swarm
+    MMA_SHA_2F product file.
+    """
+    return _load_coeff_swarm_mma_2f(
+        path, read_swarm_mma_2f_sm_internal, "gh", is_internal=False
+    )
+
+
+def load_coeff_swarm_mma_2f_sm_external(path):
+    """ Load external solar magnetic frame model coefficients from a Swarm
+    MMA_SHA_2F product file.
+    """
+    return _load_coeff_swarm_mma_2f(
+        path, read_swarm_mma_2f_sm_external, "qs", is_internal=False
+    )
+
+
+def _load_coeff_swarm_mma_2c(path, cdf_reader, variable, is_internal):
     """ Load external model coefficients from a Swarm MMA_SHA_2C product file.
     """
     with pycdf.CDF(path) as cdf:
-        data = read_swarm_mma_2c_external(cdf)
+        data = cdf_reader(cdf)
 
     return CombinedSHCoefficients(*[
         SparseSHCoefficientsTimeDependent(
-            item["nm"], item["qs"], item["t"], is_internal=False
+            item["nm"], item[variable], item["t"], is_internal=is_internal
         ) for item in data
     ])
+
+
+def _load_coeff_swarm_mma_2f(path, cdf_reader, variable, is_internal):
+    """ Load external model coefficients from a Swarm MMA_SHA_2F product file.
+    """
+    with pycdf.CDF(path) as cdf:
+        data = cdf_reader(cdf)
+
+    return SparseSHCoefficientsTimeDependent(
+        data["nm"], data[variable], data["t"], is_internal=is_internal
+    )

eoxmagmod/eoxmagmod/magnetic_model/model.py

@@ -84,21 +84,22 @@ def _eval_multi_time(self, time, coords, input_coordinate_system,
                          output_coordinate_system, **options):
         """ Evaluate spherical harmonic for multiple times. """
         result = empty(coords.shape)
-        iterator = nditer(
-            [
-                time, coords[..., 0], coords[..., 1], coords[..., 2],
-                result[..., 0], result[..., 1], result[..., 2],
-            ],
-            op_flags=[
-                ['readonly'], ['readonly'], ['readonly'], ['readonly'],
-                ['writeonly'], ['writeonly'], ['writeonly'],
-            ],
-        )
-        for time_, coord0, coord1, coord2, vect0, vect1, vect2 in iterator:
-            vect0[...], vect1[...], vect2[...] = self._eval_single_time(
-                time_, [coord0, coord1, coord2], input_coordinate_system,
-                output_coordinate_system, **options
+        if result.size > 0:
+            iterator = nditer(
+                [
+                    time, coords[..., 0], coords[..., 1], coords[..., 2],
+                    result[..., 0], result[..., 1], result[..., 2],
+                ],
+                op_flags=[
+                    ['readonly'], ['readonly'], ['readonly'], ['readonly'],
+                    ['writeonly'], ['writeonly'], ['writeonly'],
+                ],
             )
+            for time_, coord0, coord1, coord2, vect0, vect1, vect2 in iterator:
+                vect0[...], vect1[...], vect2[...] = self._eval_single_time(
+                    time_, [coord0, coord1, coord2], input_coordinate_system,
+                    output_coordinate_system, **options
+                )
         return result
 
     def _eval_single_time(self, time, coords, input_coordinate_system,