specter API

specter

A toolkit for simulating multi-object spectrographs.

The specter toolkit is centered around a generic PSF object, which describes how a given fiber (spectrum) at a particular wavelength is projected onto the detector CCD. The base class specter.psf.psf.PSF defines the interface for any PSF object. Specific instruments (e.g. BOSS or DESI) implement a subclass which provides:

xslice, yslice, pixels[ny,nx] = PSF.xypix(ispec, wavelength)

to describe how a given spectrum ispec and wavelength project onto the CCD.

specter.extract

Tools for extracting spectra from 2D images.

Notes

Initial work for adding extractions. Starting with row-by-row. This may be renamed and/or refactored into classes:

from specter.psf import load_psf
from specter.extract.ex1d import extract1d

psf = load_psf('data/boss/pixpsf-r1-00140299.fits')
img = fitsio.read('data/boss/img-r1-00140299.fits', 0)
ivar = fitsio.read('data/boss/img-r1-00140299.fits', 1)

specflux, specivar = extract1d(img, ivar, psf, specrange=(20,60), yrange=(1000,1100))

%time specflux, specivar = extract1d(img, ivar, psf, yrange=(1000,1100), specrange=(20,40) )

rows	time
10	0.32
20	0.59
50	1.75
100	4.94
200	30.7

0.59 * 25 * (4000/100) / 60 = 10 minutes. Not bad. Up to 15 minutes now after refactoring xsigma.

from specter.psf import load_psf
psf = load_psf('data/boss/pixpsf-r1-00140299.fits')
psf.xsigma(0, 7000)

### yy = N.linspace(10,4000)
yy = N.arange(10,4000,5)
ww = psf.wavelength(0, y=yy)
xsig = [psf.xsigma(0, wavelength=w) for w in ww]

plot(yy, xsig)

from specter.util import sincshift
from scipy.ndimage import center_of_mass
def xsigma(spot):
    yc, xc = center_of_mass(spot)
    xx = N.arange(spot.shape[1])
    xspot = spot.sum(axis=0)
    return N.sqrt(N.sum(xspot*(xx-xc)**2) / N.sum(xspot))

xr, yr, spot = psf.xypix(0, psf.wavelength(0, y=1000))

xx = N.arange(xr.start, xr.stop)
yy = N.arange(yr.start, yr.stop)
xspot = spot.sum(axis=0)
yspot = spot.sum(axis=1)

N.sum(xx*xspot) / N.sum(xspot)
N.sum(yy*yspot) / N.sum(yspot)

#- Variance
from scipy.ndimage import center_of_mass
yc, xc = center_of_mass(spot)
xx = N.arange(spot.shape[1])
yy = N.arange(spot.shape[0])
N.sum(xspot*(xx-xc)**2) / N.sum(xspot)
N.sum(yspot*(yy-yc)**2) / N.sum(yspot)

specter.extract.ex1d

1D Extraction like Horne 1986

Stephen Bailey, LBL Spring 2013

specter.extract.ex1d.ex1d(img, mask, psf, readnoise=2.5, specrange=None, yrange=None, nspec_per_group=20, debug=False, model=False)

Extract spectra from an image using row-by-row weighted extraction.

Parameters:

• img (array-like) – CCD image.

• mask (array-like) – Image mask, 0=good, non-zero=bad

• psf (object) – PSF object.

• readnoise (float, optional) – CCD readnoise

• specrange (list, optional) – (specmin, specmax) Spectral range to extract (default all)

• yrange (list, optional) – (ymin, ymax) CCD y (row) range to extract (default all rows) ranges are python-like, i.e. yrange=(0,100) extracts 100 rows from 0 to 99 inclusive but not row 100.

• nspec_per_group (int, optional) – Extract spectra in groups of N spectra (faster if spectra are physically separated into non-overlapping groups).

• debug (bool, optional) – If True, stop with prompt after each row

• model (bool, optional) – Unknown parameter

Returns:

Tuple containing spectra[nspec, ny] the extracted spectra, specivar[nspec, ny] the inverse variance of spectra.

Return type:

tuple

specter.extract.ex2d

2D Spectroperfectionism extractions

specter.extract.ex2d.eigen_compose(w, v, invert=False, sqr=False)

Create a matrix from its eigenvectors and eigenvalues.

Given the eigendecomposition of a matrix, recompose this into a real symmetric matrix. Optionally take the square root of the eigenvalues and / or invert the eigenvalues. The eigenvalues are regularized such that the condition number remains within machine precision for 64bit floating point values.

Parameters:

• w (array) – 1D array of eigenvalues

• v (array) – 2D array of eigenvectors.

• invert (bool) – Should the eigenvalues be inverted? (False)

• sqr (bool) – Should the square root eigenvalues be used? (False)

Returns:

A 2D numpy array which is the recomposed matrix.

specter.extract.ex2d.ex2d(image, imageivar, psf, specmin, nspec, wavelengths, xyrange=None, regularize=0.0, ndecorr=False, bundlesize=25, nsubbundles=1, wavesize=50, full_output=False, verbose=False, debug=False, psferr=None, pixpad_frac=0.8, wavepad_frac=0.2)

2D PSF extraction of flux from image patch given pixel inverse variance.

Parameters:

• image (array-like) – 2D array of pixels

• imageivar (array-like) – 2D array of inverse variance for the image

• psf (object) – PSF object

• specmin (int) – index of first spectrum to extract

• nspec (int) – number of spectra to extract

• wavelengths (array-like) – 1D array of wavelengths to extract

• xyrange (list, optional) – (xmin, xmax, ymin, ymax): treat image as a subimage cutout of this region from the full image

• regularize (float, optional) – experimental regularization factor to minimize ringing

• ndecorr (bool, optional) – if True, decorrelate the noise between fibers, at the cost of residual signal correlations between fibers.

• bundlesize (int, optional) – extract in groups of fibers of this size, assuming no correlation with fibers outside of this bundle

• nsubbundles (int, optional) – number of overlapping subbundles to use per bundle

• wavesize (int, optional) – number of wavelength steps to include per sub-extraction

• full_output (bool, optional) – Include additional outputs based upon chi2 of model projected into pixels

• verbose (bool, optional) – print more stuff

• debug (bool, optional) – if True, enter interactive ipython session before returning

• psferr (float, optional) – fractional error on the psf model. if not None, use this fractional error on the psf model instead of the value saved in the psf fits file. This is used only to compute the chi2, not to weight pixels in fit

• pixpad_frac (float, optional) – fraction of a PSF spotsize to pad in pixels when extracting

• wavepad_frac (float, optional) – fraction of a PSF spotsize to pad in wavelengths when extracting

Returns:

A tuple of (flux, ivar, Rdata):

• flux[nspec, nwave] = extracted resolution convolved flux

• ivar[nspec, nwave] = inverse variance of flux

• Rdata[nspec, 2*ndiag+1, nwave] = sparse Resolution matrix data

Return type:

tuple

Notes

• TODO: document output if full_output=True

• ex2d uses divide-and-conquer to extract many overlapping subregions and then stitches them back together. Params wavesize and bundlesize control the size of the subregions that are extracted; the necessary amount of overlap is auto-calculated based on PSF extent.

specter.extract.ex2d.ex2d_patch(image, ivar, psf, specmin, nspec, wavelengths, xyrange=None, full_output=False, regularize=0.0, ndecorr=False, use_cache=None)

2D PSF extraction of flux from image patch given pixel inverse variance.

Parameters:

• image – 2D array of pixels

• ivar – 2D array of inverse variance for the image

• psf – PSF object

• specmin – index of first spectrum to extract

• nspec – number of spectra to extract

• wavelengths – 1D array of wavelengths to extract

• xyrange – (xmin, xmax, ymin, ymax): treat image as a subimage cutout of this region from the full image

