Model spectra and synthetic photometry#

In this tutorial, we will have a look at some spectra of the DRIFT-PHOENIX atmospheric model. The spectra are first downloaded and added to the database. Then we will use the functionalities of ReadModel to extract a spectrum and calculate a photometric flux.

Getting started#

We start by importing the required Python modules.

[1]:
import species

Then we initialize species with SpeciesInit, which creates a default configuration file and the HDF5 database.

[2]:
species.SpeciesInit()
Initiating species v0.5.2... [DONE]
Creating species_config.ini... [DONE]
Database: /Users/tomasstolker/applications/species/docs/tutorials/species_database.hdf5
Data folder: /Users/tomasstolker/applications/species/docs/tutorials/data
Working folder: /Users/tomasstolker/applications/species/docs/tutorials
Grid interpolation method: linear
Creating species_database.hdf5... [DONE]
Creating data folder... [DONE]
[2]:
<species.core.init.SpeciesInit at 0x14ba61c90>

Adding model spectra to the database#

To store the spectra, we first create an instance of Database.

[3]:
database = species.Database()

Let’s check which atmospheric models are available by running the available_models method of the Database object.

[4]:
_ = database.available_models()
Available model grids:
   - AMES-Cond:
      - Label = ames-cond
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [100, 6600]
      - Wavelength range (um): [0.5, 40]
      - Resolution lambda/Dlambda: 4000
      - File size: 150 MB

   - AMES-Dusty:
      - Label = ames-dusty
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [500, 4000]
      - Wavelength range (um): [0.5, 40]
      - Resolution lambda/Dlambda: 4000
      - File size: 59 MB

   - ATMO:
      - Label = atmo
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [200, 3000]
      - Wavelength range (um): [0.4, 6000]
      - Resolution lambda/Dlambda: 10000
      - File size: 425 MB
      - Reference: Phillips et al. (2020)
      - URL: https://ui.adsabs.harvard.edu/abs/2020A%26A...637A..38P/abstract

   - blackbody:
      - Label = blackbody
      - Model parameters: ['teff']
      - Teff range (K): [10, 20000]
      - Wavelength range (um): [0.1, 5000]
      - Resolution lambda/Dlambda: 1000
      - File size: 56 MB

   - BT-Cond:
      - Label = bt-cond
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [800, 4000]
      - Wavelength range (um): [0.1, 100]
      - Resolution lambda/Dlambda: 10000
      - File size: 162 MB

   - BT-Cond [Fe/H]:
      - Label = bt-cond-feh
      - Model parameters: ['teff', 'logg', 'feh']
      - Teff range (K): [2600, 4000]
      - Wavelength range (um): [0.1, 100]
      - Resolution lambda/Dlambda: 10000
      - File size: 390 MB

   - BT-NextGen:
      - Label = bt-nextgen
      - Model parameters: ['teff', 'logg', 'feh']
      - Teff range (K): [2600, 30000]
      - Wavelength range (um): [0.1, 50]
      - Resolution lambda/Dlambda: 4000
      - File size: 368 MB

   - BT-Settl:
      - Label = bt-settl
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [400, 4000]
      - Wavelength range (um): [0.1, 100]
      - Resolution lambda/Dlambda: 10000
      - File size: 227 MB

   - BT-Settl CIFIST:
      - Label = bt-settl-cifist
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [1200, 7000]
      - Wavelength range (um): [0.1, 5000]
      - Resolution lambda/Dlambda: 20000
      - File size: 1.8 GB

   - DRIFT-PHOENIX:
      - Label = drift-phoenix
      - Model parameters: ['teff', 'logg', 'feh']
      - Teff range (K): [1000, 3000]
      - Wavelength range (um): [0.1, 50]
      - Resolution lambda/Dlambda: 4000
      - File size: 229 MB
      - Reference: Helling et al. (2008)
      - URL: https://ui.adsabs.harvard.edu/abs/2008ApJ...675L.105H/abstract

   - Exo-REM:
      - Label = exo-rem
      - Model parameters: ['teff', 'logg', 'feh', 'c_o_ratio']
      - Teff range (K): [1000, 2000]
      - Wavelength range (um): [0.9, 5]
      - Resolution lambda/Dlambda: 10000
      - File size: 706 MB
      - Reference: Charney et al. (2018)
      - URL: https://ui.adsabs.harvard.edu/abs/2018ApJ...854..172C/abstract

   - Koester WD:
      - Label = koester-wd
      - Model parameters: ['teff', 'logg']
      - Teff range (K): [5000, 80000]
      - Wavelength range (um): [0.09, 3]
      - Resolution lambda/Dlambda: 2000
      - File size: 291 MB
      - Reference: Koester (2010)
      - URL: https://ui.adsabs.harvard.edu/abs/2010MmSAI..81..921K/abstract

   - Morley et al. (2012):
      - Label = morley-2012
      - Model parameters: ['teff', 'logg', 'fsed']
      - Teff range (K): [400, 1300]
      - Wavelength range (um): [0.6, 30]
      - Resolution lambda/Dlambda: 10000
      - File size: 141 MB
      - Reference: Morley et al. (2012)
      - URL: https://ui.adsabs.harvard.edu/abs/2012ApJ...756..172M/abstract

   - petitCODE cool clear:
      - Label = petitcode-cool-clear
      - Model parameters: ['teff', 'logg', 'feh']
      - Teff range (K): [500, 1000]
      - Wavelength range (um): [0.1, 250]
      - Resolution lambda/Dlambda: 1000
      - File size: 86 MB
      - Reference: Mollière et al. (2015)
      - URL: https://ui.adsabs.harvard.edu/abs/2015ApJ...813...47M/abstract

   - petitCODE cool cloudy:
      - Label = petitcode-cool-cloudy
      - Model parameters: ['teff', 'logg', 'feh', 'fsed']
      - Teff range (K): [500, 850]
      - Wavelength range (um): [0.1, 250]
      - Resolution lambda/Dlambda: 1000
      - Extra details: log(Kzz) = 7.5
      - File size: 224 MB
      - Reference: Mollière et al. (2015)
      - URL: https://ui.adsabs.harvard.edu/abs/2015ApJ...813...47M/abstract

   - petitCODE hot clear:
      - Label = petitcode-hot-clear
      - Model parameters: ['teff', 'logg', 'feh', 'c_o_ratio']
      - Teff range (K): [1000, 2000]
      - Wavelength range (um): [0.1, 250]
      - Resolution lambda/Dlambda: 1000
      - File size: 42 MB
      - Reference: Mollière et al. (2015)
      - URL: https://ui.adsabs.harvard.edu/abs/2015ApJ...813...47M/abstract

   - petitCODE hot cloudy:
      - Label = petitcode-hot-cloudy
      - Model parameters: ['teff', 'logg', 'feh', 'c_o_ratio', 'fsed']
      - Teff range (K): [1000, 2000]
      - Wavelength range (um): [0.1, 250]
      - Resolution lambda/Dlambda: 1000
      - File size: 126 MB
      - Reference: Mollière et al. (2015)
      - URL: https://ui.adsabs.harvard.edu/abs/2015ApJ...813...47M/abstract

   - Sonora Cholla:
      - Label = sonora-cholla
      - Model parameters: ['teff', 'logg', 'log_kzz']
      - Teff range (K): [500, 1300]
      - Wavelength range (um): [1, 250]
      - Resolution lambda/Dlambda: 10000
      - Extra details: [Fe/H] = 0.0
      - File size: 497 MB
      - Reference: Karalidi et al. (2021)
      - URL: https://zenodo.org/record/4450269

   - Sonora Bobcat:
      - Label = sonora-bobcat
      - Model parameters: ['teff', 'logg', 'feh']
      - Teff range (K): [200, 2400]
      - Wavelength range (um): [0.4, 50]
      - Resolution lambda/Dlambda: 10000
      - File size: 1.1 GB
      - Reference: Marley et al. (2021)
      - URL: https://zenodo.org/record/5063476

   - Sonora Bobcat C/O:
      - Label = sonora-bobcat-co
      - Model parameters: ['teff', 'c_o_ratio']
      - Teff range (K): [200, 2400]
      - Wavelength range (um): [0.4, 50]
      - Resolution lambda/Dlambda: 10000
      - Extra details: [Fe/H] = 0.0, log(g) = 5.0
      - File size: 112 MB
      - Reference: Marley et al. (2021)
      - URL: https://zenodo.org/record/5063476

