Extracting PSF Images

Generic Code to Reconstruct the PSF

The Point-Spread Function (PSF) reconstruction can be performed without any specialized tools.

To read the PSF info from a psField file for the r band, read extension 3:

IDL> pstruct = mrdfits(psfield_file, 3)

(recall u,g,r,i,z == 0,1,2,3,4 so just add one to the index).

The resulting structure can then be used to reconstruct the image. The following IDL code, taken from the SDSSIDL routine “sdss_psfrec.pro”, demonstrates the algorithm to reconstruct the PSF at location (row,col) in the field.

nrow_b=(pstruct.nrow_b)[0]
ncol_b=(pstruct.ncol_b)[0]
;assumes they are the same for each eigen so only use the 0 one
rnrow=(pstruct.rnrow)[0]
rncol=(pstruct.rncol)[0]
nb=nrow_b*ncol_b
coeffs=fltarr(nb)
ecoeff=fltarr(3)
cmat=pstruct.c
rcs=0.001
FOR i=0L, nb-1L DO coeffs[i]=(row*rcs)^(i mod nrow_b) * (col*rcs)^(i/nrow_b)
FOR j=0,2 DO BEGIN
FOR i=0L, nb-1L DO BEGIN
ecoeff[j]=ecoeff[j]+cmat(i/nrow_b,i mod nrow_b,j)*coeffs[i]
ENDFOR
ENDFOR
psf = (pstruct.rrows)[*,0]*ecoeff[0]+$(pstruct.rrows)[*,1]*ecoeff[1]+$
(pstruct.rrows)[*,2]*ecoeff[2]


Stand Alone Code

We have also created a standalone code that serves to both read them and as a template library for inclusion in other codes. The code is available as: readAtlasImages-v5_4_11.tar.gz.

Compiling

1. make clean
2. make

If you are on a big-endian machine, remove -DSDSS_LITTLE_ENDIAN from CFLAGS in the Makefile.

Using

% read_PSF -h
Usage: read_PSF [options] input-file hdu output-file
-?      This message
-h      This message
-i      Print an ID string and exit
-v      Turn up verbosity (repeat flag for more chatter)


To reconstruct the z PSF (i.e., the 5th HDU) at the position (row, col) = (500, 600) from run 1336, column 2, field 51 you’d say:

%  read_PSF psField-001336-2-0051.fit 5 500.0 600.0 foo.fit

The desired PSF would appear as an unsigned short FITS file in foo.fit; the background level is set to the standard ‘soft bias’ of 1000. If you want a floating image, change a line in the read_PSF.c; look for

        /* create a float region */