;
 H band C/2002 C2 I-K NSFCam Data - February 23, 2002 UT
; Will make postscript file.  It will NOT output to xwindow.

; enter the name of the postscript file
psfile = 'IK_H_flat.ps'

; Change plotting scale if needed
min = -10 & max = 10000

; Read in the images?
readin = 1

; Flatten then images?
flatten = 1

; Set the platescale
platescale = 0.3

; NOTE: Data is read in with E to the left and N at the bottom

if readin eq 1 then h40a = readfits('/scrs1/cew/feb24/data0040.a')
if readin eq 1 then h41b = readfits('/scrs1/cew/feb24/data0041.b')
if readin eq 1 then h42a = readfits('/scrs1/cew/feb24/data0042.a')
if readin eq 1 then h43b = readfits('/scrs1/cew/feb24/data0043.b')
if readin eq 1 then h44a = readfits('/scrs1/cew/feb24/data0044.a')
if readin eq 1 then h45b = readfits('/scrs1/cew/feb24/data0045.b')
if readin eq 1 then h46a = readfits('/scrs1/cew/feb24/data0046.a')
if readin eq 1 then h47b = readfits('/scrs1/cew/feb24/data0047.b')
if readin eq 1 then h48a = readfits('/scrs1/cew/feb24/data0048.a')
if readin eq 1 then h49b = readfits('/scrs1/cew/feb24/data0049.b')
if readin eq 1 then h50a = readfits('/scrs1/cew/feb24/data0050.a')
if readin eq 1 then h51b = readfits('/scrs1/cew/feb24/data0051.b')

; Subtract object and sky
h0o = h40a-h41b
h1o = h42a-h43b
h2o = h44a-h45b
h3o = h46a-h47b
h4o = h48a-h49b
h5o = h50a-h51b

; Flip the images in N-S so N is on top
h0 = reverse(h0o,2)
h1 = reverse(h1o,2)
h2 = reverse(h2o,2)
h3 = reverse(h3o,2)
h4 = reverse(h4o,2)
h5 = reverse(h5o,2)

; Flatten if selected
if flatten eq 1 then restore,'/scrs1/deh/feb23/flats.idl'

if flatten eq 1 then h0 = h0/hsub_flat
if flatten eq 1 then h1 = h1/hsub_flat
if flatten eq 1 then h2 = h2/hsub_flat
if flatten eq 1 then h3 = h3/hsub_flat
if flatten eq 1 then h4 = h4/hsub_flat
if flatten eq 1 then h5 = h5/hsub_flat

; Try and make a badpixel mask
;hmask = fltarr(256,256)

;h0bad = where(h0 gt 1200)
;h1bad = where(h1 gt 1200)
;h2bad = where(h2 gt 1200)
;h3bad = where(h3 gt 1200)
;h4bad = where(h4 gt 1200)
;h5bad = where(h5 gt 1200)

;hmask(h0bad) = 1.0
;hmask(h1bad) = 1.0
;hmask(h2bad) = 1.0
;hmask(h3bad) = 1.0
;hmask(h4bad) = 1.0
;hmask(h5bad) = 1.0

; Set some plot parameters - This looks best for postscript
!x.thick=5
!y.thick=5
!p.thick=4
!p.charsize=1.5
!x.title = '!17R. A. Offset (arcsec)'
!y.title = '!17Dec. Offset (arcsec)'
!x.style=1
!y.style=1
!x.ticklen=0.04
!y.ticklen=0.04

; The center location of the comet
x0 = 148. & y0 = 125. & x1 = 123. & y1 = 151. & x2 = 123. & y2 = 97. 
x3 = 177. & y3= 96. & x4 = 180. & y4 = 150. & x5 = 154. & y5 = 122.

; Set pixel scale of array in arcsecond
nband = size(h0)
xa = platescale*findgen(256)
ya = platescale*findgen(256)

; Set up PS parameters
set_plot,'ps'
device,/portrait,/inches,filename=psfile,xoffset=1,yoffset=2,xsize=6.5,ysize=8


; Multiple images on one sheet
!p.multi = [0,2,3]

; START PLOTTING
x = x0*platescale - xa & y = -y0*platescale + ya
dispcont3,bytscl(h0,min=min,max=max),x,y,10000,aspect=1
xyouts,x(0)-2,y(0)+54*platescale,'C/2002 C1 IK 02/24/2002 UT',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+40*platescale,'NSFCam at NASA IRTF',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+27*platescale,'H,0.3"/pix, itime=1s*15,AM=2.584',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+13.5*platescale,'90" N-Nod, center = ',color=255,charsize=0.75
xyouts,x(0)-135*platescale,y(0)+13.5*platescale,fix(x1),color=255,charsize=0.75
xyouts,x(0)-176*platescale,y(0)+13.5*platescale,fix(y1),color=255,charsize=0.75
xyouts,x(0)-(x(0)-x(nband(1)-1))/2.,y(nband(2)-1)-4.,'N',color=255,charsize=0.75
xyouts,x(0)-2.,y(0)+(y(nband(2)-1)-y(0))/2.,'E',color=255,charsize=0.75

