;; 		scatter.pro 
;;
;;      creates CBT/MT2 scatter plots  
;;    from CISSCAL calibrated data files 
;;
;;    8/3/01
;;    modified 8/8/01 



device, decomposed=0

; open and read raw data files

xsize = 1024
ysize = 1024
cb2image = fltarr(xsize, ysize)
mt2image = fltarr(xsize, ysize)

openu, 1, dialog_pickfile(title = 'Select CB2 Image', $
			Filter = '*cal.raw')
readu, 1, cb2image
close, 1

openu, 1, dialog_pickfile(title = 'Select MT2 Image', $
			Filter = '*cal.raw')
readu, 1, mt2image
close, 1  


; correct image orientation (**edit accordingly**)

cb2image = rotate(cb2image, 5)
mt2image = rotate(mt2image, 5)


; display MT2 image (scaled)

mydevice = !d.name

mt2imagehalf = rebin(mt2image, xsize/2, ysize/2)
window, /free, xsize = xsize/2, ysize = ysize/2, $
		title = 'Calibrated MT2 Image'
tvscl, mt2imagehalf 
xyouts, /normal, 0.05, 0.95, 'Click on Feature to Isolate', $
			charsize = 2   

; locate feature

print, 'Click on feature to find location'
x = 0
y = 0
cursor, xcoord, ycoord, /wait
xposition = 0
yposition = 0

print, 'xcoord = ', xcoord  
print, 'ycoord = ', ycoord 
print, '________________'

xposition = round(xcoord*xsize)
yposition = round(ycoord*ysize) 
print, 'x position = ', xposition
print, 'y position = ', yposition


; create feature-only subarrays

xminmark = 0
xmaxmark = 0
yminmark = 0
ymaxmark = 0

xminmark = xposition - 50 
xmaxmark = xposition + 49
yminmark = yposition - 50
ymaxmark = yposition + 49

print, ''
print, 'xminmark = ', xminmark
print, 'xmaxmark = ', xmaxmark
print, 'yminmark = ', yminmark
print, 'ymaxmark = ', ymaxmark

cb2feature = cb2image[xminmark:xmaxmark, yminmark:ymaxmark]
mt2feature = mt2image[xminmark:xmaxmark, yminmark:ymaxmark]


; display feature-only subarrays  
 
!P.multi = [0, 2, 1]
window, /free, xsize = 210, ysize = 110, title = '      CB2 | MT2'
tvscl, cb2feature 
tvscl, mt2feature, 1 
!P.multi = 0



; realign images

window, /free, title = 'CB2 Contour for Realignment'
contour, cb2feature, /fill, nlevels = 20, xtitle = 'pixels', $
	ytitle = 'pixels', title = 'CB2'

window, /free, title = 'MT2 Contour for Realignment'
contour, mt2feature, /fill, nlevels = 20, xtitle = 'pixels', $
	ytitle = 'pixels', title = 'MT2'

mvright = 0
mvleft = 0
mvup = 0
mvdown = 0

print, ' '
print, '_____________________________________________________'
read, mvright, $
        prompt = 'Move MT2 feature RIGHT how many pixels?  '
print, '_____________________________________________________'


print, ' '
print, '_____________________________________________________'
read, mvleft, $
        prompt = 'Move MT2 feature LEFT how many pixels?  '
print, '_____________________________________________________'


print, ' '
print, '_____________________________________________________'
read, mvup, $
        prompt = 'Move MT2 feature UP how many pixels?  '
print, '_____________________________________________________'


print, ' '
print, '_____________________________________________________'
read, mvdown, $
        prompt = 'Move MT2 feature DOWN how many pixels?  '
print, '_____________________________________________________'

xminmark = xminmark - mvright + mvleft
xmaxmark = xmaxmark - mvright + mvleft  
yminmark = yminmark - mvup + mvdown  
ymaxmark = ymaxmark - mvup + mvdown  

mt2feature = mt2image[xminmark:xmaxmark, yminmark:ymaxmark]




window, /free, title = 'Realigned MT2 Contour'
contour, mt2feature, /fill, nlevels = 20, xtitle = 'pixels', $
        ytitle = 'pixels', title = 'MT2'





; Create scatter plots (for window and PS output)   

print, ' '
print, '_____________________________________________________'
date = ''
read, date, $
        prompt = 'Enter date (ex: JAN 01)  '
print, '_____________________________________________________'


print, ' '
print, '_____________________________________________________'
latitude = ''
read, latitude, $
        prompt = 'Enter latitude (ex: 23 DEG. N LAT.)  '
print, '_____________________________________________________'

print, ' '
print, '_____________________________________________________'
plottitle = ''
read, plottitle, $
        prompt = 'Enter Plot Title: '
print, '_____________________________________________________'


print, ' '
print, '_____________________________________________________'
psfilename = ''
read, psfilename, $
        prompt = 'Enter PostScript File Name (ex: 0101_23.5_f.ps)  '
print, '_____________________________________________________'


window, /free, title = 'Calibrated Scatter Plot'
plot, cb2feature, mt2feature, psym = 1, xrange = [0, 1], $
	yrange = [0, 1], xtitle = 'I/F 751 nm', $
	ytitle = 'I/F 727 nm', title = string(plottitle), $ 
	xthick = 2.0, ythick = 2.0, $
	thick = 1.5, charsize = 1.5, charthick = 1.5 
xyouts, 0.05, 0.9, string(date), charsize = 1.5, charthick = 1.5
xyouts, 0.05, 0.85, string(latitude), charsize = 1.5, charthick = 1.5

set_plot, 'ps'   
device, filename = string(psfilename)
plot, cb2feature, mt2feature, psym = 1, xrange = [0, 1], $
        yrange = [0, 1], xtitle = 'I/F 751 nm', $
	ytitle = 'I/F 727 nm', title = string(plottitle), $
	xthick = 2.0, ythick = 2.0, $
	thick = 1.5, charsize = 1.5, charthick = 1.5
xyouts, 0.05, 0.9, string(date), charsize = 1.5, charthick = 1.5
xyouts, 0.05, 0.85, string(latitude), charsize = 1.5, charthick = 1.5
device, /close
set_plot, mydevice 

           
end