Steps for reducing WFI data with IRAF:

CALIBRATION (can be done using WFIcal.cl)

1. Bias frames. cp bias frames to directory bias. Check for strange bias frames 
with imstat. Average bias frames with noao -> imred -> ccdred -> combine to produce 
bias_master.fit. 

2. Flat fields. cp flat fields to directory flats. Check for good flat fields using
imstat: average should be between 10000 and 40000 ADU. Subtract bias_master from all
flats. Median combine flats with each filter using noao -> imred -> ccdred -> combine
to produce Bflat_master.fit, etc.

3. Subtract bias from all images (noao -> imred -> ccdred -> ccdproc). Do in
bulk with list files. ( cc> epar ccdproc to set lists and approriate calc flags to yes or no.)

4. Divide all object images by flat fields (noao -> imred -> ccdred -> ccdproc). Do in
bulk with list files.

5. Trim object files. (noao -> imred -> ccdred -> ccdproc). Do in
bulk with list files. Good WFI trim section: [16:2058,23:4098]

REDUCTION (for each filter do)

FIND POTENTIAL COMPARISON STARS

1. Pick a representative image. display it. Using ds9 features determine aperture radius, 
inner sky radius and outer sky radius if doing aperture photometry.  Save radii to spreadsheet. Then:
  cl> noao
  no> digiphot
  di> daophot
  da> daoedit image 
   Use 'a' to get star stats in cl window. 
   Guesstimate sky, skysigma, fwhm and min and max ADU from stars a little dimmer
   and brighter than the object star. Sky level can also be obtained from 'h' in
   imexam: the peak of the histogram is the sky level. Record values in spreadsheet. 
  da> epar datapars (or epar daofind, :e on datapars)
   set datamin(sky-[RON x e/ADU x 3] = sky - [4.5 x 1.8 x 3] = sky - 24.2 ~= sky
   ~= 5 x skysigma), skysigma (value measured ~= (sqrt[sky x e/ADU + RON^2])/e/ADU)
   datamax (50000), readnoise (4.5 for WFI), epadu(=1.8 for WFI), 
   exposure=EXPOSURE, airmass=AIRMASS, filter=FILTER, obstime=MJD-OBS
  da> epar centerpars
   calgorithm = centroid, cbox = 2 x FWHM
  da> epar photpars
   apertures (object and sky) = ? (set several for object aperture), zmag = 25.
  *keep sky and object apertures fixed for all images*
Copy MJD from header (da> imhead l+ | page) to spreadsheet. 
   
2. Using parameters from step 1 (saved in spreadsheet, set with da> epar datapars), 
isolate stars with:
   da> daofind image
    detection threshold set high enough to include only reasonable comp stars.
   writes image.coo.1
    
3. Produce comp star map using tvmark and the coo file from step 2, produce a snapped
image for printing. Use aperture and outer sky radius circles:
  da> epar tvmark
   set mark=circle, radii=a,b (a=aperture, b=outer sky), 
   outimage = Bcomps (or Vcomps, etc.)
  da> tvmark 1 image.coo.1
print Bcomps.fits image
    
4. Trim out bad comps from coo file (must be in daophot):
   da> epar txdump 
    fields = "xcenter,ycenter,id"
   da> txdump >> image.coo
   gedit image.coo
   trim list, looking at display and printed snap image for merged stars, etc.
Record good potential comps.
Redo step 3 if desired, with image.coo in place of image.coo.1


APERTURE PHOTOMETRY

