PSF Formats
Base PSF
NOTE: we are in the process of switching the format of how x and y vs. wavelength are stored. This will create incompatibilities between old and new PSFs.
HDUs 0 and 1 contain the PSF centroid x and y location vs. wavelength stored as Legendre polynomial coefficients:
HDU 0 : XCOEFF[nspec, ncoeff]
HDU 1 : YCOEFF[nspec, ncoeff]
The header contains keywords WAVEMIN, WAVEMAX to define the wavelength range in Angstroms to map to domain [-1,1] for evaluating the Legendre polynomials. e.g. the x position for spectrum i at wavelength:
w = 2.0*(w-WAVEMIN)/(WAVEMAX-WAVEMIN) - 1.0 #- Map to [-1,1]
x[i] = Sum_j XCOEFF[i,j] L_j(w) #- evaluate Legendre series
Or in Python:
import fitsio
from numpy.polynomial.legendre import legval
h = fitsio.read_header('psf.fits', 'XCOEFF')
xcoeff = fitsio.read('psf.fits', 'XCOEFF')
w = 2.0 * (wavelength-h['WAVEMIN']) / (h['WAVEMAX'] - h['WAVEMIN']) - 1.0
x = legval(w, xcoeff[i])
Additional HDUs contain data specific to various parameterizations of the PSF as listed below. Any additional extensions should be read by EXTNAME, not by number – the order and number of other extensions is arbitrary.
HDU 0 keywords:
NPIX_X, NPIX_Y : CCD dimensions in pixels
(x,y) coordinates are centered on each pixel and start counting at 0. i.e. the lower left corner of the CCD is (-0.5, -0.5) and the center of the lower left pixel is (0,0).
Optional HDU EXTNAME=THROUGHPUT in the format as described in throughput.md . i.e. the throughput may be kept in a separate FITS file, or bundled with the instrument PSF for convenience.
If throughput isn’t available, the PSF can still be used to project photons onto a CCD, but not flux in erg/s/cm^2/A .
Old format for x and y
The old base PSF format parameterized x and y vs. wavelength as:
HDU 0 : x[nspec, nwave] EXTNAME="X"
HDU 1 : y[nspec, nwave] EXTNAME="Y"
HDU 2 : wavelength[nspec, nwave] EXTNAME="WAVELENGTH", or
loglam[nspec, nwave] EXTNAME="LOGLAM"
HDU 3+ : specific to each subformat
We may revive this format.
Spot Grid PSF
PSFTYPE = “SPOTGRID”
This PSF type provides an x,y grid of spots to be interpolated. The first interpolation dimension must monotonically increase with spectrum number, e.g. the position of a fiber along a slit head. The second interpolation dimension is wavelength.
HDU 0-2 : Same as Base PSF: X, Y, wavelength or loglam of traces:
HDU SPOTS : spot[i, j, iy, ix] #- 2D PSF spots
NAXIS1 = i = number of spot samples in the spectrum number dimension
NAXIS2 = j = number of spot samples in the wavelength direction
NAXIS3 = iy = size of spot in the CCD y dimension
NAXIS4 = ix = size of spot in the CCD x direction
HDU SPOTX : spotx[NAXIS1, NAXIS2] #- CCD X pixel loc of spot centroid
HDU SPOTY : spoty[NAXIS1, NAXIS2] #- CCD Y pixel loc of spot centroid
HDU FIBERPOS : fiberpos[nspec] #- Slit position of each fiber
HDU SPOTPOS : spotpos[NAXIS1] #- Slit positions where spots are sampled
HDU SPOTWAVE : spotwave[NAXIS2] #- Wavelengths where spots are sampled
spot[i,j] is a 2D PSF spot sampled at slit position spotpos[i] and
wavelength spotwave[j]. Its center is located on the CCD at
spotx[i,j], spoty[i,j].
Pixellated PSF PCA expansion
PSFTYPE = “PCA-PIX”
This format is a PCA-like pixelated model such that:
pix = ConstImage + x*ImageX + y*ImageY + x*y*ImageXY + ...
HDU 0-2 : Same as Base PSF: X, Y, wavelength or loglam of traces
HDU 3 : table with x and y exponents to use for each image model:
Columns IMODEL, XEXP, YEXP
One row per model image
HDU 4 : table with x and y scale factors:
Columns IGROUP, X0, XSCALE, Y0, YSCALE
One row per fiber
The model images are separately calculated for groups of fibers (could be bundles, but doesn’t have to be). Within each group, the x and y dimensions are rescaled such that:
x = xscale * (xpix - x0)
y = yscale * (ypix - y0)
HDU 5 : 4-D array with the PSF images:
Array[IGROUP, IMODEL, IY, IX]
e.g. to find the PSF for ifiber iflux:
xpix = HDU0[ifiber, iflux]
ypix = HDU1[ifiber, iflux]
x0 = HDU4.X0[ifiber]
y0 = HDU4.Y0[ifiber]
xscale = HDU4.XSCALE[ifiber]
yscale = HDU4.YSCALE[ifiber]
x = xscale*(xpix - x0)
y = yscale*(ypix - y0)
igroup = HDU4.IGROUP[ifiber]
psf = [BlankImage]
for imodel in range( len(HDU3) ):
xexp = HDU3.XEXP[imodel]
yexp = HDU3.YEXP[imodel]
psf += x^xexp * y^yexp * HDU5[igroup, imodel]
Gauss Hermite PSF
This PSF format is generated by the specex PSF package. All parameters are stored in HDU 1 as Legendre coefficients describing how each parameter varies with wavelength for each fiber.
The PSF is modeled by two Gauss-Hermite cores with different sigmas, plus a semi-Lorentzian power-law tail.
Header Keywords for HDU 1
Parameter |
Description |
|---|---|
PSFTYPE |
Must be ‘GAUSS-HERMITE2’ |
PSFVER |
Must be ‘1 ‘ |
NPIX_X |
number of columns in input CCD image |
NPIX_Y |
number of rows in input CCD image |
HSIZEX |
Half size of PSF in fit, NX=2*HSIZEX+1 |
HSIZEY |
Half size of PSF in fit, NY=2*HSIZEY+1 |
BUNDLMIN |
first bundle of fibers (starting at 0) |
BUNDLMAX |
last bundle of fibers (included) |
FIBERMIN |
first fiber (starting at 0) |
FIBERMAX |
last fiber (included) |
NPARAMS |
number of PSF parameters |
LEGDEG |
degree of Legendre pol.(wave) for parameters |
GHDEGX |
degree of Hermite polynomial along CCD columns |
GHDEGY |
degree of Hermite polynomial along CCD rows |
GHDEGX2 |
degree of Hermite polynomial along CCD columns |
GHDEGY2 |
degree of Hermite polynomial along CCD rows |
Binary table in HDU 1
Column |
Description |
|---|---|
PARAM |
Parameter name |
WAVEMIN |
Mininum wavelength |
WAVEMAX |
Maximum wavelength |
COEFF |
Legendre coefficients. Each row is size [NSPEC, LEGDEG+1] |
Each row contains the coefficients for the parameter name in the PARAM column. WAVEMIN and WAVEMAX define the wavelength extent in Angstroms that should be mapped to [-1,1] for evaluating the Legendre polynomials. The coefficients themselves are in the rows of the COEFF column. Each row of COEFF is a 2D array of shape [NSPEC, LEGDEG+1]
The parameters named in the PARAM column are:
Parameter |
Description |
|---|---|
X |
CCD column coordinate (as a function of fiber and wavelength) |
Y |
CCD row coordinate (as a function of fiber and wavelength) <BR/> (X,Y)=(0,0) means that PSF is centered on center of first pixel |
GHSIGX |
Sigma of first Gaussian along CCD columns for PSF core |
GHSIGY |
Sigma of first Gaussian along CCD rows for PSF core |
GHNSIG |
NxSigma cutoff for first Gaussian-Hermite core term |
GHSIGX2 |
Sigma of second Gaussian along CCD columns for PSF wings |
GHSIGY2 |
Sigma of second Gaussian along CCD rows for PSF wings |
GH-i-j |
Hermite pol. coefficents, i along columns, j along rows, <BR/> i is integer from 0 to GHDEGX, j is integer from 0 to GHDEGY, <BR/> there are (GHDEGX+1)*(GHDEGY+1) such coefficents. |
GH2-i-j |
Hermite pol. coefficents for sec. GH psf, i=columns, j=rows |
TAILAMP |
Amplitude of PSF tail |
TAILCORE |
Size in pixels of PSF tail saturation in PSF core |
TAILXSCA |
Scaling apply to CCD coordinate along columns for PSF tail |
TAILYSCA |
Scaling apply to CCD coordinate along rows for PSF tail |
TAILINDE |
Asymptotic power law index of PSF tail |
CONT |
Continuum flux in arc image (not part of PSF) |
These are evaluated as:
PSF_core1(X,Y) = SUM_ij (GH-i-j)*HERM(i,X/GHSIGX)*HERM(j,Y/GHSIGY)
* GAUS(X,GHSIGX)*GAUS(Y,GHSIGY)
--> Evaluated only for pixels where (X/GHSIGX)^2 + (Y/GHSIGY)^2 < GHNSIG^2
PSF_core2(X,Y) = SUM_ij (GH2-i-j)*HERM(i,X/GHSIGX2)*HERM(j,Y/GHSIGY2)
* GAUS(X,GHSIGX2)*GAUS(Y,GHSIGY2)
PSF_tail(X,Y) = TAILAMP*R^2/(TAILCORE^2+R^2)^(1+TAILINDE/2)
with R^2=(X*TAILXSCA)^2+(Y*TAILYSCA)^2
PSF_core is integrated in each pixel.
PSF_tail is not integrated, it is evaluated at the center of each pixel.
For example, if PARAM[0] = ‘X’, this means that row 0 contains the Legendre coefficiencts for the X centroids of the PSFs on the CCD.:
# Map wavelength -> [-1,1] w = 2.0*(wavelength - WAVEMIN[0]) / (WAVEMAX[0]-WAVEMIN[0]) - 1.0
# Evaluate the CCD X centroid for Fiber ispec # x = Sum_i COEFF[ispec,i] * L_i(w) from numpy.polynomial.legendre import legval xcoeff = COEFF[0] #- since PARAM[0] == ‘X’ x = legval(w, xcoeff[ispec])
Gauss Hermite PSF (Deprecated)
NOTE : This format is implemented but somewhat cumbersome. We may abandon it for a simpler format.
This format models the PSF as 2D Gauss-Hermite polynomials. The coefficients of the polynomials are modeled as smoothly varying in (x,y) using 2D Legendre polynomials. Spectra grouped in bundles with a smoothly varying solution within the bundle but independent of other bundles.
HDU 0-2 : Same as Base PSF: X, Y, wavelength or loglam of traces:
PSFTYPE = "PCA-PIX"
PSFVERS = "2.0" or above
HDU 3+ : One HDU per bundle of spectra on the CCD
Header Keywords for HDUs 3+
Keyword |
Meaning |
|---|---|
GHSIGX, GHSIGY |
Gauss-Hermite sigma in x and y directions [pixel units] |
GHDEGX, GHDEGY |
Gauss-Hermite degree in x and y directions |
FIBERMIN, FIBERMAX |
Fibers covered by this bundle min to max inclusive, 0-indexed |
LXMIN, LXMAX |
X-coordinate min/max to transform CCD x -> [-1,1] range for Legendre polynomial |
LYMIN, LYMAX |
Y-coordinate min/max to transform CCD y -> [-1,1] range for Legendre polynomial |
LDEGX, LDEGY |
Degree of Legendre polynomial in x and y directions |
Data in HDU 3+
A 4D image coeff[GHDEGY+1, GHDEGX+1, LDEGY+1, LDEGX+1] for that bundle.
Example
Find the PSF for fiber 5 at wavelength 6000 Angstroms:
Map fiber 5, wavelength 6000 A -> (x,y) on the CCD using HDUs 0-2:
x = numpy.interp(6000, WAVELENGTH[5], X[5])
y = numpy.interp(6000, WAVELENGTH[5], Y[5])
Find which bundle fiber 5 is included in (probably bundle 0 in HDU 3)
Convert x,y -> to ranges [-1,1]
xx = 2*(x-LXMIN)/(LXMAX - LXMIN) - 1
yy = 2*(y-LYMIN)/(LYMAX - LYMIN) - 1
The Gauss-Hermite coefficent c_ij = Sum_kl data[i,j,k,l] L_k(yy) L_l(xx) where L_k is the kth order Legendre Polynomial
PSF(dx, dy) = Sum_ij c_ij H_i(dy/GHSIGY) H_j(dx/GHSIGX)
Then integrate PSF(dx,dy) over the individual pixels
In practice it is better to directly integrate the functions
Notes
This is different from the original Gauss-Hermite format used by bbspec. These may be distingished by the existence of PSFVERS >= 2.0
Future versions of this format may also include additional HDUs to model the wings of the PSF and the covariance of the coefficients.
Other PSF Formats
Specter grew out of “bbspec” which includes PSF formats for:
2D rotated asymmetric Gaussian
An older format for Gauss-Hermite polynomials
These are structurally compatible with Specter PSFs but haven’t been ported yet.