Next, we will import the model spectra with the add_model method of Database. This step will automatically download the DRIFT-PHOENIX spectra (R = 2000) to the data_folder. The dowloaded data will then be unpacked and added to the database. We restrict the temperature range to 1300-1700 K, so not all spectra are added to the databse.

[5]:
database.add_model(model='drift-phoenix', teff_range=(1300., 1700.))
Downloading DRIFT-PHOENIX model spectra (229 MB)... [DONE]
Unpacking DRIFT-PHOENIX model spectra (229 MB)... [DONE]
Please cite Helling et al. (2008) when using DRIFT-PHOENIX in a publication
Reference URL: https://ui.adsabs.harvard.edu/abs/2008ApJ...675L.105H/abstract
Wavelength range (um) = 0.1 - 50
Spectral resolution = 4000
Teff range (K) = 1300.0 - 1700.0
Adding DRIFT-PHOENIX model spectra... [DONE]
Grid points stored in the database:
   - Teff = [1300. 1400. 1500. 1600. 1700.]
   - log(g) = [3.  3.5 4.  4.5 5.  5.5]
   - [Fe/H] = [-0.6 -0.3 -0.   0.3]

Two of the grid points were missing in the original data and have been added with a linear, multidimensional interpolation.

Interpolating the model grid#