1. Get sky and fwhm parameters (repeats of step 1 from FIND POTENTIAL COMPARISON STARS).
  cl> noao
  no> digiphot
  di> daophot
  da> daoedit image 
   Use 'a' to get star stats in cl window. 
   Guesstimate sky, skysigma, fwhm and min and max ADU from stars a little dimmer
   and brighter than the object star. Sky level can also be obtained from 'h' in
   imexam: the peak of the histogram is the sky level. Record values in spreadsheet. 
  da> epar datapars (or epar phot, :e on datapars)
   set datamin(sky-[RON x e/ADU x 3] = sky - [4.5 x 1.8 x 3] = sky - 24.2 ~= sky
   ~= 5 x skysigma), skysigma (value measured ~= (sqrt[sky x e/ADU + RON^2])/e/ADU)
   datamax (50000), readnoise (4.5 for WFI), epadu(=1.8 for WFI), 
   exposure=EXPOSURE, airmass=AIRMASS, filter=FILTER, obstime=MJD-OBS
  Important:
  da> epar centerpars
   calgorithm = centroid, cbox = 2 x FWHM
  da> epar photpars
   apertures (object and sky) = ? (set several for object aperture), zmag = 25.
   set datamin low enough to get sky
  da> epar fitskypars
    salgorithm = mode
    annulus = ?
    dannulus = ?
  If doing manually (recommended) set interactive to yes, coords file to nothing, 
  and maxshift in centerpars to 3.
  *keep sky and object apertures fixed for all images*
Copy MJD from header (da> imhead l+ | page) to spreadsheet. 
   
2. Do aperture photometry:
   If using coord file:
   da> phot image image.coo  (trimmed coo file from step 4 of FIND POTENTIAL COMPARISON STARS)
    output written to image.mag.1
   If manually:
   da> phot image
    Hit space bar to select star, "q" and "q" to quit and save image.mag.1
    
3. Extract data from image.mag.1 using txdump (must be in daophot). 
  da> epar txdump
   fields = "xcenter,ycenter,mag,merr,id"
  da> txdump >> objectphot
gedit objectphot
trim out bad comps from objectphot 
Copy largest aperture data (mag and merr) to spreadsheet. (check that
magnitudes have "converged".)

4. Repeat all steps (use same comps - and always use the
same apertures) for rest of stars.


PSF PHOTOMETRY

1. Get sky and fwhm parameters (repeats of step 1 from FIND POTENTIAL COMPARISON STARS).
  cl> noao
  no> digiphot
  di> daophot
  da> daoedit image 
   Use 'a' to get star stats in cl window. 
   Guesstimate:
    sky, skysigma, fwhm and min and max ADU 
   from stars a little dimmer and brighter than the object star. 
    Sky level can also be obtained from 'h' in 
   imexam: the peak of the histogram is the sky level. 
Record values in spreadsheet. 
  da> epar datapars (or epar daofind, :e on datapars)
   set datamin(sky-[RON x e/ADU x 3] = sky - [4.5 x 1.8 x 3] = sky - 24.2 ~= sky
   ~= 5 x skysigma), skysigma (value measured ~= (sqrt[sky x e/ADU + RON^2])/e/ADU)
      *[set datamin low enough to get sky]*
   datamax (50000), readnoise (4.5 for WFI), epadu(=1.8 for WFI), 
   exposure=EXPOSURE, airmass=AIRMASS, filter=FILTER, obstime=MJD-OBS
  Important:
  da> epar centerpars
   calgorithm = centroid, cbox = 2 x FWHM
   maxshif = 1 
  da> epar photpars
   aperture = 6, zmag = 25.
  da> epar fitskypars
    salgorithm = mode
    annulus = 30
    dannulus = 10
  *keep sky and object apertures fixed for all images*
Copy MJD from header (da> imhead l+ | page) to spreadsheet. 

2. Find stars for aperture photometry psf definition for use with profile 
fitting photometry. 
Using parameters from step 1 (saved in spreadsheet, set with da> epar datapars), 
find stars with:
   da> daofind image
    detection threshold 4 or 5 works well - get faint stars around object and
    planned comps.
    writes image.coo.1
Check stars found with
  da> epar tvmark
  set mark=point, make sure outimage is blank
  da> tvmark
Erase image.coo.1 and redo with higher or lower threshold if necessary.

3. Setup for aperture photometry for use with profile fitting photometry.
   *Important:
   da> epar centerpars
    calgorithm = none
   Double check from step 1: 
   da> epar photpars
   aperture = 6, zmag = 25.
  da> epar fitskypars
    salgorithm = mode
    annulus = 10
    dannulus = 10
   da> epar datapars
     set datamin low enough to get sky
    