• full_output – if True, return a dictionary of outputs including intermediate outputs such as the projection matrix.

• ndecorr – if True, decorrelate the noise between fibers, at the cost of residual signal correlations between fibers.

• use_cache – default behavior, can be turned off for testing purposes

Returns (flux, ivar, R):

flux[nspec, nwave] = extracted resolution convolved flux ivar[nspec, nwave] = inverse variance of flux R : 2D resolution matrix to convert

specter.extract.ex2d.psfabsbias(p1, p2, wave, phot, ispec=0, readnoise=3.0)

Return absolute bias from extracting with PSF p2 if the real PSF is p1.

Inputs:

p1, p2 : PSF objects wave[] : wavelengths in Angstroms phot[] : spectrum in photons

Optional Inputs:

ispec : spectrum number readnoise : CCD read out noise (optional)

Returns bias, R

bias array same length as wave R resolution matrix for PSF p1

See psfbias() for relative bias

specter.extract.ex2d.psfbias(p1, p2, wave, phot, ispec=0, readnoise=3.0)

Return bias from extracting with PSF p2 if the real PSF is p1

Inputs:

p1, p2 : PSF objects wave[] : wavelengths in Angstroms phot[] : spectrum in photons

Optional Inputs:

ispec : spectrum number readnoise : CCD read out noise (optional)

Returns:

bias array same length as wave

specter.extract.ex2d.resolution_from_icov(icov, decorr=None)

Function to generate the ‘resolution matrix’ in the simplest (no unrelated crosstalk) Bolton & Schlegel 2010 sense. Works on dense matrices. May not be suited for production-scale determination in a spectro extraction pipeline.

Parameters:

• icov (array) – real, symmetric, 2D array containing inverse covariance.

• decorr (list) – produce a resolution matrix which decorrelates signal between fibers, at the cost of correlated noise between fibers (default). This list should contain the number of elements in each spectrum, which is used to define the size of the blocks.

Returns (R, ivar):

R : resolution matrix ivar : R C R.T – decorrelated resolution convolved inverse variance

specter.extract.ex2d.split_bundle(bundlesize, n)

Partitions a bundle into subbundles for extraction

Parameters:

• bundlesize – (int) number of fibers in the bundle

• n – (int) number of subbundles to generate

Returns:

(subbundles, extract_subbundles) where

subbundles = list of arrays of indices belonging to each subbundle; extract_subbundles = list of arrays of indices to extract for each subbundle, including edge overlaps except for first and last fiber

NOTE: resulting partition is such that the lengths of the extract_subbundles differ by at most 1.

Example:

>>> split_bundle(10, 3)
([array([0, 1, 2]), array([3, 4, 5]), array([6, 7, 8, 9])],
 [array([0, 1, 2, 3]), array([2, 3, 4, 5, 6]), array([5, 6, 7, 8, 9])])

specter.io

Routines to standarize specter I/O. This may be refactored into a separate directory as this grows.

Stephen Bailey, LBL January 2013

specter.io.read_image(filename)

Read (image, ivar, header) from input filename

specter.io.read_simspec(filename)

Read an input simulated spectrum file, parse the various format options, and return a standardized spectrum dictionary for use.

Returns a dictionary with keys flux, wavelength, units, objtype

specter.io.read_simspec_image(filename)

Read an input simulated spectrum file formatted in multi-HDU FITS images.

Returns a dictionary with keys flux, wavelength, units, objtype

specter.io.read_simspec_table(filename)

Read an input simulated spectrum file formatted as a FITS binary table.

Returns a dictionary with keys flux, wavelength, units, objtype

specter.io.write_spectra(outfile, wave, flux, ivar, resolution, header)

Write spectra to outfile

Parameters:

• wave – 1D[nwave] array of wavelengths

• flux – 2D[nspec, nwave] flux

• ivar – 2D[nspec, nwave] inverse variance of flux

• resolution – 3D[nspec, ndiag, nwave] sparse resolution matrix data

• header – fits header to include in output

