Spectrum

The Spectrum class handles any 1D data representing the light from an astronomical source. A Spectrum object is created by passing the wavelength, flux density, and (optional) uncertainty with astropy.units to the class.

import astropy.units as q
import numpy as np
from sedkit import Spectrum
wavelength = np.linspace(1, 2, 100) * q.um
flux = np.ones(100) * q.erg / q.s / q.cm**2 / q.AA
unc = flux / 100.
spec = Spectrum(wavelength, flux, unc, name='My spectrum')

The Spectrum has a number of useful attributes.

spec.name           # The object's name
spec.ref            # Reference(s)
spec.header         # A header for the data
spec.wave           # Unitless wavelength array
spec.flux           # Unitless flux density array
spec.unc            # Unitless uncertainty array (optional)
spec.data           # Unitless array of [W, F] or [W, F, E]
spec.spectrum       # List of [W, F] or [W, F, E] with units
spec.wave_max       # The max wavelength
spec.wave_min       # The min wavelength
spec.wave_units     # The wavelength units
spec.flux_units     # The flux density units
spec.size           # The number of data points

After the Spectrum has been created, it be manipulated in a number of ways.

It can be trimmed by passing a list of lower and upper bounds to trim method. The include argument accepts bounds for wavelength regions to include and the exclude argument accepts bounds for regions to exclude. A list of Spectrum objects are returned unless the concat argument is set to True, which simply concatenated the trimmed segment(s) into one Spectrum.

trim_spec_include = spec.trim(include=[(1.2 * q.um, 1.6 * q.um)])
trim_spec_exclude = spec.trim(exclude=[(1.2 * q.um, 1.6 * q.um)], concat=True)

The interpolate method accepts a new wavelength array and returns a new Spectrum object interpolated to those values. The resamp method accepts the same input and resamples the spectrum onto the new wavelength array while preserving the total flux.

new_wav = np.linspace(1.2, 1.6, 50) * q.um
new_wav_interp = spec.interpolate(new_wav)
new_wav_resamp = spec.resamp(new_wav)

The integrate method integrates the curve to caluclate the area underneath using the trapezoidal rule.

area = spec.integrate()

The Spectrum can be smoothed using a Kaiser-Bessel smoothing window of narrowness beta and a given window size.

smooth_spec = spec.smooth(beta=2, window=11)

A Spectrum may be flux calibrated to a given distance by passing a distance to the flux_calibrate method.

cal_spec = spec.flux_calibrate(5.1 * q.pc, 10 * q.pc)   # Flux calibrates the spectrum from 5.1 to 10 pc

A bandpass name or svo_filters.svo.Filter object can be used to convolve the spectrum, calculate a synthetic flux or magnitude, or renormalize it to a given magnitude.

from svo_filters import Filter
jband = Filter('2MASS.J')
conv_spec = spec.convolve_filter(jband)     # Convolved spectrum
norm_spec = spec.renormalize(12.3, jband)   # Renormalized spectrum
jmag = spec.synthetic_magnitude(jband)      # Synthetic magnitude
jflux = spec.synthetic_flux(jband)          # Synthetic flux

It can also be normalized to a table of photometry weighted by the magnitude uncertainties with the norm_to_mags method. See the SED class for an example.

A Spectrum object may also interact with another Spectrum object in a number of ways. The norm_to_spec method creates a new object normalized to the input Spectrum object and the __add__ operation combines two Spectrum objects in their common wavelength region or concatenates the segments.

spec2 = Spectrum(np.linspace(1.5, 2.5, 100) * q.um, flux * 1E-11, unc * 1E-11, name='Redder spectrum')
normed_spec = spec.norm_to_spec(spec2)  # spec normalized to spec2
combined_spec = spec + spec2            # New combined spectrum

Any Spectrum may also be fit by a ModelGrid object to find the best fit model or spectrum. The best_fit_model method performs a simple goodness of fit test and returns the model with the best fit. The mcmc_fit` method performs a MCMC fit to the grid and returns the best fit parameters with uncertainties. The details of any fit are stored as a dictionary in the ``best_fit attribute.

from sedkit import BTSettl
bt_grid = BTSettl()
spec.best_fit_model(bt)             # Goodness of fit
spec.mcmc_fit(bt, params=['teff'])  # MCMC fit

For visual inspection of an interactive bokeh.plotting.figure, use the plot method.

spec.plot()

Several child classes are also available. The Blackbody class creates a blackbody spectrum for a given wavelength range and effective temperature.

from sedkit.spectrum import Blackbody
bb_spec = Blackbody(np.linspace(1.5, 2.5, 100) * q.um, teff=2456 * q.K)

The ever useful spectrum of Vega is easily created.

from sedkit.spectrum import Vega
vega = Vega()

And a Spectrum object can be created directly from a FITS or ASCII file with the FileSpectrum class.

from sedkit.spectrum import FileSpectrum
file_spec = FileSpectrum('/path/to/the/file.fits', wave_units=q.um, flux_units=q.erg/q.s/q.cm**2/q.AA)

Finally, the data can be exported to an ASCII file by passing a filepath to the export method.

spec.export('path/to/the/new/file.txt')