4. Do aperture photometry for use with profile fitting photometry:
   da> phot image
   inner sky radius = 10, width = 10, aperture radius = ~FWHM = fixed = 6
    output written to image.mag.1
   
5. Select PSF model stars.
   da> epar daopars
    function = auto
    psfrad = 30
    fitrad = 6 (~= FWHM from step 1)
   da> epar pstselect
    maxnpsf = 5
    interac = yes
   da> pstselect image
    'a' in ds9 then 'a' in plot window gets star
    [also select comp stars with 'a' then 'd' in plot window - 
    record comp star numbers in logbook]
    'q' in ds9 then 'q' in cl saves data in object.fit.pst.1
    
6. Generate model PSF
   da> psf image image.mag.1 image.pst.1 interac-
    (files can be lists)
    writes
    PSF image image.psf.1.fits
    output PSF star list image.pst.1
    output PSF star group file image.psg.1
    
7. Optional and time consuming:
Look at the "group" file image.psg.1. If the groups have more
than one member then those neighbors need to be subtracted according to
the procedure given next. Otherwise go to step 14. To subtract the neighbors
follow:
(a)
da> nstar
*reduce psf (definition) radius psfrad to 1/2 the value used in step 12.
verify group file to image.psg.1 (not default name)
verify centering to yes
psf photometry for the stars in image.psg.1 will be output to
image.nst.1
(b)
da> substar
verify psf radius to 1/2 the value used in step 12.
subtract image.nst.1 psfs from image to produce image.sub.1.fits - stars
listed in image.psg.1 will have been subtracted
If subtraction is not clean go back to step 11.
(c)
Edit image.nst.1 to remove PSF model stars (1st in each group)
(d)
da> substar
as step (b) but now image.sub.2.fits will have no neighbors around the
PSF model stars.
(e)
Repeat steps 11. and 12. on image.sub.2.fits. Make sure to increase psfrad
back to original large value. Then repeat steps (a) and (b) to check results.
  
8. Examine model PSF
    da> seepsf image.psf.1.fits image.psf.1s.fits
     input file is psf file from step 6 (or 7), output is fits of PSF
     
9. Do profile photometry.
    da> allstar image image.mag.1 image.psf.1.fits
    [allstar image, mag file from step 6 (aperture photometry), model psf 
    file from step 4]
    writes: image.als.1 and image.sub.1.fits (also writes object.fit.pst.2)
    
10. Check subtracted file image.sub.1.fits to see if fitted stars are
reasonably gone. If not redo psf steps from step 6.
    
11. Extract data from image.als.1 to test file
    da> epar txdump
     fields = "xcenter,ycenter,mag,merr,id"
    da> txdump >> imagephotpsf
gedit search for comp star numbers and transfer to spreadsheet

12. Repeat all steps (use same comps - and always use the
same apertures) for rest of stars.


CONDENSED PSF PHOTOMETRY (after parameter files are stable)

1. da> disp image
         record fwhm, MJD

2. da> daofind image
       (verify calgorithm = centroid)

   da> tvmark
     
4. da> phot image
       (verify calgorithm = none)

5. da> pstselect image

6. da> psf image image.mag.1 image.pst.1

9. da> allstar image image.mag.1 image.psf.1.fits

11. da> txdump >> imagephotpsf


*** BULK PSF PHOTOMETRY ***

1. make list of images Flist

2. da> daofind @Flist

3. da> phot @Flist

4. da> disp image 1

5. da> tvmark 1

6. da> pstselect image

7. cycle through 4 to 6, mark and record psf/comp stars for all frames

8. da> psf @Flist

9. da> allstar @Flist

10. da> txdump image.als.1 >>  imagephotpsf
       extract comps with id == 1 || id == 5 etc.
       repeat for all frames
       
11. da>  hselect *Fcal.fit $I,MJD-OBS yes > dates      
       
       
*** GOOD IRAF PARAMETERS FOR SSO 40 INCH WFI ***
notes: 
- typical FWHM of star PSF is 6 (2.25 arcsec seeing)
- aperture photometry optimized for subsequent psf fitting,
  not good for straight aperture photometry
***

 PACKAGE = daophot
   TASK = daofind

