Package ivs :: Package sed :: Module model

Module model

Interface to the SED library.

The most basic usage of this module is:

>>> wave,flux = get_table(teff=10000,logg=4.0)

This will retrieve the model SED with the specified effective temperature and logg, from the standard grid, in standard units and with zero reddening. All these things can be specified though (see below).

Section 1. Available model grids

Section 1.1 Available grids

kurucz: The Kurucz model grids, (default setting) reference: Kurucz 1993, yCat, 6039, 0
- metallicity (z): m01 is -0.1 log metal abundance relative to solar (solar abundances from Anders and Grevesse 1989)
- metallicity (z): p01 is +0.1 log metal abundance relative to solar (solar abundances from Anders and Grevesse 1989)
- alpha enhancement (alpha): True means alpha enhanced (+0.4)
- turbulent velocity (vturb): vturb in km/s
- nover= True means no overshoot
- odfnew=True means no overshoot but with better opacities and abundances
tmap: NLTE grids computed for sdB stars with the Tubingen NLTE Model Atmosphere package. No further parameters are available. Reference: Werner et al. 2003,

Section 1.2 Plotting the domains of all spectral grids

We make a plot of the domains of all spectral grids. Therefore, we first collect all the grid names

>>> grids = get_gridnames()

Preparation of the plot: set the color cycle of the current axes to the spectral color cycle.

>>> p = pl.figure(figsize=(10,8))
>>> color_cycle = [pl.cm.spectral(j) for j in np.linspace(0, 1.0, len(grids))]
>>> p = pl.gca().set_color_cycle(color_cycle)

To plot all the grid points, we run over all grid names (which are strings), and retrieve their dimensions. The dimensions are just two arrays giving the teff- and logg-coordinate of each SED in the grid. They can thus be easily plot:

>>> for grid in grids:
...    teffs,loggs = get_grid_dimensions(grid=grid)
...    p = pl.plot(np.log10(teffs),loggs,'o',ms=7,label=grid)

And we need to set some of the plotting details to make it look nicer.

>>> p = pl.xlim(pl.xlim()[::-1])
>>> p = pl.ylim(pl.ylim()[::-1])
>>> p = pl.xlabel('Effective temperature [K]')
>>> p = pl.ylabel('log( Surface gravity [cm s$^{-1}$]) [dex]')
>>> xticks = [3000,5000,7000,10000,15000,25000,35000,50000,65000]
>>> p = pl.xticks([np.log10(i) for i in xticks],['%d'%(i) for i in xticks])
>>> p = pl.legend(loc='upper left',prop=dict(size='small'))
>>> p = pl.grid()

]include figure]]ivs_sed_model_grid.png]

Section 2. Retrieval of model SEDs

Subsection 2.1 Default settings

To get information on the grid that is currently defined, you can type the following. Note that not all parameters are relevant for all grids, e.g. the convection theory parameter ct has no influence when the Kurucz grid is chosen.

>>> print defaults
{'use_scratch': False, 'Rv': 3.1, 'co': 1.05, 'c': 0.5, 'grid': 'kurucz', 'alpha': False, 'odfnew': True, 'ct': 'mlt', 'a': 0.0, 'vturb': 2, 'law': 'fitzpatrick2004', 'm': 1.0, 't': 1.0, 'z': 0.0, 'nover': False, 'He': 97}

>>> print os.path.basename(get_file())
kurucz93_z0.0_k2odfnew_sed.fits

You can change the defaults with the function set_defaults:

>>> set_defaults(z=0.5)
>>> print defaults
{'use_scratch': False, 'Rv': 3.1, 'co': 1.05, 'c': 0.5, 'grid': 'kurucz', 'alpha': False, 'odfnew': True, 'ct': 'mlt', 'a': 0.0, 'vturb': 2, 'law': 'fitzpatrick2004', 'm': 1.0, 't': 1.0, 'z': 0.5, 'nover': False, 'He': 97}

And reset the 'default' default values by calling set_defaults without arguments

>>> set_defaults()
>>> print defaults
{'use_scratch': False, 'Rv': 3.1, 'co': 1.05, 'c': 0.5, 'grid': 'kurucz', 'alpha': False, 'odfnew': True, 'ct': 'mlt', 'a': 0.0, 'vturb': 2, 'law': 'fitzpatrick2004', 'm': 1.0, 't': 1.0, 'z': 0.0, 'nover': False, 'He': 97}

Subsection 2.2 Speeding up