We will read the spectra from the database by creating an instance of ReadModel.

[6]:
read_model = species.ReadModel(model='drift-phoenix', wavel_range=(0.5, 10.))

Let’s see what the grid boundaries are from the spectra that are stored in the database.

[7]:
read_model.get_bounds()
[7]:
{'teff': (1300.0, 1700.0), 'logg': (3.0, 5.5), 'feh': (-0.6, 0.3)}

We will now interpolate the grid in the (\(T_\mathrm{eff}\), \(\log(g)\), \([\mathrm{Fe}/\mathrm{H}]\)) space at some specific parameter values, which need to be provided in a dictionary. The radius and distance are optional, otherwise the emitted flux is given at the planet surface.

[8]:
model_param = {'teff':1510., 'logg':4.1, 'feh':0.1, 'radius': 1., 'distance': 100.}

To interpolate a spectrum, we use the get_model method and provide the parameter dictionary, and also an optional spectral resolution. Together with smooth=True, the spectrum will be smoothed (but not resampeld) to the given spectral resolution.

[9]:
model_box = read_model.get_model(model_param=model_param, spec_res=200., smooth=True)

The data is stored in a ModelBox. Let’s have a look at its content.

[10]:
model_box.open_box()
Opening ModelBox...
model = drift-phoenix
type = None
wavelength = [ 0.49989727  0.50002105  0.50014486 ...  9.99710369  9.99957902
 10.00205496]
flux = [9.44339097e-20 9.38011967e-20 9.31147817e-20 ... 1.53571943e-18
 1.53924349e-18 1.54301688e-18]
parameters = {'teff': 1510.0, 'logg': 4.1, 'feh': 0.1, 'radius': 1.0, 'distance': 100.0, 'luminosity': 4.729862212008143e-05, 'mass': 4.857062223118246}
quantity = flux
contribution = None
bol_flux = None

We will now use the same atmospheric parameters but we will add some ISM extinction with the relation from Cardelli et al. (1989). Therefore, we add the V band extinction and reddening parameters to the dictionary.

[11]:
model_param = {'teff':1510., 'logg':4.1, 'feh':0.1, 'radius': 1., 'distance': 100., 'ism_ext': 5., 'ism_red': 3.}

We use again the get_model method and store the result in a different ModelBox.