image   =               @Ilist  Input image(s)
output  =              default  Output coordinate file(s) (default: image.coo.?)
(starmap=                     ) Output density enhancement image(s)
(skymap =                     ) Output sky image(s)
(datapar=                     ) Data dependent parameters
(findpar=                     ) Object detection parameters
(boundar=              nearest) Boundary extension (constant|nearest|reflect|wrap)
(constan=                   0.) Constant for boundary extension
(interac=                   no) Interactive mode ?
(icomman=                     ) Image cursor: [x y wcs] key [cmd]
(gcomman=                     ) Graphics cursor: [x y wcs] key [cmd]
(wcsout =            )_.wcsout) The output coordinate system (logical,tv,physical)
(cache  =             )_.cache) Cache the image pixels ?
(verify =            )_.verify) Verify critical daofind parameters ?
(update =            )_.update) Update critical daofind parameters ?
(verbose=           )_.verbose) Print daofind messages ?
(graphic=          )_.graphics) Graphics device
(display=           )_.display) Display device
(mode   =                   ql)

***

 PACKAGE = daophot
   TASK = datapars

(scale  =                   1.) Image scale in units per pixel
(fwhmpsf=                  7.8) FWHM of the PSF in scale units
(emissio=                  yes) Features are positive ?
(sigma  =                  30.) Standard deviation of background in counts
(datamin=                   0.) Minimum good data value
(datamax=               50000.) Maximum good data value
(noise  =              poisson) Noise model
(ccdread=                     ) CCD readout noise image header keyword
(gain   =                     ) CCD gain image header keyword
(readnoi=                  4.5) CCD readout noise in electrons
(epadu  =      1.7999999523163) Gain in electrons per count
(exposur=              EXPTIME) Exposure time image header keyword
(airmass=              AIRMASS) Airmass image header keyword
(filter =               FILTER) Filter image header keyword
(obstime=              MJD-OBS) Time of observation image header keyword
(itime  =      300.01000976562) Exposure time
(xairmas=      1.3832819461823) Airmass
(ifilter=                    B) Filter
(otime  =   5.360576979185E+04) Time of observation
(mode   =                   ql)

***

 PACKAGE = daophot
   TASK = findpars

(thresho=                   4.) Threshold in sigma for feature detection
(nsigma =                  1.5) Width of convolution kernel in sigma
(ratio  =                   1.) Ratio of minor to major axis of Gaussian kernel
(theta  =                   0.) Position angle of major axis of Gaussian kernel
(sharplo=                  0.2) Lower bound on sharpness for feature detection 
(sharphi=                   1.) Upper bound on sharpness for  feature detection
(roundlo=                  -1.) Lower bound on roundness for feature detection
(roundhi=                   1.) Upper bound on roundness for feature detection
(mkdetec=                  yes) Mark detections on the image display ?
(mode   =                   ql)

***

 PACKAGE = daophot
   TASK = phot

image   =               @Ilist  Input image(s)
coords  =              default  Input coordinate list(s) (default: image.coo.?)
output  =              default  Output photometry file(s) (default: image.mag.?)
skyfile =                       Input sky value file(s)
(plotfil=                     ) Output plot metacode file
(datapar=                     ) Data dependent parameters
(centerp=                     ) Centering parameters
(fitskyp=                     ) Sky fitting parameters
(photpar=                     ) Photometry parameters
(interac=                   no) Interactive mode ?
(radplot=                   no) Plot the radial profiles?
(icomman=                     ) Image cursor: [x y wcs] key [cmd]
(gcomman=                     ) Graphics cursor: [x y wcs] key [cmd]
(wcsin  =             )_.wcsin) The input coordinate system (logical,tv,physical,world)
(wcsout =            )_.wcsout) The output coordinate system (logical,tv,physical)
(cache  =             )_.cache) Cache the input image pixels in memory ?
(verify =            )_.verify) Verify critical phot parameters ?
(update =            )_.update) Update critical phot parameters ?
(verbose=           )_.verbose) Print phot messages ?
(graphic=          )_.graphics) Graphics device
(display=           )_.display) Display device
(mode   =                   ql)