x = x1*platescale - xa & y = -y1*platescale + ya
dispcont3,bytscl(h1,min=min,max=max),x,y,10000,aspect=1
xyouts,x(0)-2,y(0)+54*platescale,'C/2002 C1 IK 02/24/2002 UT',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+40*platescale,'NSFCam at NASA IRTF',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+27*platescale,'H,0.3"/pix, itime=1s*15,AM=2.584',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+13.5*platescale,'90" N-Nod, center = ',color=255,charsize=0.75
xyouts,x(0)-135*platescale,y(0)+13.5*platescale,fix(x1),color=255,charsize=0.75
xyouts,x(0)-176*platescale,y(0)+13.5*platescale,fix(y1),color=255,charsize=0.75
xyouts,x(0)-(x(0)-x(nband(1)-1))/2.,y(nband(2)-1)-4.,'N',color=255,charsize=0.75
xyouts,x(0)-2.,y(0)+(y(nband(2)-1)-y(0))/2.,'E',color=255,charsize=0.75


x = x2*platescale - xa & y = -y2*platescale + ya
dispcont3,bytscl(h2,min=min,max=max),x,y,10000,aspect=1
xyouts,x(0)-2,y(0)+54*platescale,'C/2002 C1 IK 02/24/2002 UT',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+40*platescale,'NSFCam at NASA IRTF',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+27*platescale,'H,0.3"/pix, itime=1s*15,AM=2.584',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+13.5*platescale,'90" N-Nod, center = ',color=255,charsize=0.75
xyouts,x(0)-135*platescale,y(0)+13.5*platescale,fix(x1),color=255,charsize=0.75
xyouts,x(0)-176*platescale,y(0)+13.5*platescale,fix(y1),color=255,charsize=0.75
xyouts,x(0)-(x(0)-x(nband(1)-1))/2.,y(nband(2)-1)-4.,'N',color=255,charsize=0.75
xyouts,x(0)-2.,y(0)+(y(nband(2)-1)-y(0))/2.,'E',color=255,charsize=0.75

x = x3*platescale - xa & y = -y3*platescale + ya
dispcont3,bytscl(h3,min=min,max=max),x,y,10000,aspect=1
xyouts,x(0)-2,y(0)+54*platescale,'C/2002 C1 IK 02/24/2002 UT',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+40*platescale,'NSFCam at NASA IRTF',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+27*platescale,'H,0.3"/pix, itime=1s*15,AM=2.584',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+13.5*platescale,'90" N-Nod, center = ',color=255,charsize=0.75
xyouts,x(0)-135*platescale,y(0)+13.5*platescale,fix(x1),color=255,charsize=0.75
xyouts,x(0)-176*platescale,y(0)+13.5*platescale,fix(y1),color=255,charsize=0.75
xyouts,x(0)-(x(0)-x(nband(1)-1))/2.,y(nband(2)-1)-4.,'N',color=255,charsize=0.75
xyouts,x(0)-2.,y(0)+(y(nband(2)-1)-y(0))/2.,'E',color=255,charsize=0.75

x = x4*platescale - xa & y = -y4*platescale + ya
dispcont3,bytscl(h4,min=min,max=max),x,y,10000,aspect=1
xyouts,x(0)-2,y(0)+54*platescale,'C/2002 C1 IK 02/24/2002 UT',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+40*platescale,'NSFCam at NASA IRTF',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+27*platescale,'H,0.3"/pix, itime=1s*15,AM=2.584',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+13.5*platescale,'90" N-Nod, center = ',color=255,charsize=0.75
xyouts,x(0)-135*platescale,y(0)+13.5*platescale,fix(x1),color=255,charsize=0.75
xyouts,x(0)-176*platescale,y(0)+13.5*platescale,fix(y1),color=255,charsize=0.75
xyouts,x(0)-(x(0)-x(nband(1)-1))/2.,y(nband(2)-1)-4.,'N',color=255,charsize=0.75
xyouts,x(0)-2.,y(0)+(y(nband(2)-1)-y(0))/2.,'E',color=255,charsize=0.75

x = x5*platescale - xa & y = -y5*platescale + ya
dispcont3,bytscl(h5,min=min,max=max),x,y,10000,aspect=1
xyouts,x(0)-2,y(0)+54*platescale,'C/2002 C1 IK 02/24/2002 UT',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+40*platescale,'NSFCam at NASA IRTF',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+27*platescale,'H,0.3"/pix, itime=1s*15,AM=2.584',color=255,charsize=0.75
xyouts,x(0)-2,y(0)+13.5*platescale,'90" N-Nod, center = ',color=255,charsize=0.75
xyouts,x(0)-135*platescale,y(0)+13.5*platescale,fix(x1),color=255,charsize=0.75
xyouts,x(0)-176*platescale,y(0)+13.5*platescale,fix(y1),color=255,charsize=0.75
xyouts,x(0)-(x(0)-x(nband(1)-1))/2.,y(nband(2)-1)-4.,'N',color=255,charsize=0.75
xyouts,x(0)-2.,y(0)+(y(nband(2)-1)-y(0))/2.,'E',color=255,charsize=0.75

device,/close

set_plot,'x'

!p.multi=0

; To find center

;window,0,xsize=500,ysize=500
;disp,h0,2,min=-100,max=30000,/c
;disp,h1,2,min=-100,max=30000,/c
;disp,h2,2,min=-100,max=30000,/c
;disp,h3,2,min=-100,max=30000,/c
;disp,h4,2,min=-100,max=30000,/c
;disp,h5,2,min=-100,max=30000,/c

end