[12]:
model_ext = read_model.get_model(model_param=model_param, spec_res=200., smooth=True)

The two boxes with the model spectra are provided to the plot_spectrum. We also include some filter profiles to indicate where the telluric windows are.

[13]:
species.plot_spectrum(boxes=[model_box, model_ext],
                      filters=['MKO/NSFCam.J', 'MKO/NSFCam.H', 'MKO/NSFCam.K', 'MKO/NSFCam.Lp', 'MKO/NSFCam.Mp'],
                      offset=(-0.15, -0.04),
                      xlim=(0.8, 5.),
                      ylim=(0., 5.5e-17),
                      legend={'loc': 'lower right', 'frameon': False, 'fontsize': 10.},
                      figsize=(7., 3.),
                      output=None)
Adding filter: MKO/NSFCam.J... [DONE]
Adding filter: MKO/NSFCam.H... [DONE]
Adding filter: MKO/NSFCam.K... [DONE]
Adding filter: MKO/NSFCam.Lp... [DONE]
Adding filter: MKO/NSFCam.Mp... [DONE]
Plotting spectrum... [DONE]
../_images/tutorials_model_spectra_31_1.png

Extracting a spectrum at a grid point#

It is also possible to extract a spectrum at one of the grid points, which doesn’t require any interpolation. Let’s check with the get_points method what parameter values are stored in the database.

[13]:
read_model.get_points()
[13]:
{'teff': array([1300., 1400., 1500., 1600., 1700.]),
 'logg': array([3. , 3.5, 4. , 4.5, 5. , 5.5]),
 'feh': array([-0.6, -0.3, -0. ,  0.3])}

We create a dictionary with values at one of the grid points.

[14]:
model_param = {'teff':1500., 'logg':4., 'feh':0.}

And now use the get_data method to extract a spectrum.

[15]:
model_full = read_model.get_data(model_param)

Let’s make another plot with plot_spectrum. The spectrum is now shown at the spectral resolution as stored in the database (\(R = 2000\)).

[17]:
species.plot_spectrum(boxes=[model_full],
                      filters=None,
                      offset=(-0.12, -0.05),
                      xlim=(0.8, 5.),
                      ylim=(0., 1e5),
                      legend={'loc': 'upper right', 'fontsize': 10.},
                      figsize=(7., 3.),
                      output=None)
Plotting spectrum... [DONE]
../_images/tutorials_model_spectra_40_1.png

Calculating synthetic photometry#

The ReadModel class can also be used for calculating photometric fluxes and magnitudes. To do so, we create a new instance and set the filter_name argument to the VLT/NACO M’ filter. This will automatically downloadd and addd the filter profile.

[17]:
read_model = species.ReadModel(model='drift-phoenix', filter_name='Paranal/NACO.Mp')
Adding filter: Paranal/NACO.Mp... [DONE]

We create again a dictionary with the parameters but now run the get_flux method, which returns the flux in W m-2 um-1.

[18]:
model_param = {'teff':1510., 'logg':4.1, 'feh':0.1, 'radius': 1., 'distance': 100.}
flux = read_model.get_flux(model_param)
print(f'Flux (W m-2 um-1) = {flux[0]:.2e}')
Flux (W m-2 um-1) = 1.33e-17

Since we provided a radius and distance, the emitted flux at the planet surface has been scaled by \((\mathrm{radius}/\mathrm{distance})^2\).

Similarly, we can use the get_magnitude method to calculate the magnitude for the NACO M’ filter. Note that the returned absolute magnitude is set to None if the parameter dictionary does not contain a radius and distance.

[19]:
app_mag, abs_mag = read_model.get_magnitude(model_param)
print(f'Apparent magnitude = {app_mag:.2f}')
print(f'Absolute magnitude = {abs_mag:.2f}')
Downloading Vega spectrum (270 kB)... [DONE]
Adding Vega spectrum... [DONE]
Apparent magnitude = 15.53
Absolute magnitude = 10.53

As expected, at a distance of 100 pc, the difference between the apparent and absolute magnitude is 5.