Writes data to outfile in 4 HDUs with EXTNAME FLUX, IVAR, WAVE, RESOLUTION

specter.psf

specter.psf.load_psf(filename, psftype=None)

Load a PSF fits file, using the “PSFTYPE” header keyword to determine which specific subclass to use.

psftype : override fits header keyword PSFTYPE

specter.psf.gausshermite

GaussHermitePSF - PSF modeled with 2D Gauss-Hermite polynomials as generated by the specex package at https://github.com/julienguy/specex

Stephen Bailey December 2013

class specter.psf.gausshermite.GaussHermitePSF(filename)

Model PSF with two central Gauss-Hermite cores with different sigmas plus power law wings.

_gh(x, m=0, xc=0.0, sigma=1.0)

return Gauss-Hermite function value, NOT integrated, for display of PSF. :param x: coordinates baseline array :param m: order of Hermite polynomial multiplying Gaussian core :param xc: sub-pixel position of Gaussian centroid relative to x baseline :param sigma: sigma parameter of Gaussian core in units of pixels

_value(x, y, ispec, wavelength)

return PSF 2D array (shape from x and y), NOT integrated, for display of PSF.

Parameters:

• x – x-coordinates baseline array (1D or 2D)

• y – y-coordinates baseline array (1D or 2D; if 2D same shape as x)

• ispec – fiber

• wavelength – wavelength

Returns: 2D psf

_xypix(ispec, wavelength, ispec_cache=None, iwave_cache=None)

Two branches of this function which does legendre series fitting First branch = no caching, will recompute legval every time eval is used (slow) Second branch = yes caching, will look up precomputed values of legval (faster)

Parameters:

• ispec – the index of each spectrum

• wavelength – the wavelength at which to evaluate the legendre series

Optional arguments:

ispec_cache: the index of each spectrum which starts again at 0 for each patch iwave_cache: the index of each wavelength which starts again at 0 for each patch

cache_params(specrange, wavelengths)

this is implemented in specter.psf.gausshermite, everywhere else just an empty function

xsigma(ispec, wavelength)

Return Gaussian sigma of PSF spot in cross-dispersion direction in CCD pixel units.

ispec : spectrum index wavelength : scalar or vector wavelength(s) to evaluate spot sigmas

ysigma(ispec, wavelength)

Return Gaussian sigma of PSF spot in dispersion direction in CCD pixel units.

ispec : spectrum index wavelength : scalar or vector wavelength(s) to evaluate spot sigmas

specter.psf.gausshermite.generate_core(degx1, degy1, npfx, npfy, spot1, core1, c1_array)

Funtion to speed up some of the expensive parts of _xypix by using numba to JIT-compile.

Parameters:

• degx1 – self._polyparams[‘GHDEGX’], gauss-hermite x degree

• degy1 – self._polyparams[‘GHDEGY’], gauss-hermite y degree

• npfx – xfunc1 converted to numpy array (numba doesn’t like lists of numpy arrays)

• npfy – yfunc1 converted to numpy array (numba doesn’t like lists of numpy arrays)

• spot1 – a preallocated empty 2d array that is the same size as core1

• core1 – a 2d array that is modified in place and then returned as the final function output

• c1_array – the legval values for degx1, degy1, ispec_cache, and iwave_cache which cannot be fed in directly because numba will not handle dictionaries

specter.psf.gausshermite.pgh(x, m=0, xc=0.0, sigma=1.0)

Pixel-integrated (probabilist) Gauss-Hermite function. :param x: pixel-center baseline array :param m: order of Hermite polynomial multiplying Gaussian core :param xc: sub-pixel position of Gaussian centroid relative to x baseline :param sigma: sigma parameter of Gaussian core in units of pixels

Uses the relationship Integral{ H_k(x) exp(-0.5 x^2) dx} = -H_{k-1}(x) exp(-0.5 x^2) + const Written: Adam S. Bolton, U. of Utah, fall 2010 Adapted for efficiency by S. Bailey while dropping generality modified from the orig _pgh to enable jit-compiling –> no longer passing in self, calling custom numba functions in util