***  
 
PACKAGE = daophot
   TASK = centerpars

(calgori=             centroid) Centering algorithm
(cbox   =                  12.) Centering box width in scale units
(cthresh=                   0.) Centering threshold in sigma above background
(minsnra=                   1.) Minimum signal-to-noise ratio for centering algorithim
(cmaxite=                   10) Maximum iterations for centering algorithm
(maxshif=                   1.) Maximum center shift in scale units
(clean  =                   no) Symmetry clean before centering
(rclean =                   1.) Cleaning radius in scale units
(rclip  =                   2.) Clipping radius in scale units
(kclean =                   3.) K-sigma rejection criterion in skysigma
(mkcente=                   no) Mark the computed center
(mode   =                   ql)

***

PACKAGE = daophot
   TASK = fitskypars

(salgori=                 mode) Sky fitting algorithm
(annulus=                  30.) Inner radius of sky annulus in scale units
(dannulu=                  10.) Width of sky annulus in scale units
(skyvalu=                  50.) User sky value
(smaxite=                   10) Maximum number of sky fitting iterations
(sloclip=                   0.) Lower clipping factor in percent
(shiclip=                   0.) Upper clipping factor in percent
(snrejec=                   50) Maximum number of sky fitting rejection iterations
(sloreje=                   3.) Lower K-sigma rejection limit in sky sigma
(shireje=                   3.) Upper K-sigma rejection limit in sky sigma
(khist  =                   3.) Half width of histogram in sky sigma
(binsize=     0.10000000149012) Binsize of histogram in sky sigma
(smooth =                   no) Boxcar smooth the histogram
(rgrow  =                   0.) Region growing radius in scale units
(mksky  =                  yes) Mark sky annuli on the display
(mode   =                   ql)

***

PACKAGE = daophot
   TASK = photpars

(weighti=             constant) Photometric weighting scheme
(apertur=                    6) List of aperture radii in scale units
(zmag   =                  25.) Zero point of magnitude scale
(mkapert=                  yes) Draw apertures on the display
(mode   =                   ql)

***

PACKAGE = daophot
   TASK = pstselect

image   =          ccd0126Ical  Image for which to build psf star list
photfile=              default  Photometry file (default: image.mag.?)
pstfile =              default  Output psf star list file (default: image.pst.?)
maxnpsf =                   20  Maximum number of psf stars
(mkstars=                   no) Mark deleted and accepted psf stars
(plotfil=                     ) Output plot metacode file
(datapar=                     ) Data dependent parameters
(daopars=                     ) Psf fitting parameters
(interac=                  yes) Select psf stars interactively ?
(plottyp=                 mesh) Default plot type (mesh|contour|radial)
(icomman=                     ) Image cursor: [x y wcs] key [cmd]
(gcomman=                     ) Graphics cursor: [x y wcs] key [cmd]
(wcsin  =             )_.wcsin) The input coordinate system (logical,tv,physical,world)
(wcsout =            )_.wcsout) The output coordinate system (logical,tv,physical)
(cache  =             )_.cache) Cache the input image pixels in memory ?
(verify =            )_.verify) Verify critical pstselect parameters?
(update =            )_.update) Update critical pstselect parameters?
(verbose=           )_.verbose) Print pstselect messages?
(graphic=          )_.graphics) Graphics device
(display=           )_.display) Image display device
(mode   =                   ql)

***

PACKAGE = daophot
   TASK = daopars

(functio=                 auto) Form of analytic component of psf model
(varorde=                    0) Order of empirical component of psf model
(nclean =                    0) Number of cleaning iterations for computing psf model
(saturat=                   no) Use wings of saturated stars in psf model computation ?
(matchra=                   3.) Object matching radius in scale units
(psfrad =                  30.) Radius of psf model in scale units
(fitrad =                   10.) Fitting radius in scale units
(recente=                  yes) Recenter stars during fit ?
(fitsky =                   no) Recompute group sky value during fit ?
(groupsk=                  yes) Use group rather than individual sky values ?
(sannulu=                   0.) Inner radius of sky fitting annulus in scale units
(wsannul=                  11.) Width of sky fitting annulus in scale units
(flaterr=                 0.75) Flat field error in percent
(proferr=                   5.) Profile error in percent
(maxiter=                   50) Maximum number of fitting iterations
(clipexp=                    6) Bad data clipping exponent
(clipran=                  2.5) Bad data clipping range in sigma
(mergera=                INDEF) Critical object merging radius in scale units
(critsnr=                   1.) Critical S/N ratio for group membership
(maxnsta=                10000) Maximum number of stars to fit
(maxgrou=                   60) Maximum number of stars to fit per group
(mode   =                   ql)