When fitting an sed using the builder class, or repeatedly reading model seds, or integrated photometry, the main bottleneck on the speed will be the disk access This can be circumvented by using the scratch disk. To do this, call the function copy2scratch() after setting the default settings as explained above. f.x.:

>>> set_defaults(grid='kurucz', z=0.5)
>>> copy2scratch()

You have to do this every time you change a grid setting. This function creates a directory named 'your_username' on the scratch disk and works from there. So you won`t disturbed other users.

After the fitting process use the function

>>> clean_scratch()

to remove the models that you used from the scratch disk. Be carefull with this, because it will remove the models without checking if there is another process using them. So if you have multiple scripts running that are using the same models, only clean the scratch disk after the last process is finished.

The gain in speed can be up to 70% in single sed fitting, and up to 40% in binary and multiple sed fitting.

For the sake of the examples, we'll set the defaults back to z=0.0:

>>> set_defaults()

Subsection 2.3 Model SEDs

Be careful when you supply parameters: e.g., not all grids are calculated for the same range of metallicities. In get_table, only the effective temperature and logg are 'interpolatable' quantities. You have to set the metallicity to a grid value. The reddening value can take any value: it is not interpolated but calculated. You can thus also specify the type of reddening law (see reddening).

>>> wave,flux = get_table(teff=12345,logg=4.321,ebv=0.12345,z=0.5)

but

>>> try:
...     wave,flux = get_table(teff=12345,logg=4.321,ebv=0.12345,z=0.6)
... except IOError,msg:
...     print msg
File sedtables/modelgrids/kurucz93_z0.6_k2odfnew_sed.fits not found in any of the specified data directories /STER/pieterd/IVSDATA/, /STER/kristofs/IVSdata

Since the Kurucz model atmospheres have not been calculated for the value of z=0.6.

Instead of changing the defaults of this module with set_defaults, you can also give extra arguments to get_table to specify the grid you want to use. The default settings will not change in this case.

>>> wave,flux = get_table(teff=16321,logg=4.321,ebv=0.12345,z=0.3,grid='tlusty')

The default units of the SEDs are angstrom and erg/s/cm2/AA/sr. To change them, do:

>>> wave,flux = get_table(teff=16321,logg=4.321,wave_units='micron',flux_units='Jy/sr')

To remove the steradian from the units when you know the angular diameter of your star in milliarcseconds, you can do (we have to convert diameter to surface):

>>> ang_diam = 3.21 # mas
>>> scale = conversions.convert('mas','sr',ang_diam/2.)
>>> wave,flux = get_table(teff=9602,logg=4.1,ebv=0.0,z=0.0,grid='kurucz')
>>> flux *= scale

The example above is representative for the case of Vega. So, if we now calculate the synthetic flux in the GENEVA.V band, we should end up with the zeropoint magnitude of this band, which is close to zero:

>>> flam = synthetic_flux(wave,flux,photbands=['GENEVA.V'])
>>> print '%.3f'%(conversions.convert('erg/s/cm2/AA','mag',flam,photband='GENEVA.V')[0])
0.063

Compare this with the calibrated value

>>> print filters.get_info(['GENEVA.V'])['vegamag'][0]
0.061

Section 3. Retrieval of integrated photometry

Instead of retrieving a model SED, you can immediately retrieve pre-calculated integrated photometry. The benefit of this approach is that it is much faster than retrieving the model SED and then calculating the synthetic flux. Also, you can supply arbitrary metallicities within the grid boundaries, as interpolation is done in effective temperature, surface gravity, reddening and metallicity.

Note that also the reddening law is fixed now, you need to recalculate the tables for different parameters if you need them.

The massive speed-up is accomplished the following way: it may take a few tens of seconds to retrieve the first pre-integrated SED, because all available files from the specified grid will be loaded into memory, and a `markerarray' will be made allowing a binary search in the grid. This makes it easy to retrieve all models around the speficied point in N-dimensional space. Next, a linear interpolation method is applied to predict the photometric values of the specified point.

All defaults set for the retrieval of model SEDs are applicable for the integrated photometry tables as well.

When retrieving integrated photometry, you also get the absolute luminosity (integration of total SED) as a bonus. This is the absolute luminosity assuming the star has a radius of 1Rsol. Multiply by Rstar**2 to get the true luminosity.

Because photometric filters cannot trivially be assigned a wavelength to (see filters.eff_wave), by default, no wavelength information is retrieved. If you want to retrieve the effective wavelengths of the filters themselves (not taking into account the model atmospheres), you can give an extra keyword argument wave_units. If you want to take into account the model atmosphere, use filters.eff_wave.

