Fixes #830 (#832)

* Update tardis_model_reader; reads all columns now

* Split CSV and CMGFEN composition parsers

* Update tests to work with new parsers

* Update and add test data files

* Fix some issues about photospheric composition

* Update documentation about CSV composition

* Update test config

* Update scheme to work with cmfgen density file type

* Fix bug in model_reader

docs/examples/abundancecust.rst

@@ -92,6 +92,13 @@ The example file shown here has three simple layers:
     :code:`Ni` and :code:`Ni56`. 
     Here, Ni will refer to the stable Nickel, i.e. (Z=26, A=58).
 
+
+.. note::
+
+    As with the custom density file, the first row will be ignored. It is
+    supposed to give the composition below the photosphere. Thus, the first row
+    (after the header) can be filled with dummy values
+
 .. warning::
 
    The example given here is to show the format only. It is not a
@@ -118,4 +125,4 @@ If you want to convert an ASCII abundance file(say "abund.dat") to CSV format, y
 
     from tardis.util import convert_abundances_format
     df = convert_abundances_format('abund.dat')
-    df.to_csv('converted_abund.csv', index=False, sep=' ')  
+    df.to_csv('converted_abund.csv', index=False, sep=' ')  

docs/examples/tardis_configv1_isotope_abundance_cust_example.yml

@@ -2,4 +2,4 @@ model:
     abundances:
         type: file
         filename: tardis_model_abund.csv
-        filetype: tardis_model
+        filetype: custom_composition

tardis/io/model_reader.py

@@ -1,5 +1,6 @@
 #reading different model files
 
+import warnings
 import numpy as np
 from numpy import recfromtxt, genfromtxt
 import pandas as pd
@@ -43,11 +44,11 @@ def read_density_file(filename, filetype):
     """
     file_parsers = {'artis': read_artis_density,
                     'simple_ascii': read_simple_ascii_density,
-                    'tardis_model': read_cmfgen_density}
+                    'cmfgen_model': read_cmfgen_density}
 
     electron_densities = None
     temperature = None
-    if filetype == 'tardis_model':
+    if filetype == 'cmfgen_model':
         (time_of_model, velocity,
          unscaled_mean_densities, electron_densities, temperature) = read_cmfgen_density(filename)
     else:
@@ -95,11 +96,12 @@ def read_abundances_file(abundance_filename, abundance_filetype,
 
     file_parsers = {'simple_ascii': read_simple_ascii_abundances,
                     'artis': read_simple_ascii_abundances,
-                    'tardis_model': read_simple_isotope_abundances}
+                    'cmfgen_model': read_cmfgen_composition,
+                    'custom_composition': read_csv_composition}
 
     isotope_abundance = pd.DataFrame()
-    if abundance_filetype == 'tardis_model':
-        index, abundances, isotope_abundance = read_simple_isotope_abundances(
+    if abundance_filetype in ["cmfgen_model", "custom_composition"]:
+        index, abundances, isotope_abundance = file_parsers[abundance_filetype](
             abundance_filename)
     else:
         index, abundances = file_parsers[abundance_filetype](
@@ -275,6 +277,9 @@ def read_cmfgen_density(fname):
     temperature: ~np.ndarray
 
     """
+    warnings.warn("The current CMFGEN model parser is deprecated",
+                  DeprecationWarning)
+
     df = pd.read_csv(fname, comment='#', delimiter='\s+', skiprows=[0, 2])
 
     with open(fname) as fh:
@@ -323,18 +328,65 @@ def read_simple_ascii_abundances(fname):
     return index, abundances
 
 
-def read_simple_isotope_abundances(fname, delimiter='\s+'):
+def read_cmfgen_composition(fname, delimiter='\s+'):
+    """Read composition from a CMFGEN model file
+
+    The CMFGEN file format contains information about the ejecta state in the
+    first four columns and the following ones contain elemental and isotopic
+    abundances.
+
+    WARNING: deprecated
+
+    fname: str
+        filename of the csv file
     """