***

PACKAGE = daophot
   TASK = allstar

image   =          ccd0126Ical  Image corresponding to photometry
photfile=              default  Input photometry file (default: image.mag.?)
psfimage=              default  PSF image (default: image.psf.?)
allstarf=              default  Output photometry file (default: image.als.?)
rejfile =              default  Output rejections file (default: image.arj.?)
subimage=              default  Subtracted image (default: image.sub.?)
(datapar=                     ) Data dependent parameters
(daopars=                     ) Psf fitting parameters
(wcsin  =             )_.wcsin) The input coordinate system (logical,tv,physical,world)
(wcsout =            )_.wcsout) The output coordinate system (logical,tv,physical)
(wcspsf =            )_.wcspsf) The psf coordinate system (logical,tv,physical)
(cache  =                  yes) Cache the data in memory?
(verify =            )_.verify) Verify critical allstar parameters?
(update =            )_.update) Update critical allstar parameters?
(verbose=           )_.verbose) Print allstar messages?
(mode   =                   ql)

***

PACKAGE = daophot
   TASK = txdump

textfile=    ccd0125Ical.als.1  Input apphot/daophot text database(s)
fields  = xcenter,ycenter,mag,merr,id  Fields to be extracted
expr    = id == 608 || id == 749 || id == 784 || id == 916 || id == 1969 || id == 1239  Boolea
(headers=                   no) Print the field headers ?
(paramet=                  yes) Print the parameters if headers is yes ?
(mode   =                   ql) Mode of task


*** Inserting HJD into fits header ***

Use setjd in noao.astutil with:

as> lpar setjd
       images = "@list"         Images
 (observatory = "sso")          Observatory of observation
        (date = "date-obs")     Date of observation keyword
        (time = "ut")           Time of observation keyword
    (exposure = "exptime")      Exposure time keyword
          (ra = "ra")           Right ascension (hours) keyword
         (dec = "dec")          Declination (degrees) keyword
       (epoch = "equinox")      Epoch (years) keyword\n
          (jd = "jd")           Output Julian date keyword
         (hjd = "hjd")          Output Helocentric Julian date keyword
         (ljd = "ljd")          Output local Julian date keyword\n
      (utdate = yes)            Is observation date UT?
      (uttime = yes)            Is observation time UT?
    (listonly = no)             List only without modifying images?
        (mode = "ql")           

where list is a file with all the fits image names.

The output written to the screen contains:
imagename  jd  hjd  integerpart

Note: The jd time stored originally in the SSO header fits file is the
time at the start of the exposure if 0.5 is added to it.

Extract hjd using hselect:

da> hselect *cal.fit $I,hjd yes > hjdates

or do for each filter:
da> hselect *Ical.fit $I,hjd yes > Ihjdates
da> hselect *Rcal.fit $I,hjd yes > Rhjdates
da> hselect *Vcal.fit $I,hjd yes > Vhjdates
da> hselect *Bcal.fit $I,hjd yes > Bhjdates

this makes copying to spreadsheets easier.


*** Process for doing inhomogeneous ensemble photometry ***

Notes.

1. The inhomogeneous photometry (IEP) software is available from:
http://spiff.rit.edu/ensemble/

2. Install the dump.cl file in /iraf (or wherever you want). Here's dump.cl:

string *inlist, *outlist
inlist = "inlist"
outlist = "outlist"
while (fscan(inlist,s1) != EOF) {
  if(fscan(outlist,s2) == 1) {
  print(s1, " -> ", s2)
  txdump (s1, "xcenter,ycenter,mag,merr,id", yes, >> s2)
  }
  ;
}