>>> photbands = ['GENEVA.U','2MASS.J']
>>> fluxes,Labs = get_itable(teff=16321,logg=4.321,ebv=0.12345,z=0.123,photbands=photbands)
>>> waves,fluxes,Labs = get_itable(teff=16321,logg=4.321,ebv=0.12345,z=0.123,photbands=photbands,wave_units='AA')

Note that the integration only gives you fluxes, and is thus independent from the zeropoints of the filters (but dependent on the transmission curves). To get the synthetic magnitudes, you can do

>>> mymags = [conversions.convert('erg/s/cm2/AA','mag',fluxes[i],photband=photbands[i]) for i in range(len(photbands))]

The mags don't mean anything in this case because they have not been corrected for the distance to the star.

The retrieval of integrated photometry can go much faster if you want to do it for a whole set of parameters. The get_itable_pix function has a much more flexible, reliable and fast interpolation scheme. It is possible to interpolate also over doppler shift and interstellar Rv, as long as the grids have been computed before. See get_itable_pix for more information.

Subsection 3. Full example

We build an SED of Vega and compute synthetic magnitudes in the GENEVA and 2MASS bands.

These are the relevant parameters of Vega and photometric passbands

>>> ang_diam = 3.21 # mas
>>> teff = 9602
>>> logg = 4.1
>>> ebv = 0.0
>>> z = 0.0
>>> photbands = ['GENEVA.U','GENEVA.G','2MASS.J','2MASS.H','2MASS.KS']

We can compute (R/d) to scale the synthetic flux as

>>> scale = conversions.convert('mas','sr',ang_diam/2.)

We retrieve the SED

>>> wave,flux = get_table(teff=teff,logg=logg,ebv=ebv,z=z,grid='kurucz')
>>> flux *= scale

Then compute the synthetic fluxes, and compare them with the synthetic fluxes as retrieved from the pre-calculated tables

>>> fluxes_calc = synthetic_flux(wave,flux,photbands)
>>> wave_int,fluxes_int,Labs = get_itable(teff=teff,logg=logg,ebv=ebv,z=z,photbands=photbands,wave_units='AA')
>>> fluxes_int *= scale

Convert to magnitudes:

>>> m1 = [conversions.convert('erg/s/cm2/AA','mag',fluxes_calc[i],photband=photbands[i]) for i in range(len(photbands))]
>>> m2 = [conversions.convert('erg/s/cm2/AA','mag',fluxes_int[i],photband=photbands[i]) for i in range(len(photbands))]

And make a nice plot

>>> p = pl.figure()
>>> p = pl.loglog(wave,flux,'k-',label='Kurucz model')
>>> p = pl.plot(wave_int,fluxes_calc,'ro',label='Calculated')
>>> p = pl.plot(wave_int,fluxes_int,'bx',ms=10,mew=2,label='Pre-calculated')
>>> p = [pl.annotate('%s: %.3f'%(b,m),(w,f),color='r') for b,m,w,f in zip(photbands,m1,wave_int,fluxes_calc)]
>>> p = [pl.annotate('%s: %.3f'%(b,m),(w-1000,0.8*f),color='b') for b,m,w,f in zip(photbands,m2,wave_int,fluxes_int)]
>>> p = pl.xlabel('Wavelength [Angstrom]')
>>> p = pl.ylabel('Flux [erg/s/cm2/AA]')

]include figure]]ivs_sed_model_example.png]

Functions

[hide private]

_get_itable_markers(photbands, teffrange=(-np.inf,np.inf), loggrange=(-np.inf,np.inf), ebvrange=(-np.inf,np.inf), zrange=(-np.inf,np.inf), include_Labs=True, clear_memory=True, **kwargs)
Get a list of markers to more easily retrieve integrated fluxes.

source code

_get_pix_grid(photbands, teffrange=(-np.inf,np.inf), loggrange=(-np.inf,np.inf), ebvrange=(-np.inf,np.inf), zrange=(-np.inf,np.inf), rvrange=(-np.inf,np.inf), vradrange=(-np.inf,np.inf), include_Labs=True, clear_memory=True, variables=['teff','logg','ebv','z','rv','vrad'], **kwargs)
Prepare the pixalted grid.

source code

_get_flux_from_table(fits_ext, photbands, index=None, include_Labs=True)
Retrieve flux and flux ratios from an integrated SED table.

source code

Interface to library

set_defaults(*args, **kwargs)
Set defaults of this module