-    Reading an abundance file of the following structure (example; lines starting with hash will be ignored):
-    The first line of abundances describe the abundances in the center of the model and are not used.
-    First 4 columns contain values related to velocity, density, electron_density and temperature.
-    From 5th column onwards, abundances of elements and isotopes begin.
-    The file consists of a header row and next row contains unit of the respective attributes
-    Since abundance fractions are unitless , its unit row is filled with ones
-    Example 
-    velocity...temperature C O Ni56
-    km/s.........K         1 1 1
-    ...................... 0.4 0.3 0.2
+
+    warnings.warn("The current CMFGEN model parser is deprecated",
+                  DeprecationWarning)
+
+    return read_csv_isotope_abundances(fname, delimiter=delimiter,
+                                       skip_columns=4, skip_rows=[0, 2])
+
+def read_csv_composition(fname, delimiter='\s+'):
+    """Read composition from a simple CSV file
+
+    The CSV file can contain specific isotopes or elemental abundances in the
+    different columns. The first row must contain the header in which the
+    contents of each column is specified by the elemental symbol (for elemental
+    abundances) or by the symbol plus mass number (for isotopic abundances).
+
+    Example: C O Fe Ni56 Co
+
+    The i-th row specifies the composition in the i-th shell
+
+    fname: str
+        filename of the csv file
+    """
+
+    return read_csv_isotope_abundances(fname, delimiter=delimiter,
+                                       skip_columns=0, skip_rows=[1])
+
+
+def read_csv_isotope_abundances(fname, delimiter='\s+', skip_columns=0,
+                                skip_rows=[1]):
+    """
+    A generic parser for a TARDIS composition stored as a CSV file
+
+    The parser can read in both elemental and isotopic abundances. It also
+    allows for additional information to be stored in the first skip_columns
+    columns. These will be ignored if skip_columns > 0. Specific header lines
+    can be skipped by the skip_rows keyword argument
+
+    It is expected that the first row of the date block (after skipping the
+    rows specified in skip_rows) specifies the different elements and isotopes.
+    Each row after contains the composition in the corresponding grid shell.
+    The first composition row describes the composition of the photosphere and
+    is essentially ignored (for the default value of skip_rows).
+
+    Example
+    C O Ni56
+    1 1 1
+    0.4 0.3 0.2
 
     Parameters
     ----------
@@ -347,10 +399,10 @@ def read_simple_isotope_abundances(fname, delimiter='\s+'):
 
     index: ~np.ndarray
     abundances: ~pandas.DataFrame
-    isotope_abundance: ~pandas.MultiIndex    
+    isotope_abundance: ~pandas.MultiIndex
     """
     df = pd.read_csv(fname, comment='#',
-                     delimiter=delimiter, skiprows=[0, 2])
+                     delimiter=delimiter, skiprows=skip_rows)
     df = df.transpose()
 
     abundance = pd.DataFrame(columns=np.arange(df.shape[1] - 1),
@@ -363,9 +415,8 @@ def read_simple_isotope_abundances(fname, delimiter='\s+'):
     isotope_abundance = pd.DataFrame(columns=np.arange(df.shape[1] - 1),
                                      index=isotope_index,
                                      dtype=np.float64)
-
-    #First 4 columns related to density parser (e.g. velocity)
-    for element_symbol_string in df.index[4:]:
+
+    for element_symbol_string in df.index[skip_columns:]:
         if element_symbol_string in nucname.name_zz:
             z = nucname.name_zz[element_symbol_string]
             abundance.loc[z, :] = df.loc[element_symbol_string].tolist()[1:]

tardis/io/schemas/model.yml

@@ -104,7 +104,7 @@ definitions:
           enum:
           - simple_ascii
           - artis
-          - tardis_model
+          - cmfgen_model
           description: file type
         v_inner_boundary:
           type: quantity

tardis/io/tests/data/csv_composition.csv

@@ -0,0 +1,12 @@
+C O Mg Si Ni56 Ni58
+1 1 1 1 1 1
+0 0 0 0.6 0.4 0
+0 0 0 0.1 0.5 0.4
+0 0 0 0.3 0 0.7
+0 0.2 0.8 0 0 0
+0 0.3 0.7 0 0 0
+0 0.2 0.8 0 0 0
+0 0.2 0.8 0 0 0
+0 0.2 0.8 0 0 0
+0.5 0.5 0 0 0 0
+0.5 0.5 0 0 0 0

tardis/io/tests/data/tardis_configv1_tardis_model_format.yml

@@ -9,13 +9,13 @@ model:
 
     structure:
         type: file
-        filename: tardis_model_format.csv
-        filetype: tardis_model
+        filename: cmfgen_model.csv
+        filetype: cmfgen_model
 
     abundances:
         type: file
-        filename: tardis_model_format.csv
-        filetype: tardis_model
+        filename: cmfgen_model.csv
+        filetype: cmfgen_model
 
 plasma:
     ionization: lte

tardis/io/tests/test_model_reader.py

@@ -7,7 +7,7 @@
 import tardis
 from tardis.io.config_reader import Configuration
 from tardis.io.model_reader import (
-    read_artis_density, read_simple_ascii_abundances, read_simple_isotope_abundances, read_uniform_abundances, read_cmfgen_density)
+    read_artis_density, read_simple_ascii_abundances, read_csv_composition, read_uniform_abundances, read_cmfgen_density, read_cmfgen_composition)
 
 data_path = os.path.join(tardis.__path__[0], 'io', 'tests', 'data')
 