(another script I have written, extract.cl. will extract the dates)


Procedure.

A. If hjd dates haven't been put in the fits header:

1. In the image directory make the file "inlist" with the names
of the als files listed:
% ls *.als* > inlist
(make sure that multiple versions of als files are not present)

2. % cp inlist outlist

3. edit image.als.1 to imageall in outlist

4. goto noao.digiphot.daophot in IRAF and run:
da> cl < /iraf/dump.cl

This produces, for each image, a file (imageall) with columns
x   y   instmag   err  id

4b.
goto noao.digiphot.daophot in IRAF and run:
da> cl < /iraf/extract.cl

This produces, for each image, a file (imageall) with columns
x   y   instmag   err  id 	hjd


5. Sort the output imageall files by filter to make input lists for
multipht; use the names:
Bmplist
Vmplist
Rmplist
Implist
% ls *Bcalall > Bmplist
% ls *Vcalall > Vmplist
% ls *Rcalall > Rmplist
% ls *Icalall > Implist

6. Further sort by pierside if necessary, use:
da> hselect *Bcal.fit $I,pierside yes > Bpierinfo
da> hselect *Vcal.fit $I,pierside yes > Vpierinfo
da> hselect *Rcal.fit $I,pierside yes > Rpierinfo
da> hselect *Ical.fit $I,pierside yes > Ipierinfo

Then split Bmplist into (use gedit Bpierinfo)
BmplistE
BmplistW
% cp Bmplist BmplistE
% cp Bmplist BmplistW

split Vmplist into (use gedit Vpierinfo)
VmplistE
VmplistW
% cp Vmplist VmplistE
% cp Vmplist VmplistW

split Rmplist into (use gedit Rpierinfo)
RmplistE
RmplistW
% cp Rmplist RmplistE
% cp Rmplist RmplistW

split Implist into (use gedit Ipierinfo)
ImplistE
ImplistW
% cp Implist ImplistE
% cp Implist ImplistW

If pier sort was necessary then follow:

7. Run multipht:
% multipht print x=0 y=1 mag=2 err=3 auto list=BmplistE outfile=BsplistE star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=BmplistW outfile=BsplistW star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=VmplistE outfile=VsplistE star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=VmplistW outfile=VsplistW star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=RmplistE outfile=RsplistE star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=RmplistW outfile=RsplistW star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=ImplistE outfile=IsplistE star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=ImplistW outfile=IsplistW star_cut=0.1 image_cut=3

Adjust star_cut to capture 500 to 1000 stars. Follow each multipht with solvepht (step 8).

8. Run solvepht
% solvepht infile=BsplistE outfile=BsplistE.out imfile=BsplistE.img sigfile=BsplistE.sig niter=1
% solvepht infile=BsplistW outfile=BsplistW.out imfile=BsplistW.img sigfile=BsplistW.sig niter=1
% solvepht infile=VsplistE outfile=VsplistE.out imfile=VsplistE.img sigfile=VsplistE.sig niter=1
% solvepht infile=VsplistW outfile=VsplistW.out imfile=VsplistW.img sigfile=VsplistW.sig niter=1
% solvepht infile=RsplistE outfile=RsplistE.out imfile=RsplistE.img sigfile=RsplistE.sig niter=1
% solvepht infile=RsplistW outfile=RsplistW.out imfile=RsplistW.img sigfile=RsplistW.sig niter=1
% solvepht infile=IsplistE outfile=IsplistE.out imfile=IsplistE.img sigfile=IsplistE.sig niter=1
% solvepht infile=IsplistW outfile=IsplistW.out imfile=IsplistW.img sigfile=IsplistW.sig niter=1

(be sure to look at the .img file to check for bad frames)


If pier sort was not necessary then follow:

7. Run multipht:
% multipht print x=0 y=1 mag=2 err=3 auto list=Bmplist outfile=Bsplist star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=Vmplist outfile=Vsplist star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=Rmplist outfile=Rsplist star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 auto list=Implist outfile=Isplist star_cut=0.1 image_cut=3

