Add .splups file reader (#278)

.github/workflows/ci.yml

@@ -38,6 +38,14 @@ jobs:
         - linux: py311
       cache-path: ~/.chianti
       cache-key: chianti-${{ github.event.number }}
+  test_database_v7:
+    needs: [test]
+    uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@v1
+    with:
+      posargs: '--ascii-dbase-root ~/.chianti --ascii-dbase-url http://download.chiantidatabase.org/CHIANTI_v7.1.4_database.tar.gz --disable-file-hash --skip-version-check'
+      toxdeps: "'tox<4' tox-pypi-filter"
+      envs: |
+        - linux: py311
   test_database_v9:
     needs: [test]
     uses: OpenAstronomy/github-actions-workflows/.github/workflows/tox.yml@v1

fiasco/conftest.py

@@ -82,27 +82,36 @@ def dbase_version(ascii_dbase_root):
 
 @pytest.fixture(autouse=True)
 def requires_dbase_version(request, dbase_version):
+    """
+    Skip a test if not all version requirements are met. Multiple requirements are joined by AND.
+    """
     # NOTE: Fixtures that depend on other fixtures are awkward to implement.
     # See this SO answer: https://stackoverflow.com/a/28198398
     if marker := request.node.get_closest_marker('requires_dbase_version'):
-        # NOTE: This has to have a space between the operator and the target
-        if  len(marker.args) != 2:
-            raise ValueError("Arguments must contain a condition and a version number, e.g. '<', '8.0.7'")
-        operator, target_version = marker.args
         op_dict = {'<': np.less,
                    '<=': np.less_equal,
                    '>': np.greater,
                    '>=': np.greater_equal,
                    '=': np.equal,
                    '==': np.equal,
                    '!=': np.not_equal}
-        if operator not in op_dict:
-            raise ValueError(f'''{operator} is not a supported comparison operation.
-                                 Must be one of {list(op_dict.keys())}.''')
-        target_version = Version(target_version)
-        allowed_dbase_version = op_dict[operator](dbase_version, target_version)
-        if not allowed_dbase_version:
-            pytest.skip(f'Skip because database version {dbase_version} is not {operator} {target_version}.')
+
+        def _evaluate_condtion(condition_string):
+            condition_array = condition_string.split()
+            if  len(condition_array) != 2:
+                raise ValueError("Arguments must contain a condition and a version number with a space, e.g. '< 8.0.7'")
+            operator, target_version = condition_array
+            if operator not in op_dict:
+                raise ValueError(f'''{operator} is not a supported comparison operation.
+                                    Must be one of {list(op_dict.keys())}.''')
+            target_version = Version(target_version)
+            allowed_dbase_version = op_dict[operator](dbase_version, target_version)
+            return allowed_dbase_version, operator, target_version
+
+        conditions = np.atleast_1d(marker.args)
+        for is_met, operator, target_version in list(map(_evaluate_condtion, conditions)):
+            if not is_met:
+                pytest.skip(f'Skipping because database version {dbase_version} is not {operator} {target_version}.')
 
 
 def pytest_configure(config):

fiasco/io/sources/ion_sources.py

@@ -11,6 +11,7 @@
     'ElvlcParser',
     'FblvlParser',
     'ScupsParser',
+    'SplupsParser',
     'PsplupsParser',
     'EasplomParser',
     'EasplupsParser',
@@ -132,6 +133,70 @@ def postprocessor(self, df):
         return df
 
 
+class SplupsParser(GenericIonParser):
+    """
+    Spline fits to scaled collisions strengths (denoted by upsilon) between energy levels as described
+    in :cite:t:`burgess_analysis_1992`.  These files were used in CHIANTI versions prior to 8.0, and
+    were replaced by ``.scups`` files in versions after that.
+
+    Notes
+    -----
+    * The number of spline points for the rates depends on the fit type, 5 points for type 6
+      fits and 9 points for type 2.
+    """
+    filetype = 'splups'
+    dtypes = [int, int, int, int, int, float, float, float, 'object']
+    units = [
+        None,
+        None,
+        None,
+        None,
+        None,
+        u.dimensionless_unscaled,
+        u.Ry,
+        u.dimensionless_unscaled,
+        u.dimensionless_unscaled,
+    ]
+    headings = [
+        'Z',
+        'ion',
+        'lower_level',
+        'upper_level',
+        'bt_type',
+        'gf',
+        'delta_energy',
+        'bt_c',
+        'bt_upsilon',
+    ]
+    descriptions = [
+        'atomic number',
+        'ionization state',
+        'lower level index',
+        'upper level index',
+        'Burgess-Tully scaling type',
+        'oscillator strength',
+        'delta energy',
+        'Burgess-Tully scaling parameter',
+        'Burgess-Tully scaled effective collision strengths',
+    ]
+
+    def preprocessor(self, table, line, index):
+        n_spline =  9  # Max number of spline points
+        fformat = fortranformat.FortranRecordReader(f'(5I3,{3+n_spline}E10.3)')
+        line = fformat.read(line)
+        # NOTE: The first eight entries are fixed. The last entry is the scaled
+        # spline fit to the array and can vary in length.
+        # NOTE: Some spline fits only have 5 points and the scaling type is not
+        # a reliable way to determine this so we have to filter these manually.
+        # When fortranformat has missing entries, it fills them in as None. We
+        # remove them here to avoid the undefined behavior of None in a ragged
+        # array within an astropy Table.
+        spline_fit = line[8:]
+        spline_fit = [sf for sf in spline_fit if sf is not None]
+        row = line[:8] + [np.array(spline_fit)]
+        table.append(row)
+
+
 class PsplupsParser(ScupsParser):
     """
     Spline fits to scaled collision rates for protons. These files are discussed in

fiasco/io/sources/tests/test_sources.py

@@ -11,7 +11,7 @@
 @pytest.mark.parametrize('filename', [
     'h_1.elvlc',
     'h_1.fblvl',
-    'h_1.scups',
+    pytest.param('h_1.scups', marks=pytest.mark.requires_dbase_version('>= 8')),
     'c_2.psplups',
     'be_2.easplom',
     'al_3.easplups',
@@ -22,8 +22,10 @@
     'fe_2.trparams',
     'fe_12.drparams',
     'al_3.diparams',
-    pytest.param('fe_23.auto', marks=pytest.mark.requires_dbase_version('>=', '9')),
-    pytest.param('fe_23.rrlvl', marks=pytest.mark.requires_dbase_version('>=', '9')),
+    pytest.param('fe_23.auto', marks=pytest.mark.requires_dbase_version('>= 9')),
+    pytest.param('fe_23.rrlvl', marks=pytest.mark.requires_dbase_version('>= 9')),
+    pytest.param('c_5.splups', marks=pytest.mark.requires_dbase_version('< 8')),
+    pytest.param('c_6.splups', marks=pytest.mark.requires_dbase_version('< 8')),
 ])
 def test_ion_sources(ascii_dbase_root, filename,):
     parser = fiasco.io.Parser(filename, ascii_dbase_root=ascii_dbase_root)

fiasco/tests/test_collections.py

@@ -106,6 +106,7 @@ def test_free_bound(another_collection, wavelength):
     index_t = 24  # This is approximately where the ioneq for Fe V peaks
     assert u.allclose(fb[index_t, index_w], 3.057781475607237e-36 * u.Unit('erg cm3 s-1 Angstrom-1'))
 
+@pytest.mark.requires_dbase_version('>= 8')
 @pytest.mark.parametrize('wavelength', [wavelength, wavelength[450]])
 def test_two_photon(collection, wavelength, hdf5_dbase_root):
     # add Li III to the test to include an ion that throws a MissingDatasetException
@@ -120,12 +121,13 @@ def test_two_photon(collection, wavelength, hdf5_dbase_root):
     # This value has not been checked for correctness
     assert u.allclose(tp[index_w, index_t, 0], 3.48580645e-27 * u.Unit('erg cm3 s-1 Angstrom-1'))
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_radiative_loss(collection):
     rl = collection.radiative_loss(1e9*u.cm**(-3))
     # This value has not been checked for correctness
     assert u.allclose(rl[0,0], 3.90235371e-24*u.Unit('erg cm3 s-1'))
 
-
+@pytest.mark.requires_dbase_version('>= 8')
 def test_spectrum(hdf5_dbase_root):
     i1 = fiasco.Ion('H 1', 1 * u.MK, hdf5_dbase_root=hdf5_dbase_root)
     i2 = fiasco.Ion('Fe 5', 1 * u.MK, hdf5_dbase_root=hdf5_dbase_root)

fiasco/tests/test_ion.py

@@ -85,6 +85,7 @@ def test_level_indexing(ion):
     assert levels[2].__repr__() == fiasco.Level(10, ion._elvlc).__repr__()
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_level(ion):
     level = ion[0]
     assert isinstance(level, fiasco.Level)
@@ -157,6 +158,7 @@ def test_abundance(ion):
     assert u.allclose(ion.abundance, 0.0001258925411794166)
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_proton_collision(fe10):
     rate = fe10.proton_collision_excitation_rate
     assert u.allclose(rate[0, 0], 4.69587161e-13 * u.cm**3 / u.s)
@@ -184,6 +186,7 @@ def test_missing_ip(hdf5_dbase_root):
         _ = ion.ip
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_level_populations(ion):
     pop = ion.level_populations(1e8 * u.cm**-3)
     assert pop.shape == ion.temperature.shape + (1,) + ion._elvlc['level'].shape
@@ -193,19 +196,22 @@ def test_level_populations(ion):
     assert u.allclose(pop.squeeze().sum(axis=1), 1, atol=None, rtol=1e-15)
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_level_populations_proton_data_toggle(ion):
     # Fe V has no psplups data so the toggle should have no effect
     lp_protons = ion.level_populations(1e9*u.cm**(-3), include_protons=True)
     lp_no_protons = ion.level_populations(1e9*u.cm**(-3), include_protons=False)
     assert u.allclose(lp_protons, lp_no_protons, atol=0, rtol=0)
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_contribution_function(ion):
     cont_func = ion.contribution_function(1e7 * u.cm**-3)
     assert cont_func.shape == ion.temperature.shape + (1, ) + ion._wgfa['wavelength'].shape
     # This value has not been tested for correctness
     assert u.allclose(cont_func[0, 0, 0], 2.51408088e-30 * u.cm**3 * u.erg / u.s)
 
+@pytest.mark.requires_dbase_version('>= 8')
 @pytest.mark.parametrize(('density', 'use_coupling'), [
     ([1e9,] * u.cm**(-3), False),
     ([1e8, 1e9, 1e10] * u.cm**(-3), False),
@@ -231,6 +237,7 @@ def test_emissivity_shape(c6, density, use_coupling):
     assert emiss.shape == c6.temperature.shape + density_shape + wavelength.shape
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_coupling_unequal_dimensions_exception(ion):
     with pytest.raises(ValueError, match='Temperature and density must be of equal length'):
         _ = ion.level_populations([1e7, 1e8]*u.cm**(-3), couple_density_to_temperature=True)
@@ -247,12 +254,13 @@ def pops_no_correction(fe20):
                                   include_level_resolved_rate_correction=False).squeeze()
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_level_populations_normalized(pops_no_correction, pops_with_correction):
     assert u.allclose(pops_with_correction.sum(axis=1), 1, atol=None, rtol=1e-15)
     assert u.allclose(pops_no_correction.sum(axis=1), 1, atol=None, rtol=1e-15)
 
 
-@pytest.mark.requires_dbase_version('<=','8.0.7')
+@pytest.mark.requires_dbase_version('>= 8', '<= 8.0.7')
 def test_level_populations_correction(fe20, pops_no_correction, pops_with_correction):
     # Test level-resolved correction applied to correct levels
     i_corrected = np.unique(np.concatenate([fe20._cilvl['upper_level'], fe20._reclvl['upper_level']]))
@@ -270,13 +278,15 @@ def test_level_populations_correction(fe20, pops_no_correction, pops_with_correc
                       atol=0.0, rtol=1e-5)
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_emissivity(ion):
     emm = ion.emissivity(1e7 * u.cm**-3)
     assert emm.shape == ion.temperature.shape + (1, ) + ion._wgfa['wavelength'].shape
     # This value has not been tested for correctness
     assert u.allclose(emm[0, 0, 0], 2.18000422e-16 * u.erg / u.cm**3 / u.s)
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 @pytest.mark.parametrize('em', [
     1e29 * u.cm**-5,
     [1e29] * u.cm**-5,
@@ -334,7 +344,7 @@ def test_free_bound(ion):
     assert u.allclose(emission[0, 0], 9.7902609e-26 * u.cm**3 * u.erg / u.Angstrom / u.s)
 
 # The two-photon test currently fails for dbase_version >= 9 because it is missing c_5.reclvl
-@pytest.mark.requires_dbase_version('<=','8.0.7')
+@pytest.mark.requires_dbase_version('>= 8','<= 8.0.7')
 def test_two_photon(c4, c5, c6):
     # test both H-like and He-like ions, and one that doesn't have two-photon emission
     c4_emission = c4.two_photon(200 * u.Angstrom, electron_density = 1e10* u.cm**(-3))
@@ -374,6 +384,7 @@ def test_radd_ions(ion, another_ion):
     assert collection[0] == another_ion
 
 
+@pytest.mark.requires_dbase_version('>= 8')
 def test_transitions(ion):
     trans = ion.transitions
     assert isinstance(trans, fiasco.Transitions)