@@ -19,10 +19,13 @@ def artis_density_fname():
 def artis_abundances_fname():
     return os.path.join(data_path, 'artis_abundances.dat')
 
+@pytest.fixture
+def cmfgen_fname():
+    return os.path.join(data_path, 'cmfgen_model.csv')
 
 @pytest.fixture
-def tardis_model_fname():
-    return os.path.join(data_path, 'tardis_model_format.csv')
+def csv_composition_fname():
+    return os.path.join(data_path, 'csv_composition.csv')
 
 
 @pytest.fixture
@@ -48,16 +51,27 @@ def test_read_simple_ascii_abundances(artis_abundances_fname):
     assert np.isclose(abundances[23].ix[2], 2.672351e-08 , atol=1.e-12)
 
 
-def test_read_simple_isotope_abundances(tardis_model_fname):
-    index, abundances, isotope_abundance = read_simple_isotope_abundances(
-        tardis_model_fname)
+def test_read_simple_isotope_abundances(csv_composition_fname):
+    index, abundances, isotope_abundance = read_csv_composition(
+        csv_composition_fname)
     assert np.isclose(abundances.loc[6, 8], 0.5, atol=1.e-12)
     assert np.isclose(abundances.loc[12, 5], 0.8, atol=1.e-12)
     assert np.isclose(abundances.loc[14, 1], 0.3, atol=1.e-12)
     assert np.isclose(isotope_abundance.loc[(28, 56), 0], 0.5, atol=1.e-12)
     assert np.isclose(isotope_abundance.loc[(28, 58), 1], 0.7, atol=1.e-12)
 
 
+def test_read_cmfgen_isotope_abundances(cmfgen_fname):
+    index, abundances, isotope_abundance = read_cmfgen_composition(
+        cmfgen_fname)
+    assert np.isclose(abundances.loc[6, 8], 0.5, atol=1.e-12)
+    assert np.isclose(abundances.loc[12, 5], 0.8, atol=1.e-12)
+    assert np.isclose(abundances.loc[14, 1], 0.3, atol=1.e-12)
+    assert np.isclose(isotope_abundance.loc[(28, 56), 0], 0.5, atol=1.e-12)
+    assert np.isclose(isotope_abundance.loc[(28, 58), 1], 0.7, atol=1.e-12)
+
+
+
 def test_read_uniform_abundances(isotope_uniform_abundance):
     abundances, isotope_abundance = read_uniform_abundances(
         isotope_uniform_abundance, 20)
@@ -67,9 +81,9 @@ def test_read_uniform_abundances(isotope_uniform_abundance):
     assert np.isclose(isotope_abundance.loc[(28, 58), 2], 0.05, atol=1.e-12)
 
 
-def test_simple_read_cmfgen_density(tardis_model_fname):
+def test_simple_read_cmfgen_density(cmfgen_fname):
     time_of_model, velocity, mean_density, electron_densities, temperature = read_cmfgen_density(
-        tardis_model_fname)
+        cmfgen_fname)
 
     assert np.isclose(0.976 * u.day, time_of_model, atol=1e-7 * u.day)
     assert np.isclose(mean_density[4], 4.2539537e-09 * u.g / u.cm**3, atol=1.e-6