source code

set_defaults_multiple(*args)
Set defaults for multiple stars

source code

copy2scratch(**kwargs)
Copy the grids to the scratch directory to speed up the fitting process.

source code

clean_scratch(**kwargs)
Remove the grids that were copied to the scratch directory by using the function copy2scratch().

source code

defaults2str()
Convert the defaults to a string, e.g.

source code

defaults_multiple2str()
Convert the defaults to a string, e.g.

source code

list of str

get_gridnames(grid=None)
Return a list of available grid names.

source code

str

get_file(integrated=False, **kwargs)
Retrieve the filename containing the specified SED grid.

source code

_blackbody_input(fctn)
Prepare input and output for blackbody-like functions.

source code

array

blackbody(x, T, wave_units='AA', flux_units='erg/s/cm2/AA', disc_integrated=True, ang_diam=None)
Definition of black body curve.

source code

array

rayleigh_jeans(x, T, wave_units='AA', flux_units='erg/s/cm2/AA', disc_integrated=True, ang_diam=None)
Rayleigh-Jeans approximation of a black body.

source code

array

wien(x, T, wave_units='AA', flux_units='erg/s/cm2/AA', disc_integrated=True, ang_diam=None)
Wien approximation of a black body.

source code

(ndarray,ndarray)

get_table_single(teff=None, logg=None, ebv=None, rad=None, star=None, wave_units='AA', flux_units='erg/s/cm2/AA/sr', **kwargs)
Retrieve the spectral energy distribution of a model atmosphere.

source code

(ndarray,)ndarray,float

get_itable_single(teff=None, logg=None, ebv=0, z=0, rad=None, photbands=None, wave_units=None, flux_units='erg/s/cm2/AA/sr', **kwargs)
Retrieve the spectral energy distribution of a model atmosphere in photometric passbands.

source code

(ndarray,ndarray)

get_itable(photbands=None, wave_units=None, flux_units='erg/s/cm2/AA/sr', grids=None, **kwargs)
Retrieve the integrated spectral energy distribution of a combined model atmosphere.

source code

get_itable_single_pix(teff=None, logg=None, ebv=None, z=0, rv=3.1, vrad=0, photbands=None, wave_units=None, flux_units='erg/s/cm2/AA/sr', **kwargs)
Super fast grid interpolator.

source code

get_itable_pix(photbands=None, wave_units=None, flux_units='erg/s/cm2/AA/sr', grids=None, **kwargs)
Super fast grid interpolator for multiple tables, completely based on get_itable_pix.

source code

(ndarray,ndarray)

get_table(wave_units='AA', flux_units='erg/cm2/s/AA/sr', grids=None, full_output=False, **kwargs)
Retrieve the spectral energy distribution of a combined model atmosphere.

source code

(ndarray,ndarray)

get_grid_dimensions(**kwargs)
Retrieve possible effective temperatures and gravities from a grid.

source code

(ndarray,ndarray,ndarray)

get_igrid_dimensions(**kwargs)
Retrieve possible effective temperatures, surface gravities and reddenings from an integrated grid.

source code

(3x1Darray,3Darray,interp_func)

get_grid_mesh(wave=None, teffrange=None, loggrange=None, **kwargs)
Return InterpolatingFunction spanning the available grid of atmosphere models.

source code

Calibration

list of str

list_calibrators(library='calspec')
Print and return the list of calibrators

source code

(ndarray,ndarray)

get_calibrator(name='alpha_lyr', version=None, wave_units=None, flux_units=None, library='calspec')
Retrieve a calibration SED

source code

read_calibrator_info(library='ngsl')

source code

calibrate()
Calibrate photometry.

source code

Synthetic photometry

ndarray

synthetic_flux(wave, flux, photbands, units=None)
Extract flux measurements from a synthetic SED (Fnu or Flambda).

source code

ndarray

synthetic_color(wave, flux, colors, units=None)
Construct colors from a synthetic SED.

source code

float

luminosity(wave, flux, radius=1.)
Calculate the bolometric luminosity of a model SED.

source code

Variables

[hide private]

new_scipy = False

logger = logging.getLogger("SED.MODEL")

caldir = os.sep.join(['sedtables', 'calibrators'])

__defaults__ = dict(grid= 'kurucz', odfnew= True, z=+ 0.0, vtu...

defaults = __defaults__.copy()

defaults_multiple = [defaults.copy(), defaults.copy()]

basedir = 'sedtables/modelgrids/'

scratchdir = None
hash(x)

Function Details