7b. If dates were extracted using extract.cl then use:

% multipht print x=0 y=1 mag=2 err=3 jd=5 auto list=Bmplist outfile=Bsplist star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 jd=5 auto list=Vmplist outfile=Vsplist star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 jd=5 auto list=Rmplist outfile=Rsplist star_cut=0.1 image_cut=3
% multipht print x=0 y=1 mag=2 err=3 jd=5 auto list=Implist outfile=Isplist star_cut=0.1 image_cut=3

Adjust star_cut to capture 500 to 1000 stars. Follow each multipht with solvepht (step 8).

8. Run solvepht
% solvepht infile=Bsplist outfile=Bsplist.out imfile=Bsplist.img sigfile=Bsplist.sig niter=1 varstar=XXX
% solvepht infile=Vsplist outfile=Vsplist.out imfile=Vsplist.img sigfile=Vsplist.sig niter=1 varstar=XXX
% solvepht infile=Rsplist outfile=Rsplist.out imfile=Rsplist.img sigfile=Rsplist.sig niter=1 varstar=XXX
% solvepht infile=Isplist outfile=Isplist.out imfile=Isplist.img sigfile=Isplist.sig niter=1 varstar=XXX

(be sure to look at the .img file to check for bad frames)


*** 

Saskatoon physics roof information
in /home/iraf/iraf/noao/lib/obsdb.dat

observatory = "saskatoon"
	name = "U of S Physics Building"
	longitude = 106.683
	latitude = 52.167
	altitude = 500
	timezone = 6

*** imalign, etc. ***

1. run daofind

Note: for U of S data use    
(threshold = 20.)            Threshold in sigma for feature detection
in findpars and
(datamin = 500.)           Minimum good data value
in datapars

To use imalign:

1. Pick a reference image.

2. epar imexam to set
(logfile=           coords.txt) logfile
(keeplog=                  yes) log output results

3. Using imexam, use , on stars of interest,
q to finish

4. turn logfile in epars imexam to no to prevent future writing to
coords.txt

1. Put input images into list:
% ls *.fit > flist

2. % cp flist slist
modify for output image names

3. set imalign parameters as (this worked once - bigbox and number of iterations
need trial and error adjustment):

da> lpar imalign
        input = "@flist"        Input images
    reference = "LPH1001-20060505E100-001I" Reference image
       coords = "coords.txt"    Reference coordinates file
       output = "@slist"        Output images
      (shifts = "")             Initial shifts file
     (boxsize = 7)              Size of the small centering box
      (bigbox = 101)            Size of the big centering box
    (negative = no)             Are the features negative ?
  (background = INDEF)          Reference background level
       (lower = INDEF)          Lower threshold for data
       (upper = INDEF)          Upper threshold for data
    (niterate = 20)             Maximum number of iterations
   (tolerance = 0)              Tolerance for convergence
    (maxshift = INDEF)          Maximum acceptable pixel shift
 (shiftimages = yes)            Shift the images ?
 (interp_type = "linear")       Interpolant
(boundary_typ = "nearest")      Boundary type
    (constant = 0.)             Constant for constant boundary extension
  (trimimages = yes)            Trim the shifted images ?
     (verbose = yes)            Print the centers, shifts, and trim section ?
        (list = "")             
        (mode = "ql")           


*** Procedure for U of S data ***

Assuming all data have been calibrated in MaxIm:

1. Set hjd in the images using setjd (see above). Sometimes ra, dec and equinox
are missing - these have to be added to the image headers with hedit before setjd will
work.

2. Pick a reference image and align images using imalign (see above).

3. Find the stars in the reference image using daofind. This will make the reference
coo file.

4. Make a file with repeated names of the reference coo file. This list will be used for
phot next.

5. Run phot to produce the reference aperture photometry.

6. Run pstselect on the reference image (disp refimage 1, tvmark 1 refimage.coo.1,
pstselect refimage).

7. Make a file with repeated names of the reference pst file. This list will be used for
psf next.

8. Run psf.

9. Run allstar.

10. Extract hjd from headers as given above.

11. Run IEP as given above.