To begin, we need a spectral model library that we will use for our fitting. One common example are the PHOENIX models, most recently computed by T.O. Husser. We provide many interfaces directly with different libraries, which can be viewed in Raw Grid Interfaces.
As a convenience, we provide a helper to download PHOENIX models from the Goettingen servers. Note this will skip any files already on disk.
In [1]:
import numpy as np
from Starfish.grid_tools import download_PHOENIX_models
ranges = [[5700, 8600], [4.0, 6.0], [-0.5, 0.5]] # T, logg, Z
download_PHOENIX_models(path="PHOENIX", ranges=ranges)
Now that we have the files downloaded, let's set up a grid interface
In [2]:
from Starfish.grid_tools import PHOENIXGridInterfaceNoAlpha
grid = PHOENIXGridInterfaceNoAlpha(path="PHOENIX")
From here, we will want to set up our HDF5 interface that will allow us to go on to using the spectral emulator, but first we need to determine our model subset and instrument.
We set up an HDF5 interface in order to allow much quicker reading and writing than compared to loading FITS files over and over again. In addition, when considering the application to our likelihood methods, we know that for a given dataset, any effects characteristic of the instrument can be pre-applied to our models, saving on computation time during the maximum likelihood estimation.
Looking towards our fitting examples, we know we will try fitting some data from TRES Spectrograph. This instrument is available in our grid tools, but if yours isn't, you can always supply the FWHM in km/s. The FWHM ($\Gamma$) can be found using the resolving power, $R$
$$ \Gamma = \frac{c}{R} $$with $c$ in km/s. Let’s also say that, for a given dataset, we want to only use a reasonable subset of our original model grid. The data provided in future examples is a ~F3V star, so we will limit our model parameter ranges appropriately.
In [3]:
from Starfish.grid_tools.instruments import SPEX
from Starfish.grid_tools import HDF5Creator
creator = HDF5Creator(
grid, "F_SPEX_grid.hdf5", instrument=SPEX(), wl_range=(0.9e4, np.inf), ranges=ranges
)
creator.process_grid()
In [4]:
from Starfish.emulator import Emulator
# can load from string or HDF5Interface
emu = Emulator.from_grid("F_SPEX_grid.hdf5")
emu
Out[4]:
In [5]:
%time emu.train(options=dict(maxiter=1e5))
emu
Out[5]:
We can do a sanity check on the optimization by looking at slice of the emulator's parameter space and the corresponding Gaussian process fit. We should see a smooth line connecting all the parameter values with some uncertainty that grows with large gaps or turbulent weights.
In [6]:
%matplotlib inline
from Starfish.emulator.plotting import plot_emulator
plot_emulator(emu)
If we are satisfied, let's save this emulator and move on to fitting some data.
In [7]:
emu.save("F_SPEX_emu.hdf5")