specter.psf.gausshermite2

GaussHermitePSF - PSF modeled with 2D Gauss-Hermite polynomials as generated by the specex package at https://github.com/julienguy/specex

Stephen Bailey December 2013

TODO: GaussHermitePSF (no 2) was copied and pasted from this and then

modified. Can they be refactored to share code?

class specter.psf.gausshermite2.GaussHermite2PSF(filename)

Model PSF with two central Gauss-Hermite cores with different sigmas plus power law wings.

_pgh(x, m=0, xc=0.0, sigma=1.0)

Pixel-integrated (probabilist) Gauss-Hermite function.

Parameters:

• x – pixel-center baseline array

• m – order of Hermite polynomial multiplying Gaussian core

• xc – sub-pixel position of Gaussian centroid relative to x baseline

• sigma – sigma parameter of Gaussian core in units of pixels

Uses the relationship Integral{ H_k(x) exp(-0.5 x^2) dx} = -H_{k-1}(x) exp(-0.5 x^2) + const

Written: Adam S. Bolton, U. of Utah, fall 2010 Adapted for efficiency by S. Bailey while dropping generality

specter.psf.monospot

MonoSpotPSF - …

class specter.psf.monospot.MonoSpotPSF(filename, spot=None, scale=1.0)

_xypix(ispec, wavelength, ispec_cache=None, iwave_cache=None)

Return xslice, yslice, pix for PSF at spectrum ispec, wavelength

specter.psf.pixpsf

Pixelated 2D PSF

David Schlegel & Stephen Bailey, Summer 2011

class specter.psf.pixpsf.PixPSF(filename)

Pixelated PSF[ny, nx] = Ic[ny,nx] + x*Ix[ny, nx] + y*Iy[ny,nx] + …

_xypix(ispec, wavelength, ispec_cache=None, iwave_cache=None)

Evaluate PSF for a given spectrum and wavelength

returns xslice, yslice, pixels[yslice, xslice]

specter.psf.psf

Base class for 2D PSFs

Provides PSF base class which defines the interface for other code using PSFs. Subclasses implement specific models of the PSF and override/extend the __init__ and xypix(ispec, wavelength) methods, while allowing interchangeable use of different PSF models through the interface defined in this base class.

Stephen Bailey, Fall 2012

class specter.psf.psf.PSF(filename)

Base class for 2D PSFs

Subclasses need to extend __init__ to load format-specific items from the input fits file and implement _xypix(ispec, wavelength) to return xslice, yslice, pixels[y,x] for the PSF evaluated at spectrum ispec at the given wavelength. All interactions with PSF classes should be via the methods defined here, allowing interchangeable use of different PSF models.

_fit_spot_sigma(ispec, axis=0, npoly=5)

Fit the cross-sectional Gaussian sigma of PSF spots vs. wavelength. Return callable Legendre object.

Parameters:

• ispec – spectrum number

• axis – 0 or ‘x’ for cross dispersion sigma; 1 or ‘y’ or ‘w’ for wavelength dispersion

• npoly – order of Legendre poly to fit to sigma vs. wavelength

Returns:

legfit such that legfit(w) returns fit at wavelengths w

_value(x, y, ispec, wavelength)

this is implemented in specter.psf.gausshermite and specter.psf.spotgrid, everywhere else just an empty function

_xypix(ispec, wavelength, ispec_cache=None, iwave_cache=None)

Subclasses of PSF should implement this to return xslice, yslice, pixels[iy,ix] for their particular models. Don’t worry about edge effects – PSF.xypix will take care of that.

angstroms_per_pixel(ispec, wavelength)

Return CCD pixel width in Angstroms for spectrum ispec at given wavlength(s). Wavelength may be scalar or array.

cache_params(spec_range, wavelengths)

this is implemented in specter.psf.gausshermite, everywhere else just an empty function

