Merge pull request #3 from wkerzendorf/plasma/ab_zeta_fix

Plasma/ab zeta fix

tardis/__init__.py

@@ -9,9 +9,11 @@
 
 from tardis.base import run_tardis
 
+logging.captureWarnings(True)
 logger = logging.getLogger('tardis')
 logger.setLevel(logging.INFO)
 console_handler = logging.StreamHandler()
 console_formatter = logging.Formatter('%(name)s - %(levelname)s - %(message)s')
 console_handler.setFormatter(console_formatter)
-logger.addHandler(console_handler)
+logger.addHandler(console_handler)
+logging.getLogger('py.warnings').addHandler(console_handler)

tardis/plasma/properties/atomic.py

@@ -153,8 +153,6 @@ class ZetaData(BaseAtomicDataProperty):
     outputs = ('zeta_data',)
 
     def _filter_atomic_property(self, zeta_data, selected_atoms):
-        for temperature in range(42000, 82000, 2000):
-            zeta_data[temperature] = zeta_data[40000]
         zeta_data['atomic_number'] = zeta_data.index.labels[0] + 1
         zeta_data['ion_number'] = zeta_data.index.labels[1] + 1
         zeta_data = zeta_data[zeta_data.atomic_number.isin(selected_atoms)]

tardis/plasma/properties/ion_population.py

@@ -1,4 +1,5 @@
 import logging
+import warnings
 
 import numpy as np
 import pandas as pd
@@ -99,16 +100,20 @@ def calculate(t_rad, w, zeta_data, t_electrons, delta,
 
     @staticmethod
     def get_zeta_values(zeta_data, general_phi, t_rad):
-        try:
-            zeta = interpolate.interp1d(zeta_data.columns.values, zeta_data.ix[
-                general_phi.index].values)(t_rad)
-            zeta = zeta.astype(float)
-        except ValueError:
-            raise ValueError('t_rads outside of zeta factor interpolation'
-                             ' zeta_min={0:.2f} zeta_max={1:.2f} '
-                             '- requested {2}'.format(
+        zeta_t_rad = zeta_data.columns.values.astype(np.float64)
+        zeta_values = zeta_data.ix[general_phi.index].values.astype(np.float64)
+        zeta = interpolate.interp1d(zeta_t_rad, zeta_values, bounds_error=False,
+                                    fill_value=np.nan)(t_rad)
+        zeta = zeta.astype(float)
+
+        if np.any(np.isnan(zeta)):
+            warnings.warn('t_rads outside of zeta factor interpolation'
+                                 ' zeta_min={0:.2f} zeta_max={1:.2f} '
+                                 '- replacing with 1s'.format(
                 zeta_data.columns.values.min(), zeta_data.columns.values.max(),
                 t_rad))
+            zeta[np.isnan(zeta)] = 1.0
+
         return zeta
 
 class RadiationFieldCorrection(ProcessingPlasmaProperty):

tardis/tests/test_plasma_simple.py

@@ -4,6 +4,7 @@
 import tardis
 from tardis import atomic
 import pytest
+import warnings
 from tardis.plasma.standard_plasmas import LegacyPlasmaArray
 from tardis.io.util import parse_abundance_dict_to_dataframe
 #from numpy.testing import assert_allclose
@@ -29,9 +30,8 @@ def setup(self):
 #            ionization_mode='nebular', excitation_mode='dilute-lte')
 
     def test_high_temperature(self):
-        with pytest.raises(ValueError) as excinfo:
+        with warnings.catch_warnings(record=True) as w:
             self.plasma.update_radiationfield(t_rad=[100000.],
                 ws=[0.5], j_blues=None, nlte_config=None)
-
-        assert str(excinfo.value).startswith('t_rads outside of zeta '
-                                                'factor interpolation')
+        assert str(w[0].message).startswith('t_rads outside of zeta factor'
+                                               ' interpolation')