YALO IR - April 2001

YALO IR channel notes for SN2000cx data (7-9/2000)
written 4/2001
Basic calibration for your data 2

SN2000cx data reductions
Data taken June-Sept 2001

April 2001

DATA REDUCTION

Gain 6.6e-/ADU
read noise 14.1e-

PRELIMIARIES:

0. Copy all the images from fits to imh.

cpimh ir*.fits del+

1. The YALO FITS headers have some features which I change.

equinox ==> epoch
observat ==> "ctio"
and move the jd to JD-OLD

I run the script "yalohead" to convert the FITS headers into something
more standard for IRAF.

yalohead ir*.imh
setjd ir*.imh date="UTDATE" time="UT" exposure="EXPTIME" epoch="EQUINOX"
setairmass ir*.imh

2. I have written a small task to put a dither parameter in the header.

If you have standards that were taken with dithers, you may want to use this.

Check the tilt parameters first as:

hsel *.imh $I,tilt1,tilt2,tilt3 yes
dtilt *.imh

dtilt:

          images = "a*.imh" input images
  (dither = 40 Tilt step: 10,20,20, etc
  (tilt1 = 1320 Tilt position 1
  (tilt2 = 2225 Tilt position 2
  (tilt3 = 1860 Tilt position 3
  (imglist = "tmp$tmp.562ga")  
  (mode = "ql"  

 

 

BASIC CCDRED STUFF.

Making the biases:

The IR detector has a numerical bias of 400 units. On top of that, the dark frame at the same exptime as an object frame has warm pixels that are similar to biases. The biases we will use are in order of preference:

1. An averaged dark taken at the same time as the object frame.

Check to see if the darks look okay. Sometimes the first one is bad.

zerocomb nickdark.????.imh out=irdark45 comb=med rej=minmax nlow=1 nhigh=1
displ irdark45 1

2. The DOME_OFF frame

Note the the K DOME_OFF actually has some light on it. You must do a getsky on this image, see what the sky value is, and subtract a constant to bring it to 400.

imar ir000323.K_OFF - 460 ir000323.K_OFF

3. A numerical bias frame with 400. in all the pixels. If you have to make a numerical bias, then:

imcopy ir991121.flatj zero400
imrep zero400 400. lower=INDEF upper=INDEF
hedit zero400 IMAGETYP zero up+
hedit zero400 title "Numerical bias of 400" up+


4 IMPORTANT! Whatever bias you are using, you must declare the image as a ZERO image.

hedit irdark45 IMAGETYP zero up+ ver-

 

 

MAKING THE FLATS AND CORRECTING THE VIGNETTING


1. The flats are calculated as DOME_ON-DOME_OFF.

YOU MUST EDIT THE HEADER OF YOUR FLAT FRAMES TO SHOW THAT THESE FLATS ARE ZERO-CORRECTED. (They are already zero-corrected because they were calculated as DOME_ON-DOME_OFF). IF YOU DON'T DO THIS, THE FLATS WILL BE ZERO-CORRECTED BY CCDPR, AND THIS IS VERY WRONG!

hedit irflat?.????.imh ZEROCOR "Corrected by DOME_OFF" add+

2. The flats may have 0 value pixels, which will cause the ccdpr to crash.

The low pixels should be replaced by a large number (here I chose saturation) so that in the division, they will be close to 0. You may want to run the flats through imrep as:

imreplace irflat?.????.imh 15000 lower=INDEF upper=1

3. Now we flatten the data

The data are now [ZF].

ccdr:

(pixeltype = "real real") Output and calculation pixel datatypes
(verbose = yes) Print log information to the standard output?
(logfile = "logfile") Text log file
(plotfile = "") Log metacode plot file
(backup = "") Backup directory or prefiz
(instrument = "myiraf$/yalo_ir.dat") CCD instrument file
(ssfile = "myiraf$/yalo_ir.sub") Subset translation file
(graphics = "stdgraph") Interactive graphics output device
(cursor = "") Graphics cursor input
(version = "2:October 1987")  
(mode = "")  
($nargs = 0)  

 

ccdpr:

images = "a*.imh" List od CCD images to correct
(output = "") List of output CCD images
(ccdtype = "") CCD image type to correct
(max_cache = 0) Maximum image caching memory (in Mbytes)
(noproc = no) List processing steps only?\n
(fixpix = no) Fix bad CCD lines and columns?
(overscan = no) Apply overscan strip correction?
(trim = no) Trim the image
(zerocor = yes) Apply zero level correction?
(darkcor = no) Apply dark count correction?
(flatcor = no) Apply flat field correction?
(illumcor = no) Apply illumination correction?
(fringecor = no) Apply fringe correction?
(readcor = no) Convert zero level image to readout correction?
(scancor = no) Convert flat field image to scan correction?\n
(readaxis = "line") Read out axis (column|line)
(fixfile = "") File describing the bad lines and columns
(biassec = "") Overscan strip image section
(trimsec = "") Trim data section
(zero = "irdark45") Zero level calibration image
(dark = "") Dark count calibration image
(flat = "irflat?.imh") Flat field images
(illum = "") Illumination correction images
(fringe = "") Fringe correction images
(minreplace = 1.) Minimum flat field value
(scantype = "shortscan") Scan type (shortscan|longscan)
(nscan = 1) Number of short scan lines\n
(interactive = yes) Fit overscan interactively?
(function = "legendre") Fitting function
(order = 4) Number of polynomial terms or spline pieces
(sample = "*") Sample points to fit
(naverage = 1) Number of sample points to combine
(niterate = 3 Number of rejection iterations
(low_reject = 2.5) Low sigma rejection factor
(high_reject = 2.5 High sigma rejection factor
(grow = 0.) Rejection growing radius
(mode = "ql")  

      

myiraf$/yalo_ir.dat:
exptime exptime
imagetyp imagetyp
subset IRFLTID

 

        OBJECT   object
        DARK   zero
        FLAT flat  
        BIAS   zero
        MASK   other

 

myiraf$/yalo_ir.sub

        'H' H
        'J' J
        'K' K

 You can rename the reduced data to something simple, like the following. If you do this command, make sure you don't make typos!

imren ir000729.0*.imh %ir000729.0%r%*.imh

The r*.imh data are [ZF]

 

MAKE THE MASK

Make a mask image as follows. Here we use the dome flats corrected for DOME_OFF for the mask. Note that there are very many warm pixels with the detector and about 10% of these change flux during the night. If the warm pixels change flux between the ON and OFF images, they will be flagged as bad pixels here.

The philosophy of the masks is that all pixels in a normalize image that are less than some value like 0.7 are probably bad, and will be marked as a bad pixel.

mask1.cl:
# to make the mask, use imhist and look for the limits
# first flatten the flats and remove the edge defects
#
string img
img = "ir000729.flath"
#
imdel("temp*.imh", >>& "dev$null")
imdel("mask.imh,maskdao.imh,mask.pl", >>& "dev$null")
imtrans(img,"temp1")
fmed("temp1","temp2", xwin=201, ywin=1, boundary="wrap")
imtrans("temp2","temp3")
imar(img, "/", "temp3", "mask")
imdel("temp*.imh", >>& "dev$null")
imrep mask[*,1:10] 0 lower=INDEF upper=INDEF
imrep mask[*,1020:1024] 0 lower=INDEF upper=INDEF
imrep mask[1:1,*] 0 lower=INDEF upper=INDEF
imrep mask[1021:1024,*] 0 lower=INDEF upper=INDEF
#
# now check the historgram and change the limits if needed.
#
imhist mask z1=0.4 z2=1.4 nbins=100
displ mask 1 zs- zr- z1=0.4 z2=1.4


mask2.cl
#
# good pix are 0, bad are 1 for IRAF mask
# the values 0.65 and 1.25 need to be checked on the histogram
# each time you make the mask.
#
real lll,uuu
real hist1,hist2,hist3,xjunk,histsum,nax1,nax2,npixx,ratio
lll = 0.7
uuu = 1.19
#
imhist('mask',z1=lll,z2=uuu,list+,nbin=1) | scan(xjunk,hist1)
imhist('mask',z1=INDEF,z2=lll,list+,nbin=1) | scan(xjunk,hist2)
imhist('mask',z1=uuu,z2=INDEF,list+,nbin=1) | scan(xjunk,hist3)
histsum= hist1+hist2+hist3
hsel('mask','naxis1','yes') | scan(nax1)
hsel('mask','naxis2','yes') | scan(nax2)
npixx=nax1*nax2
ratio=(hist2+hist3)/npixx
printf("Fraction rejected=%9.3f\n",ratio)
#
imhist('mask',z1=lll,z2=uuu,list+,nbin=1)
imdel("temp*.imh")
imcopy mask temp
displ mask 1 zs- zr- z1=0.4 z2=1.4
imrep("mask", lower=INDEF, upper=lll, val=-1 )
imrep("mask", lower=uuu, upper=INDEF, val=-1)
imrep("mask", lower=lll, upper=uuu, val=0)
imar mask * mask mask
imcopy mask.imh mask.pl
# make DAOPHOT mask where bad pix are 0 and good are 1
imrename mask.imh maskdao
imar maskdao - 1 maskdao
imar maskdao * -1 maskdao
#
displ mask.pl 2 zs- zr- z1=0 z2=1

You can check frames 1,2 to see if the mask looks good.

 

VIGNETTING CORRECTION

1. Make a directory called "vig" and copy the K data from r*.imh which is [ZF] data.

2. Fixpix r*.imh using

fixpix r*.imh mask=mask.pl

3. Edit out the star and galaxy.

We edit the stars out with the "b" aperture and an radius of 18. Edit out the stars first, then change the radius to 40 by typing ":rad 40", and then edit out the galaxy. You may have to use the "c" key instead if the galaxy is sitting on a large gradient. In this case you mark the lower left and upper right corners. The rectangle should be big and long in the column direction!

imedit r264 a264 aper="circular" radius=15 buf=10 wid=10
imedit r265 a265 aper="circular" radius=15 buf=10 wid=10
etc.

2. Divide the images by the first image in dither ==> b???.imh

imar a264 / a264 b264
imar a265 / a265 b265
etc.

3. You must make the b???.imh roughly equal to 1 before doing the filtering in the next steps.

Pick a statsec on the images where there is no bright star or galaxy. Figure out the minimum and maximum vignetting in the divided images. For dither=30, zlo=0.9 and zhi=1.1. These are very important! For larger dithers, the numbers must be more like 0.8 and 1.2. Run;

task normrat = home$scripts/normrat.cl
normrat b*.imh

normrat:

images = "b*.imh" input images
(statsec = "[50:600,25:1000]") Stat sec
(sigma = 2.5) sigma clip for stats
(niter = 10) interations for sigma clipping
(pre1 = "b") input prefix
(pre2 = "c") output prefix
(zlo = 0.8) Low cutoff for fmed
(zhi = 1.2) High cutoff for fmed
(xwin = 1) xwin for fmed
(ywin = 351) ywin for fmed
(outfile1 = "temp1.cl") output file for norm script1
(outfile2 = "temp2.cl") output file for norm script2
(outfile3 = "temp3.cl") output file for norm script3
(imglist1 = "tmp$tmp.8179fa")  
(mode = "ql")  

 

This will produce 3 scripts, temp1.cl, temp2.cl, and temp3.cl.
To normalize the b???.imh, do

cl < temp1.cl

You can run getsky to check that the norm looks good.

NOTE THAT THE FIRST B???.IMH IMAGE IS 1.0 AND CAN BE IGNORED.

To produce the medianed flat, run

cl < temp2.cl

which will output c*.imh images.


5. Run prog50 to force 1.0 on the good part of the chip.

Measure xlo and xhi for the good column limits. xlo is where the vignetting starts on the first row, and xhi is where the vignetting starts on the last row. All pixels to the left of this will be set to 1.0, and all pixels to the right will have the value as in the c???.imh image. The script temp3.cl has the basic commands, but YOU MUST EDIT THE XLO AND
XHI VALUES.

Note that the dithers:
1,4,7
2,6
3,5
should be pretty much the same vignetting.

Find the xlo,xhi:

yalocen c*.imh iraf- zz1=0.9 zz2=1.1

# dither 2
!$myprog/prog50 c265.imh c265.imh 730 640
# dither 3
!$myprog/prog50 c266.imh c266.imh 770 680
# dither 4
!$myprog/prog50 c267.imh c267.imh 750 660
# dither 5
!$myprog/prog50 c268.imh c268.imh 770 680
# dither 6
!$myprog/prog50 c269.imh c269.imh 730 640
# dither 7
!$myprog/prog50 c270.imh c270.imh 750 660

730 640
770 680
750 660
770 680
730 640
750 660

These values really should not be changing!

The c???.imh data now represent dither=x/dither=1 corrections. We need to divide these data into the original r*.imh data in the upper directory. To keep the bookkeepping straight, first lets make an image for dither=1.

imcopy r264 vig01
imrep vig01 1 lower=INDEF upper=INDEF

Now rename the c*.imh data to vig01,vig02, ... vig07. Make sure it looks okay as:

hsel vig*.imh $I,dither yes

6. Now correct the r*.imh data. Go to the upper directory and:

imren vig/vig*.imh .
hsel r*.imh $I,dither yes > in2

Edit in2 as:
imar r250 / vig01 f250
hedit f250 VIGCOR "Corrected for vignetting by vig01" add+ up+
imar r251 / vig02 f251
hedit f251 VIGCOR "Corrected for vignetting by vig02" add+ up+
imar r252 / vig03 f252
hedit f252 VIGCOR "Corrected for vignetting by vig03" add+ up+
imar r253 / vig04 f253
hedit f253 VIGCOR "Corrected for vignetting by vig04" add+ up+
imar r254 / vig05 f254
hedit f254 VIGCOR "Corrected for vignetting by vig05" add+ up+
etc.

 

SKY SUBTRACTION

Make inj,inh,ink files for all the SN data. These will be used to make the sky.

imdel in*
files f*.imh > in1
hsel @in1 $I,irfltid yes | grep "J" - | fields - 1 > inj
hsel @in1 $I,irfltid yes | grep "H" - | fields - 1 > inh
hsel @in1 $I,irfltid yes | grep "K" - | fields - 1 > ink

Run irsky. MAKE SURE THAT THE INSUF AND OUTSUF ARE CORRECTLY SET.

irsky:

  images = "@inh" input images
  (statsec = "[25:600,25:1000]") Stat sec
  (sigma = 2.5) sigma clip for stats
  (niter = 9) interations for sigma clipping
  (irfltid = "IRFLTID") keyword for filter
  (outimage = "Sky") Output root for sky image
  (nlow = 0) number of low pixels to reject
  (nhigh = 1) number of high pixels to reject
  (combine = "median") type of combine function
  (reject = "minmax") type of rejection
==> (insuf = "f") Root suffix for input image
==> (outsuf = "s") Root suffix for output image
  (imglist1 = "t1.jnk")  
  (mode = "al")  


You may have to play with the nhigh to reduce the print-through.

 

This program outputs a file called sub.cl which you run to do the sky subtractions.

cl < sub.cl


This is now sky subtracted data. All the data should be near 0 sky. You can check this with getsky.

task getsky = home$scripts/getsky.cl

 

FIX BAD PIXELS

For the final mosaic, you should set the bad pixels to a large number. Since saturation is 12000, 20000ADU is a good value.

imar s*.imh / maskdao s*.imh divzero=20000

If you want to fix the bad pixels for pretty pictures:

fixpix s???.imh mask=mask.pl

The data will be [BZF] now.

 

FINAL MOSAIC

The final mosaic is a piece of art, and I don't have the killer technique yet. The following does an excellent job if the night is photometric. The main problem we face is to detect and remove the warm pixels/cr's without removing flux from the stars.

The first step is to shift the data. If the seeing is >3 pix or so, use integer shifts.

We will now operate on the s*.imh images. Run:
!mv inj temp ; sed s/f/s/ temp > inj ; rm temp
!mv inh temp ; sed s/f/s/ temp > inh ; rm temp
!mv ink temp ; sed s/f/s/ temp > ink ; rm temp

del junk.dat
yalocenter @inj
!$myprog/prog48 junk.dat
cl < shift.cl
etc.

This will produce integer shift image called temp*.imh. You can modify prog48 if you want real-valued shifts but I would not recommend it.

The final combine can be made as follows.

Use stsdas.hst_calib.wfpc package and run noisemodel on your data. Converge on the read noise and scalenoise. You will see a plot with a bunch of points at the lower left and two paralllel sequences to the upper right. Fudge the read noise until it passes thought the lower left. Then fugde the scalenoise (in units of percent) until it passes through the LOWER sequence. These are the stars. The upper sequence are warm pixels.

stsdas
hst_calib
wfpc

noisemodel s111 xb=10 yb=10

Input these parameter to imcomb, remembering to convert from percent to fraction. For instance, I found:

imdel t.imh,t.pl
# H
#imcomb temp??.imh t plf=t.pl comb=ave reject=ccd lth=-100 hth=13000 \\
# gain=6.5 rdn=49 snoise=0.30 lsig=4 hsig=4
# K
#imcomb temp??.imh t plf=t.pl comb=ave reject=ccd lth=-200 hth=13000 \\
# gain=6.5 rdn=104 snoise=0.30 lsig=4 hsig=4
# J
#imcomb temp??.imh t plf=t.pl comb=ave reject=ccd lth=-50 hth=13000 \\
# gain=6.5 rdn=22 snoise=0.25 lsig=4 hsig=4
displ t.imh 1 zs- zr- z1=-20 z2=100
displ t.pl 2

When the detector had lots of warm pixels, I usedL

imdel t.imh,t.pl
# H
#imcomb temp??.imh t plf=t.pl comb=ave reject=ccd lth=-50 hth=13000 \\
# gain=6.5 rdn=72 snoise=0.60 lsig=6 hsig=5
# K
#imcomb temp??.imh t plf=t.pl comb=ave reject=ccd lth=-500 hth=13000 \\
# gain=6.5 rdn=140 snoise=0.55 lsig=7 hsig=6
# J
imcomb temp??.imh t plf=t.pl comb=ave reject=ccd lth=-50 hth=13000 \\
gain=6.5 rdn=55 snoise=0.70 lsig=7 hsig=5
displ t.imh 1 zs- zr- z1=-20 z2=100
displ t.pl 2

Move the t.imh and t.pl to J.imh J.pl etc.

If the night was not photometric, we have to estimate a scale
factor. I have not figured this out yet but it will require scaling on
the galaxy or some stars, but doing the calculation throwing out bad
pixels.

 

SOME UNRESOLVED QUESTIONS

1. Why do the warm pixels come and go through the night? Sometimes I have very clean data and other times 2-3% of the pixels are >1000ADU.

2. Why is the Sky? frame not very flat? I see maybe 1-3% variations in the non-vignetted part of the detector. The sky has beeen divided by the domes, and should be flat. There are not enough data to create a proper sky flat for division, so I am stuck with the present technique. Note that this can produce systematic errors (especially as a function of x) since the slope is ususally in the same direction in x.

The dome white spot is canted about 8deg wrt to the telescope ring. The flat is illuminated from a lamp in the middle of the ring. This can't be a very good illumination and the general slope (and vignetting problem) may be due in part to this.

3. Why does some of the data show a 2% discontinuity across the two detector halves *after* flatfielding? Is this a large variable bias?

back to top