pix(ispec, wavelength)

Evaluate PSF for spectrum[ispec] at given wavelength

returns 2D array pixels[iy,ix]

also see xypix(ispec, wavelength)

project(wavelength, phot, specmin=0, xyrange=None, verbose=False)

Returns 2D image or 3D images of spectra projected onto the CCD

Required inputs:

phot[nwave] or phot[nspec, nwave] or phot[nimage, nspec, nwave]

as photons on CCD per bin

wavelength[nwave] or wavelength[nspec, nwave] in Angstroms

if wavelength is 1D and spectra is 2D or 3D, then wavelength[] applies to all phot[i]

Optional inputs:

specmin : starting spectrum number xyrange : (xmin, xmax, ymin, ymax) range of CCD pixels

if phot is 1D or 2D, output is a single 2D[ny,nx] image if phot is 3D[nimage,nspec,nwave], output is 3D[nimage,ny,nx]

projection_matrix(spec_range, wavelengths, xyrange, use_cache=None)

Returns sparse projection matrix from flux to pixels

Inputs:

spec_range = (ispecmin, ispecmax) or scalar ispec wavelengths = array_like wavelengths xyrange = (xmin, xmax, ymin, ymax)

Optional inputs:

use_cache= default True, legval values will be precomputed

Usage:

xyrange = xmin, xmax, ymin, ymax A = psf.projection_matrix(spec_range, wavelengths, xyrange) nx = xmax-xmin ny = ymax-ymin img = A.dot(phot.ravel()).reshape((ny,nx))

shift_xy(dx, dy)

Shift the x,y trace locations of this PSF while preserving wavelength grid: xnew = x + dx, ynew = y + dy

wavelength(ispec=None, y=None)

Return wavelength of spectrum[ispec] evaluated at y.

ispec can be None, scalar, or vector y can be None, scalar, or vector

May return a view of the underlying array; do not modify unless specifying copy=True to get a copy of the data.

wdisp(ispec, wavelength)

Return Gaussian sigma of PSF spot in wavelength-dispersion direction in units of Angstroms.

Also see ysigma(…) which returns sigmas in units of pixels.

ispec : spectrum index wavelength : scalar or vector wavelength(s) to evaluate spot sigmas

See notes in xsigma(…) about caching of Legendre fit coefficients.

property wmax

Maximum wavelength seen by any spectrum

property wmax_all

Maximum wavelength seen by all spectra

property wmin

Minimum wavelength seen by any spectrum

property wmin_all

Minimum wavelength seen by all spectra

x(ispec=None, wavelength=None)

Return CCD X centroid of spectrum ispec at given wavelength(s).

ispec can be None, scalar, or vector wavelength can be None, scalar or a vector

ispec wavelength returns +——-+———–+—— None None array[nspec, npix_y] None scalar vector[nspec] None vector array[nspec, nwave] scalar None array[npix_y] scalar scalar scalar scalar vector vector[nwave] vector None array[nspec, npix_y] vector scalar vector[nspec] vector vector array[nspec, nwave]

xsigma(ispec, wavelength)

Return Gaussian sigma of PSF spot in cross-dispersion direction in CCD pixel units.

ispec : spectrum index wavelength : scalar or vector wavelength(s) to evaluate spot sigmas

The first time this is called for a spectrum, the PSF is sampled at 20 wavelengths and the variation is fit with a 5th order Legendre polynomial and the coefficients are cached. The actual value (and subsequent calls) use these cached Legendre fits to interpolate the sigma value. If this is not fast enough and/or accurate enough, PSF subtypes may override this function to provide a more accurate xsigma measurement.

xy(ispec=None, wavelength=None)

Utility function to return self.x(…) and self.y(…) in one call

xypix(ispec, wavelength, xmin=0, xmax=None, ymin=0, ymax=None, ispec_cache=None, iwave_cache=None)

Evaluate PSF for spectrum[ispec] at given wavelength

