spectrum¶
The spectrum module consists of the Spectrum class, with child classes ObsBlock and Spectra for multi-spectrum analysis of different types.
The Spectrum class is the main functional code. ObsBlocks are containers of multiple spectra of different objects The Spectra class is a container of multiple spectra of the same object at different wavelengths/frequencies
- class pyspeckit.spectrum.classes.Spectrum(filename=None, filetype=None, xarr=None, data=None, error=None, header=None, doplot=False, maskdata=True, unit=None, plotkwargs={}, xarrkwargs={}, **kwargs) [github] [bitbucket]¶
Bases: object
The core class for the spectroscopic toolkit. Contains the data and error arrays along with wavelength / frequency / velocity information in various formats.
Create a Spectrum object.
Must either pass in a filename or ALL of xarr, data, and header, plus optionally error.
kwargs are passed to the file reader
Parameters: filename : string
The file to read the spectrum from. If data, xarr, and error are specified, leave filename blank.
filetype : string
Specify the file type (only needed if it cannot be automatically determined from the filename)
xarr : units.SpectroscopicAxis or np.ndarray
The X-axis of the data. If it is an np.ndarray, you must pass xarrkwargs or a valid header if you want to use any of the unit functionality.
data : np.ndarray
The data array (must have same length as xarr)
error : np.ndarray
The error array (must have same length as the data and xarr arrays)
header : pyfits.Header or dict
The header from which to read unit information. Needs to be a pyfits.Header instance or another dictionary-like object with the appropriate information
maskdata : boolean
turn the array into a masked array with all nan and inf values masked
doplot : boolean
Plot the spectrum after loading it?
plotkwargs : dict
keyword arguments to pass to the plotter
xarrkwargs : dict
keyword arguments to pass to the SpectroscopicAxis initialization (can be used in place of a header)
unit : str
The data unit
Examples
>>> sp = pyspeckit.Spectrum(data=np.random.randn(100), xarr=np.linspace(-50, 50, 100), error=np.ones(100)*0.1, xarrkwargs={'unit':'km/s', 'refX':4.829, 'refX_unit':'GHz', 'xtype':'VLSR-RAD'}, header={})
>>> xarr = pyspeckit.units.SpectroscopicAxis(np.linspace(-50,50,100), units='km/s', refX=6562.83, refX_unit='angstroms') >>> data = np.random.randn(100)*5 + np.random.rand(100)*100 >>> err = np.sqrt(data/5.)*5. # Poisson noise >>> sp = pyspeckit.Spectrum(data=data, error=err, xarr=xarr, header={})
>>> # if you already have a simple fits file >>> sp = pyspeckit.Spectrum('test.fits')
- copy(deep=True) [github] [bitbucket]¶
Create a copy of the spectrum with its own plotter, fitter, etc. Useful for, e.g., comparing smoothed to unsmoothed data
- crop(x1, x2, unit=None, **kwargs) [github] [bitbucket]¶
Replace the current spectrum with a subset from x1 to x2 in current units
Fixes CRPIX1 and baseline and model spectra to match cropped data spectrum
- downsample(dsfactor) [github] [bitbucket]¶
Downsample the spectrum (and all of its subsidiaries) without smoothing
Parameters: dsfactor : int
Downsampling Factor
- flux¶
The data in the spectrum (flux = data, for compatibility with astropy’s Spectrum1D object).
- classmethod from_hdu(hdu) [github] [bitbucket]¶
Create a pyspeckit Spectrum object from an HDU
- classmethod from_spectrum1d(spec1d) [github] [bitbucket]¶
Tool to load a pyspeckit Spectrum from a specutils object
(this is intended to be temporary; long-term the pyspeckit Spectrum object will inherit from a specutils Spectrum1D object)
- getlines(linetype='radio', **kwargs) [github] [bitbucket]¶
Access a registered database of spectral lines. Will add an attribute with the name linetype, which then has properties defined by the speclines module (most likely, a table and a “show” function to display the lines)
- interpnans(spec) [github] [bitbucket]¶
Interpolate over NAN values, replacing them with values interpolated from their neighbors using linear interpolation.
- measure(z=None, d=None, fluxnorm=None, miscline=None, misctol=10.0, ignore=None, derive=True, **kwargs) [github] [bitbucket]¶
Initialize the measurements class - only do this after you have run a fitter otherwise pyspeckit will be angry!
- moments(unit='km/s', **kwargs) [github] [bitbucket]¶
Return the moments of the spectrum. In order to assure that the 1st and 2nd moments are meaningful, a ‘default’ unit is set. If unit is not set, will use current unit.
Documentation imported from the moments module:
Returns the gaussian parameters of a 1D distribution by calculating its moments. Depending on the input parameters, will only output a subset of the above.
Theory, from first principles (in the absence of noise): integral(gaussian) = sqrt(2*pi*sigma^2) * amp sigma = integral / amp / sqrt(2*pi)
In the presence of noise, this gets much more complicated. The noisy approach is inspired by mpfit
Parameters: height: :
is the background level
amplitude: :
is the maximum (or minimum) of the data after background subtraction
x: :
is the first moment
width_x: :
is the second moment
estimator: function :
A function to estimate the “height” or “background level” of the data, e.g. mean or median. If masked arrays are being used, use the np.ma versions of the numpy functions
negamp: bool or None :
Force the peak negative (True), positive (False), or the sign of the peak will be “autodetected” (negamp=None)
nsigcut: float or None :
If specified, the code will attempt to estimate the noise and only use data above/below n-sigma above the noise. The noise will be estimated from the data unless the noise is specified with noise_estimate
noise_estimate: float or None :
Guess for the noise value. Only matters if nsigcut is specified.
Returns: (height, amplitude, x, width_x) :
- parse_hdf5_header(hdr) [github] [bitbucket]¶
HDF5 reader will create a hdr dictionary from HDF5 dataset attributes if they exist. This routine will convert that dict to a pyfits header instance.
- parse_header(hdr, specname=None) [github] [bitbucket]¶
Parse parameters from a .fits header into required spectrum structure parameters
- parse_text_header(Table) [github] [bitbucket]¶
Grab relevant parameters from a table header (xaxis type, etc)
This function should only exist for Spectrum objects created from .txt or other atpy table type objects
- shape¶
Return the data shape (a property of the Spectrum)
- slice(start=None, stop=None, unit='pixel', copy=True, preserve_fits=False) [github] [bitbucket]¶
Slicing the spectrum
Warning
this is the same as cropping right now, but it returns a copy instead of cropping inplace
Parameters: start : numpy.float or int
start of slice
stop : numpy.float or int
stop of slice
unit : str
allowed values are any supported physical unit, ‘pixel’
copy : bool
Return a ‘view’ of the data or a copy?
preserve_fits : bool
Save the fitted parameters from self.fitter?
- smooth(smooth, downsample=True, **kwargs) [github] [bitbucket]¶
Smooth the spectrum by factor smooth.
Documentation from the smooth module:
Parameters: smooth : float
Number of pixels to smooth by
smoothtype : [ ‘gaussian’,’hanning’, or ‘boxcar’ ]
type of smoothing kernel to use
downsample : bool
Downsample the data?
downsample_factor : int
Downsample by the smoothing factor, or something else?
convmode : [ ‘full’,’valid’,’same’ ]
see numpy.convolve. ‘same’ returns an array of the same length as ‘data’ (assuming data is larger than the kernel)
- stats(statrange=(), interactive=False) [github] [bitbucket]¶
Return some statistical measures in a dictionary (somewhat self-explanatory)
Parameters: statrange : 2-element tuple
X-range over which to perform measures
interactive : bool
specify range interactively in plotter
- unit¶
- units¶
- write(filename, type=None, **kwargs) [github] [bitbucket]¶
Write the spectrum to a file. The available file types are listed in spectrum.writers.writers
type - what type of file to write to? If not specified, will attempt to determine type from suffix
- class pyspeckit.spectrum.classes.Spectra(speclist, xunit='GHz', **kwargs) [github] [bitbucket]¶
Bases: pyspeckit.spectrum.classes.Spectrum
A list of individual Spectrum objects. Intended to be used for concatenating different wavelength observations of the SAME OBJECT. Can be operated on just like any Spectrum object, incuding fitting. Useful for fitting multiple lines on non-continguous axes simultaneously. Be wary of plotting these though...
Can be indexed like python lists.
X array is forcibly sorted in increasing order
- fiteach(**kwargs) [github] [bitbucket]¶
Fit each spectrum within the Spectra object
- ploteach(xunit=None, inherit_fit=False, plot_fit=True, plotfitkwargs={}, **plotkwargs) [github] [bitbucket]¶
Plot each spectrum in its own window inherit_fit - if specified, will grab the fitter & fitter properties from Spectra
- smooth(smooth, **kwargs) [github] [bitbucket]¶
Smooth the spectrum by factor “smooth”. Options are defined in sm.smooth
because ‘Spectra’ does not have a header attribute, don’t do anything to it...
- class pyspeckit.spectrum.classes.ObsBlock(speclist, xtype='frequency', xarr=None, force=False, **kwargs) [github] [bitbucket]¶
Bases: pyspeckit.spectrum.classes.Spectra
An Observation Block
Consists of multiple spectra with a shared X-axis. Intended to hold groups of observations of the same object in the same setup for later averaging.
ObsBlocks can be indexed like python lists.
- average(weight=None, inverse_weight=False, error='erravgrtn', debug=False) [github] [bitbucket]¶
Average all scans in an ObsBlock. Returns a single Spectrum object
Parameters: weight : string
a header keyword to weight by. If not specified, the spectra will be averaged without weighting
inverse_weight : bool
Is the header keyword an inverse-weight (e.g., a variance?)
error : [‘scanrms’,’erravg’,’erravgrtn’]
estimate the error spectrum by one of three methods. ‘scanrms’ : the standard deviation of each pixel across all scans ‘erravg’ : the average of all input error spectra ‘erravgrtn’ : the average of all input error spectra divided by sqrt(n_obs)
- smooth(smooth, **kwargs) [github] [bitbucket]¶
Smooth the spectrum by factor “smooth”. Options are defined in sm.smooth
- pyspeckit.spectrum.register_reader(filetype, function, suffix, default=False) [github] [bitbucket]¶
Register a reader function.
Parameters: filetype: str :
The file type name
function: function :
The reader function. Should take a filename as input and return an X-axis object (see units.py), a spectrum, an error spectrum (initialize it to 0’s if empty), and a pyfits header instance
suffix: int :
What suffix should the file have?
- pyspeckit.spectrum.register_writer(filetype, function, suffix, default=False) [github] [bitbucket]¶
Register a writer function.
Parameters: filetype: string :
The file type name
function: function :
The writer function. Will be an attribute of Spectrum object, and called as spectrum.Spectrum.write_hdf5(), for example.
suffix: int :
What suffix should the file have?