Added support for the T/Y dwarf model spectra from Morley et al. (2012)

docs/overview.rst

@@ -27,6 +27,7 @@ The following data and models are currently supported:
 - `DRIFT-PHOENIX <http://svo2.cab.inta-csic.es/theory/newov2/index.php?models=drift>`_
 - `Exo-REM <https://ui.adsabs.harvard.edu/abs/2018ApJ...854..172C/abstract>`_
 - `petitCODE <http://www.mpia.de/~molliere/#petitcode>`_
+- `T/Y dwarf spectra from Morley et al. (2012) <http://svo2.cab.inta-csic.es/svo/theory/newov2/index.php?models=morley12>`_
 
 **Spectral libraries**
 

species/data/database.py

@@ -21,7 +21,7 @@
 from species.data import ames_cond, ames_dusty, atmo, blackbody, btcond, btcond_feh, btnextgen, \
                          btsettl, btsettl_cifist, companions, drift_phoenix, dust, exo_rem, \
                          filters, irtf, isochrones, leggett, petitcode, spex, vega, vlm_plx, \
-                         kesseli2017
+                         kesseli2017, morley2012
 from species.read import read_calibration, read_filter, read_model, read_object, read_planck
 from species.util import data_util, dust_util, read_util
 
@@ -325,27 +325,31 @@ def add_model(self,
                   teff_range: Optional[Tuple[float, float]] = None,
                   data_folder: Optional[str] = None) -> None:
         """
+        Method for adding a grid of model spectra to the database. All spectra have been resampled
+        to a lower, constant spectral resolution (typically :math:`R = 5000`).
+
         Parameters
         ----------
         model : str
             Model name ('ames-cond', 'ames-dusty', 'atmo', 'bt-settl', 'bt-settl-cifist',
             'bt-nextgen', 'drift-phoenix', 'petitcode-cool-clear', 'petitcode-cool-cloudy',
             'petitcode-hot-clear', 'petitcode-hot-cloudy', 'exo-rem', 'blackbody', bt-cond',
-            or 'bt-cond-feh).
+            'bt-cond-feh, 'morley-2012').
         wavel_range : tuple(float, float), None
-            Wavelength range (um). Optional for the DRIFT-PHOENIX and petitCODE models. For
-            these models, the original wavelength points are used if set to ``None``.
-            which case the argument can be set to ``None``.
+            Wavelength range (um) for adding a subset of the spectra. The full wavelength range
+            is used if the argument is set to ``None``.
         spec_res : float, None
-            Spectral resolution. The parameter is optional for the DRIFT-PHOENIX, petitCODE,
-            BT-Settl, and Exo-REM models. The argument is only used if ``wavel_range`` is not
+            Spectral resolution to which the spectra will be resampled. This parameter is optional
+            since the spectra have already been resampled to a lower, constant resolution
+            (typically :math:`R = 5000`). The argument is only used if ``wavel_range`` is not
             ``None``.
         teff_range : tuple(float, float), None
-            Effective temperature range (K). Setting the value to None for will add all available
-            temperatures.
+            Effective temperature range (K) for adding a subset of the model grid. The full
+            parameter grid will be added if the argument is set to ``None``.
         data_folder : str, None
-            Folder with input data. This parameter has been deprecated since all model spectra
-            are publicly available.
+            DEPRECATED: Folder where the input data is located. This parameter is no longer in use
+            since all model spectra are publicly available. The parameter will be removed in a
+            future release.
 
         Returns
         -------
@@ -473,6 +477,15 @@ def add_model(self,
 
             data_util.add_missing(model, ['teff', 'logg', 'feh'], h5_file)
 
+        elif model == 'morley-2012':
+            morley2012.add_morley2012(self.input_path,
+                                      h5_file,
+                                      wavel_range,
+                                      teff_range,
+                                      spec_res)
+
+            data_util.add_missing(model, ['teff', 'logg', 'fsed'], h5_file)
+
         elif model == 'petitcode-cool-clear':
             petitcode.add_petitcode_cool_clear(self.input_path,
                                                h5_file,
@@ -908,7 +921,7 @@ def add_object(self,
 
                                 corr_warn = f'The covariance matrix from {value[1]} contains ' \
                                             f'ones along the diagonal. Converting this ' \
-                                            f'correlation  matrix into a covariance matrix.'
+                                            f'correlation matrix into a covariance matrix.'
 
                                 if data.ndim == 2 and data.shape[0] == data.shape[1]:
                                     if key not in read_cov:

species/data/morley2012.py

@@ -0,0 +1,157 @@
+"""
+Module for T/Y dwarf model spectra from Morley et al. (2012).
+"""
+
+import os
+import tarfile
+import urllib.request
+
+from typing import Optional, Tuple
+
+import h5py
+import spectres
+import numpy as np
+
+from typeguard import typechecked
+
+from species.util import data_util, read_util
+
+
+@typechecked
+def add_morley2012(input_path: str,
+                   database: h5py._hl.files.File,
+                   wavel_range: Optional[Tuple[float, float]],
+                   teff_range: Optional[Tuple[float, float]],
+                   spec_res: Optional[float]) -> None:
+    """
+    Function for adding the T/Y dwarf model spectra from Morley et al. (2012) to the database.
+    The spectra have been at solar metallicity in the Teff range from 500 to 1300 K. The spectra
+    have been downloaded from the Theoretical spectra web server
+    (http://svo2.cab.inta-csic.es/svo/theory/newov2/index.php?models=morley12) and resampled
+    to a spectral resolution of 5000 from 0.6 to 30 um.
+
+    Parameters
+    ----------
+    input_path : str
+        Folder where the data is located.
+    database : h5py._hl.files.File
+        Database.
+    wavel_range : tuple(float, float), None
+        Wavelength range (um). The original wavelength points are used if set to ``None``.
+    teff_range : tuple(float, float), None
+        Effective temperature range (K). All temperatures are selected if set to ``None``.
+    spec_res : float, None
+        Spectral resolution for resampling. Not used if ``wavel_range`` is set to ``None``.
+
+    Returns
+    -------
+    NoneType
+        None
+    """
+
+    if not os.path.exists(input_path):
+        os.makedirs(input_path)
+
+    input_file = 'morley2012.tgz'
+
+    data_folder = os.path.join(input_path, 'morley2012/')
+    data_file = os.path.join(input_path, input_file)
+
+    if not os.path.exists(data_folder):
+        os.makedirs(data_folder)
+
+    url = 'https://home.strw.leidenuniv.nl/~stolker/species/morley2012.tgz'
+
+    if not os.path.isfile(data_file):
+        print('Downloading Morley et al. (2012) model spectra (141 MB)...', end='', flush=True)
+        urllib.request.urlretrieve(url, data_file)
+        print(' [DONE]')
+
+    print('Unpacking Morley et al. (2012) model spectra (141 MB)...', end='', flush=True)
+    tar = tarfile.open(data_file)
+    tar.extractall(data_folder)
+    tar.close()
+    print(' [DONE]')
+
+    teff = []
+    logg = []
+    fsed = []
+    flux = []
+
+    if wavel_range is not None and spec_res is not None:
+        wavelength = read_util.create_wavelengths(wavel_range, spec_res)
+    else:
+        wavelength = None
+
+    print_message = ''
+
+    for _, _, file_list in os.walk(data_folder):
+        for filename in sorted(file_list):
+            if filename[:11] == 'morley2012_':
+                file_split = filename.split('_')
+
+                teff_val = float(file_split[2])
+                logg_val = float(file_split[4])
+                fsed_val = float(file_split[6])
+
+                if teff_range is not None:
+                    if teff_val < teff_range[0] or teff_val > teff_range[1]:
+                        continue
+
+                empty_message = len(print_message) * ' '
+                print(f'\r{empty_message}', end='')
+
+                print_message = f'Adding Morley et al. (2012) model spectra... {filename}'
+                print(f'\r{print_message}', end='')
+
+                data_wavel, data_flux = np.loadtxt(os.path.join(data_folder, filename), unpack=True)
+                print(data_wavel)
+
+                teff.append(teff_val)
+                logg.append(logg_val)
+                fsed.append(fsed_val)
+
+                if wavel_range is None or spec_res is None:
+                    if wavelength is None:
+                        wavelength = np.copy(data_wavel)  # (um)
+
+                    if np.all(np.diff(wavelength) < 0):
+                        raise ValueError('The wavelengths are not all sorted by increasing value.')
+
+                    flux.append(data_flux)  # (W m-2 um-1)
+
+                else:
+                    flux_resample = spectres.spectres(wavelength,
+                                                      data_wavel,
+                                                      data_flux,
+                                                      spec_errs=None,
+                                                      fill=np.nan,
+                                                      verbose=False)
+
+                    if np.isnan(np.sum(flux_resample)):
+                        raise ValueError(f'Resampling is only possible if the new wavelength '
+                                         f'range ({wavelength[0]} - {wavelength[-1]} um) falls '
+                                         f'sufficiently far within the wavelength range '
+                                         f'({data_wavel[0]} - {data_wavel[-1]} um) of the input '
+                                         f'spectra.')
+
+                    flux.append(flux_resample)  # (W m-2 um-1)
+
+    empty_message = len(print_message) * ' '
+    print(f'\r{empty_message}', end='')
+
+    print_message = 'Adding Morley et al. (2012) model spectra... [DONE]'
+    print(f'\r{print_message}', end='')
+
+    data_sorted = data_util.sort_data(np.asarray(teff),
+                                      np.asarray(logg),
+                                      None,
+                                      None,
+                                      np.asarray(fsed),
+                                      wavelength,
+                                      np.asarray(flux))
+
+    data_util.write_data('morley-2012',
+                         ['teff', 'logg', 'fsed'],
+                         database,
+                         data_sorted)