returns xslice, yslice, pixels[iy,ix] such that image[yslice,xslice] += photons*pixels adds the contribution from spectrum ispec at that wavelength.

if xmin or ymin are set, the slices are relative to those minima (useful for simulating subimages)

Optional inputs:

ispec_cache = an index into the spectrum number that starts again at 0 for each patch iwave_cache = an index into the wavelength number that starts again at 0 for each patch

xyrange(spec_range, wavelengths)

Return recommended range of pixels which cover these spectra/fluxes: (xmin, xmax, ymin, ymax)

spec_range = indices specmin,specmax (python style indexing),

or scalar for single spectrum index

wavelengths = wavelength range wavemin,wavemax inclusive

or sorted array of wavelengths

BUG: will fail if asking for a range where one of the spectra is completely off the CCD

y(ispec=None, wavelength=None)

Return CCD Y centroid of spectrum ispec at given wavelength(s).

ispec can be None, scalar, or vector wavelength can be scalar or a vector (but not None)

ispec wavelength returns +——-+———–+—— None scalar vector[nspec] None vector array[nspec,nwave] scalar scalar scalar scalar vector vector[nwave] vector scalar vector[nspec] vector vector array[nspec, nwave]

ysigma(ispec, wavelength)

Return Gaussian sigma of PSF spot in wavelength-dispersion direction in units of pixels.

Also see wdisp(…) which returns sigmas in units of Angstroms.

ispec : spectrum index wavelength : scalar or vector wavelength(s) to evaluate spot sigmas

See notes in xsigma(…) about caching of Legendre fit coefficients.

specter.psf.spotgrid

SpotGridPSF - Linear interpolate hi-res sampled spots to model PSF

Stephen Bailey Fall 2012

class specter.psf.spotgrid.SpotGridPSF(filename)

Model PSF with a linear interpolation of high resolution sampled spots

_value(x, y, ispec, wavelength)

return PSF 2D array (shape from x and y), NOT integrated, for display of PSF.

Parameters:

• x – x-coordinates baseline array (1D or 2D)

• y – y-coordinates baseline array (1D or 2D; if 2D same shape as x)

• ispec – fiber

• wavelength – wavelength

Returns: 2D psf

_xypix(ispec, wavelength, ispec_cache=None, iwave_cache=None)

Return xslice, yslice, pix for PSF at spectrum ispec, wavelength

_xypix_interp(ispec, wavelength)

Return xslice, yslice, pix for PSF at spectrum ispec, wavelength

specter.psf.spotgrid.new_pixshift(xc, yc, pix_spot_values, rebin)

Does a new pixel shift.

Parameters:

• xc (array-like) – Center of the PSF in CCD coordinates.

• yc (array-like) – Center of the PSF in CCD coordinates.

• pix_spot_values (array-like) – 2D array of interpolated values.

• rebin (array-like) – Ratio of spot to CCD pixel size.

Returns:

The resampled and scaled 2D array of pix_spot_values.

Return type:

array-like

specter.throughput

A class for tracking throughput.

Hacks:

• Doesn’t support spatial variation of input sources.

• Doesn’t support per-fiber throughput.

• Do I really want to impose a clumsy ObjType ENUM?

How to handle fiber size and sky units?

specter.throughput.load_throughput(filename)

Create Throughput object from FITS file with EXTNAME=THROUGHPUT HDU

specter.util

specter.util.cachedict

Implement a dictionary that caches the last N entries and throws the rest away

class specter.util.cachedict.CacheDict(n, d=None)

A dictionary that keeps only the last n items added

specter.util.pixspline

Pixel-integrated spline utilities.

Written by A. Bolton, U. of Utah, 2010-2013.

exception specter.util.pixspline.PixSplineError(value)

class specter.util.pixspline.PixelSpline(pixbound, flux)

Pixel Spline object class.

Initialize as follows:

PS = PixelSpline(pixbound, flux)

where

pixbound = array of pixel boundaries in baseline units

(if same len as flux, treat as centers instead of edges)

and

flux = array of specific flux values in baseline units.

Assumptions:

‘pixbound’ should have one more element than ‘flux’, and units of ‘flux’ are -per-unit-baseline, for the baseline units in which pixbound is expressed, averaged over the extent of each pixel.

find_extrema(minima=False)

Return array of extrema of underlying pixelspline

if minima is False, return only maxima

point_evaluate(xnew, missing=0.0)

Evaluate underlying pixel spline at array of points

Also see: PixelSpline.resample()

resample(pb_new)

Resample a pixelspline analytically onto a new set of pixel boundaries

specter.util.pixspline.cen2bound(pixelcen)

Convenience function to do the obvious thing to transform pixel centers to pixel boundaries.

specter.util.traceset

Handle sets of Legendre coefficients

specter.util.traceset.fit_traces(x, yy, deg=5, domain=None)

returns TraceSet object modeling y[i] vs. x

Parameters:

• x – 1D array

• y – 2D array[nspec, nx]

• deg – optional Legendre degree

specter.util.util

Utility functions and classes for specter

Stephen Bailey Fall 2012

class specter.util.util.LinearInterp2D(x, y, data)

Linear interpolation on a 2D grid. Allows values to be interpolated to be multi-dimensional.

specter.util.util._sincfunc(x, dx, dampfac=3.25)

sinc helper function for sincshift()

specter.util.util.custom_erf(y)

Custom implementation of scipy.special.erf() to enable jit-compiling with Numba (which as of 10/2018 does not support scipy). This functionality is equilvalent to:

scipy.special.erf(y)

with the exception that scalar values of y are not supported.

Parameters:

y (array-like) – Points at which the error function will be evaluated.

Returns:

The value of the error function at points in array y.

Return type:

array-like

Notes

This function has been translated from the original fortran function to Python. The original scipy erf function can be found at: https://github.com/scipy/scipy/blob/8dba340293fe20e62e173bdf2c10ae208286692f/scipy/special/cdflib/erf.f Note that this new function introduces a small amount of machine-precision numerical error as compared to the original scipy function.

specter.util.util.custom_hermitenorm(n, u)

Custom implementation of scipy.special.hermitenorm to enable jit-compiling

with Numba (which as of 10/2018 does not support scipy). This functionality is equivalent to: fn = scipy.special.hermitenorm(n) return fn(u) with the exception that scalar values of u are not supported.

Inputs:

n: the degree of the hermite polynomial u: (requires array) points at which the polynomial will be evaulated.

Outputs:

res: the value of the hermite polynomial at array points(u)

specter.util.util.gaussint(x, mean=0.0, sigma=1.0)

Return integral from -inf to x of normalized Gaussian with mean and sigma

specter.util.util.gausspix(x, mean=0.0, sigma=1.0)

Return Gaussian(mean,sigma) integrated over unit width pixels centered at x[].

specter.util.util.rebin_image(image, n)

rebin 2D array pix into bins of size n x n

New binsize must be evenly divisible into original pix image

specter.util.util.resample(x, xp, yp, xedges=False, xpedges=False)

IN PROGRESS. Resample a spectrum to a new binning using PixelSpline

1 <= x.ndim <= xp.ndim <= yp.ndim <= 2

specter.util.util.sincshift(image, dx, dy, sincrad=10, dampfac=3.25)

Return image shifted by dx, dy using sinc interpolation.

For speed, do each dimension independently which can introduce edge effects. Also see sincshift2d().

specter.util.util.sincshift2d(image, dx, dy, sincrad=10, dampfac=3.25)

Return image shifted by dx, dy using full 2D sinc interpolation

specter.util.util.trapz(edges, xp, yp)

Perform trapezoidal integration between edges using sampled function yp vs. xp. Returns array of length len(edges)-1.

Input xp array must be sorted in ascending order.

See also numpy.trapz, which integrates a single array

specter.util.util.weighted_solve(A, b, w)

Solve A x = b with weights w on b Returns x, inverseCovarance(x)