; Copyright (c) 1998-2001 A.P. Hitchcock  All rights reserved
; ORIGINAL from CXRO, LBNL - adapted Feb-01 by aph
;+
;
;NAME:
;		XMLOAD
;
;LAST CHANGED: ----------------------------------- 12-may-02 (aph)
;
; PURPOSE:
; This function loads an XM-1 image, removes camera artifacts,
; and rotates the image to correspond with the Winview orientation.
;
; CONTENTS OF THIS FILE (in order)
;
;	LOADCMSA
;	COMPRESSPIXELS
;	NOJUNK
;	EXP_TIME
;	XMLOAD
;
; CATEGORY:
;	Image read in
;
; CALLING SEQUENCE:
;	Result = XMLOAD(FNAME)
;
; INPUTS:
; 	FNAME - name of the file (optional)
;
; KEYWORDS:
;	ZOOM
;	NOSHOW		DO NOT DISPLAY (if not set, uses LOOK in standalone mode)
;	MZPSHIFT
;	MAG
;	BINX
;	BINY
;	ENERGY
;	WAVEL
;	NOPRINT
;	VERSION
;	IMTYPE
;
; OUTPUTS:
;	No explicit outputs.  In stand alone mode, displays corrected image
;
; COMMON BLOCKS:
;	AXIS_COM	standard set of common blocks; locked into running through AXIS
;	CMSA
;
; SIDE EFFECTS: none
;
; RESTRICTIONS:
;	LOOK - not yet functioning
;
; PROCEDURE:
;	Image is read in using loadcmsa.pro
;	Camera artefacts corrected etc using nojunk.pro
;	Display in stand alone mode using look.pro
;
; MODIFICATIONS
; 8/4/95: John Heck first version
; 1-9-98: modified to return 0 if file not found by loadcmsa, Angelic
; 5-1-97: modified for filenames bigger than ymmdd999, Werner
; (21-feb-01 aph) adapting for windows environment & reading into AXIS
; (23-feb-01 aph) implement calibration and background corrections
; (12-may-02 aph) correct header name format
;-

; *******************************************************************
; -------------------------------------------------------------------
; *******************************************************************
;
;------------------------------------
; Ken Goldberg  9/13/94
; IDL procedure:  loadcmsa.pro
; this procedure is meant to read cmsa data files, from the winview
; program and display them in IDL.
;
; The CMSA format is discussed in the CMSA User's Manual, pages 181-3.
; Essentially, there are 4100 bytes of "header" preceeding the data.  The
; header is filled with specific information scattered about the 4100 bytes.
;
;  Changes: (mi) 1-3-95:
;	Usage out=loadcmsa(file,zoom=f,noshow=1,noload=1,notype=1,square=1)
;	(everything optional)
;	if the file is not found in, /cxro/xm1disk/, it is looked for in,
;	/cxro/xm1data/, then in, /cxro/mss/xm1/
;
;	zoom is a new keyword that if given, displays the image
;	using the zoomfactor f
;
;	12-9-98 (mi) added mzpx and mzpy parameter
;
;       1-22-98 (ael) fixed minor error
;
;       7-7-99 (ael) added skip_mss keyword.  when set, will not load image if found in mss
;
; (20-feb-01 aph) incorporating XM-1 read in into AXIS2000
;		changed first line due to length limit
;		got rid of all the XM-1 specific location checking
;-------------------------------------
;

function loadcmsa, fName, zoom=zoom, noload=noload, noshow=noshow, $
notype=notype,square=square,exptime=exptime,mzpx=mzpx,mzpy=mzpy,mzpz=mzpz, $
ccdz=ccdz,mag=mag,binx=binx,biny=biny,skip_mss=skip_mss,wavel=wavel, $
noprint=noprint,version=version, imtype=imtype,mzpvoltages=mzpvoltages

common cmsa, d
@axis_com

; @exp_time

if NOT KEYWORD_SET(skip_mss) then skip_mss=0
if (N_Elements(fName) EQ 0) then begin
	fName=dialog_pickfile(FILT='*.spe', GET_PATH=P)
	if (strlen(fName) EQ 0) then return,0
endif

openr, unit, fname, /GET_LUN
b=bytarr(4100)
xSize = 0
readu, unit, b
xSize = (fix(b,0,1))(0)
if !version.os EQ 'sunos' then byteorder, xSize			; byteorder  needed only for unix system
dataType = (fix(b,108,1))(0)
;      byteorder, dataType
DT=["Float","Long Integer","Integer","Unsigned Integer","String/Char", $
    "Double","Byte","Unsigned Byte"]
if not Keyword_Set(notype) then print, 'dataType = ',DT(dataType)
b1 = b(110:655)
b(110:655) = b1
ySize = (fix(b,656,1))(0)
if !version.os EQ 'sunos' then byteorder, ySize
; print the comments stored in 200-599 as 5x80 bytes.
for j = 200,599,80 do begin
      st = string(b(j:j+80))
      if  (strlen(st) GT 0) then print, st
endfor

exp_time, header=b, exptime=exptime ,noprint=noprint    ; print exposure length

if Keyword_Set(noload) then begin
   free_lun, unit
   return, 0   ; exp_time doesn't need the data
endif

if not Keyword_Set(noprint) then $
     print, FORMAT='("Image Size:          ",i4," x ",i4)', xSize, ySize

mzpx=(swap_endian(double(b,1492,1)))(0)
mzpy=(swap_endian(double(b,1500,1)))(0)
mzpz=(swap_endian(double(b,848,1)))(0)

ccdz=(swap_endian(double(b,856,1)))(0)

mag=(swap_endian(double(b,814,1)))(0)
binx=(swap_endian(fix(b,788,1)))(0)
biny=(swap_endian(fix(b,790,1)))(0)

wavel=(swap_endian(double(b,806,1)))(0)
imtype=String(BYTE((swap_endian(double(b,1508,1)))(0)))
reads,string(b(688:704)),version

mzpvoltages=fltarr(4)
mzpvoltages(0)=(swap_endian(double(b,1600,1)))(0)
mzpvoltages(1)=(swap_endian(double(b,1608,1)))(0)
mzpvoltages(2)=(swap_endian(double(b,1616,1)))(0)
mzpvoltages(3)=(swap_endian(double(b,1624,1)))(0)

; here's the manipulation which takes the data from cmsa to idl format
; print, 'data type is ' , dataType
case dataType of
    0: begin  ; float (4 byte)
         d = fltarr(xSize, ySize)
         readu, unit, d
		if !version.os EQ 'sunos' then byteorder, d, /LSWAP
       end
    1: begin  ; long integer (4 byte)
         d = lonarr(xSize, ySize)
         readu, unit, d
		if !version.os EQ 'sunos' then byteorder, d, /LSWAP
       end
    2: begin  ; integer (2 byte)
         d = intarr(xSize, ySize)
         readu, unit, d
		if !version.os EQ 'sunos' then byteorder, d
       end
    3: begin  ; unsigned integer (2 bytes)
         d = intarr(xSize, ySize)
         readu, unit, d
		if !version.os EQ 'sunos' then byteorder, d
         d = long(d) + long(65536)*(d LT 0)
       end
    4: begin ; string/char (1 byte)     {not used}
         d = bytarr(xSize, ySize)
         readu, unit, d
         d = string(d)
       end
    5: begin  ; double (8 byte)         {actually, not implemented}
         d = dblarr(xSize, ySize)
         readu, unit, d
       end
    6: begin ; byte (1 byte)            {not used}
         d = bytarr(xSize, ySize)
         readu, unit, d
       end
    7: begin ; unsigned byte (1 byte)   {not used}
         d = bytarr(xSize, ySize)
         readu, unit, d
         print, "Data is in Unsigned Byte format, which I will treat as byte"
       end
    else: print, 'dataType error!'
endcase
free_lun, unit
NP = xSize < ySize
if ((xSize NE ySize) and Keyword_Set(square)) then begin $
    print, 'Image has been cropped square.'
    d = float( d(0:NP-1,0:NP-1) )  ; make sure it's square by cropping excess
    end
  	if (not Keyword_Set(noshow)) then begin
	    if (Keyword_set( zoom )) then begin
		look,zm(d,zoom)
  	endif else begin
		window, xsize = xSize+25, ysize = ySize+25
		tvscl, d
  	endelse
endif
if not Keyword_Set(noprint) then print, 'Loaded file: ',fname
return, d
end

; *******************************************************************
; -------------------------------------------------------------------
; *******************************************************************

; Angelic replaced nohot.pro with code to read file name from the calibration
; file which contains information about the current CCD's hot spots.
; Also included function for converting 1-D pixel values to account for binning
; 7-6-99
;
; Werner and Ryohei commented out the dark current limit, as it makes no
; sense with the now defect camera, also should include binning!
;
; John Heck, 3/29/96 (changed dark-test to 6025 + 100 exptime Werner 4-3-96)
; and restricted output image to > 0.1
;
; This procedure creates a new image, with the CCD junk subtracted, from an
; existing image.
;
; It also replaces the camera's hot spot with averaged values
;
; Input old image (icounts), exposure time, and filename so the correct
; calibration set can be chosen.
;
; Output is a new image, 'iphotons'.

; -----------

FUNCTION compresspixels, pixarr, bin
;equation for converting 1-D pixels when 2-D image has binning
RETURN, (pixarr/(1024*bin))*1024/bin + (pixarr mod 1024)/bin
END

; -----------

pro nojunk,icounts,iphotons,exptime,fname,COUNTS_ONLY=counts_only,noprint=noprint

@axis_com

icounts=float(icounts)
ims=size(icounts)
imsiz=ims(1)
;print,'Size of Image...',imsiz

if NOT( KEYWORD_SET( COUNTS_ONLY ) ) then begin
if (imsiz ne 1024) and (imsiz ne 512) and (imsiz ne 256) then begin
	iphotons=1.*icounts
	print,"Not calibrated, counts replaced"
	return
endif
endif

; use sunos as the trigger for identifying use of XM1 Unix computers
if !version.os EQ 'sunos' then IDL_XM1_DIR = !IDL_XM1_DIR else IDL_XM1_DIR = ''

; Determine appropriate calibration set and read parameters...
calibFileXM1 = IDL_XM1_DIR + 'pros/calibrat.dat'
nsets=0
nlines=0
calbeg=0		  ; 60000000  date variable ??
calend=0          ; 60000000
c3=''
offav=0.
gfname='bgnd'
gmin=0.
gmf=1.
hot_file='hot'
set = 0
getout = 0
on_ioerror, SKIP_CALIB
calibFile = dialog_pickfile(title = 'XM-1 calibration', file='calibFileXM1', path = DefPath)
openr,unit,calibFile,/GET_LUN
readf,unit,nsets
readf,unit,nlines
while getout eq 0 do begin
   readf,unit,calbeg
   readf,unit,calend
   readf,unit,c3
   readf,unit,c4
   readf,unit,offav
   readf,unit,gfname
   readf,unit,gmin
   readf,unit,gmf
   readf,unit,c9
   readf,unit,c10
   readf,unit,hot_file
   readf,unit,c12
   readf,unit,c13
   readf,unit,c14
   readf,unit,c15
   readf,unit,c16
   set = set + 1
   if calend lt 0 then getout=1
   if fname gt calbeg and fname lt calend then getout=1
endwhile
close,unit
FREE_LUN,unit

SKIP_CALIB:
if ( not Keyword_Set(noprint)) then begin
if (calend lt 0) then print,strcompress('Using calibration set # ' + $
   string(set) + ' for images after ' + string(calbeg)) else $
   print,strcompress('Using calibration set # ' + string(set) + $
   ' for images between ' + string(calbeg) + ' and ' + string(calend))
 endif

;print,'nsets ',nsets
;print,'nlines ',nlines
;print,'calbeg ',calbeg
;print,'calend ',calend
;print,'offav ',offav
;print,'gfname ',gfname
;print,'gmin ',gmin
;print,'gmf ',gmf

; Read background data
gain = make_array(imsiz,imsiz,/byte,value=1)
on_ioerror, SKIP_Background
if IDL_XM1_DIR NE '' then t = IDL_XM1_DIR +'pros/' else t = ''
BgndFileXM1 = strcompress(t + gfname + '.' + $
   string(imsiz),/remove_all)
BgndFile = dialog_pickfile(title='XM-1 background file', file='BgndFileXM1', path = DefPath)
openr,unit,BgndFile,/GET_LUN
readu,unit,gain
close,unit
FREE_LUN,unit

SKIP_Background:

; --------- Make iphotons ---------------
; smooth the image (!!)
bg = median(icounts(imsiz-9:imsiz-5,imsiz-9:imsiz-5))

if NOT( KEYWORD_SET( COUNTS_ONLY ) ) then begin
	iphotons = ((float(icounts) - bg) / (float(gain) + gmin) * gmf / 24.4) > 0.1
	if not Keyword_Set(noprint) then print,'Divided by 24.4 (counts/photon @ 2.4 nm) and made positive'
endif else begin
	iphotons = (((float(icounts) - bg) / (float(gain) + gmin) * gmf ) + bg) > 0.1
endelse

; Remove hot spot of image
if (hot_file NE '0') then begin
    iphotons=rotate(iphotons, 7)
    bin=1024/imsiz
	GOTO, SKIP_Hot
	if IDL_XM1_DIR NE '' then t = IDL_XM1_DIR +'pros/' else t = ''
   	HotFileXM1 =strcompress(t+hot_file, /remove_all)
	HotFile = dialog_pickfile(file=HotFileXM1)
;print, 'opening hot file'
    openr, unit,HotFile, /GET_LUN
    numspots=0
    numdefect=0L
    numrim=0L
    readf, unit, numspots
    for i=1, numspots do begin
;read in pixels for defect and rim
        readf, unit, numdefect
        defect=lonarr(numdefect)
        readf, unit, defect
        readf, unit, numrim
        rim=lonarr(numrim)
        readf, unit, rim
;change pixels if image was binned
        defect=compresspixels(defect, bin)
        rim=compresspixels(rim, bin)
;replace each pixel in defect with a value near the median of the rim
;some pixels will be replaced multiple times if binning occurs, not a problem
;seed variable in randomn call is intentionally undefined
        iphotons(defect) = median(iphotons(rim)) + stdev(iphotons(rim)) * randomn(seed, numdefect)
    endfor
    free_lun, unit

    SKIP_Hot:

    iphotons=rotate(iphotons, 7)
endif


finito:
end

; *******************************************************************
; -------------------------------------------------------------------
; *******************************************************************

pro exp_time, header=header,exptime=tt,noprint=noprint
;
; ------------------------------------
; Ken Goldberg  9/13/94
; IDL procedure:  exp_time.pro
;
; This procedure reads the header to a cmsa file and prints the
; exposure time on the screen.
;
; The CMSA format is discussed in the CMSA User's Manual, pages 181-3.
; Essentially, there are 4100 bytes of "header" preceeding the data.  The
; header is filled with specific information scattered about the 4100 bytes.
; -------------------------------------
;


  if Keyword_Set(header) then begin     ; use the header supplied
    t = (long(header, 660, 1))(0)
    byteorder, t, /LSWAP
   if not Keyword_Set(noprint) then begin
    if (t LT 1000) then print, FORMAT='("Exposure Time: ",i8," msec")', t $
                   else print, FORMAT='("Exposure Time: ",f8.3," sec")',t/1000.
    endif
    endif else begin
    q = loadcmsa(/noload, /noshow, /notype)
  endelse
tt = t/1000.

return
end

; *******************************************************************
; -------------------------------------------------------------------
; *******************************************************************

function xmload,fname,zoom=zoom,noshow=noshow,mzpshift=mzpshift, $
     mag=mag,binx=binx,biny=biny,energy=energy,wavel=wavel,noprint=noprint, $
     version=version,imtype=imtype

@axis_com

icounts = loadcmsa(fname,nosh=1,exptime=exptime,notype=1,mzpx=mzpx,mzpy=mzpy,mag=mag,binx=binx,biny=biny,wavel=wavel,noprint=noprint,version=version,imtype=imtype)

;check if fname found by function loadcmsa, if not, then return
s=SIZE(icounts)
IF (s(0) EQ 0) THEN RETURN, 0

;print,'Exposure Time: ',exptime
;read,'Enter Exposure Time: ',exptime

; (aph) assume that fname is meant to be name only
; typically 10206123.spe YMMDDnnn
; so convert the path&name to name only
t=ax_name(fname)
fname = t(1)+'.' + t(2)

;print,'Removing CCD junk... '
namebase=1000*long(strmid(fname,0,5))
 nojunk,icounts,iphotons,exptime,namebase,noprint=noprint

;print,'Rotating image... '
iphotons = rotate(iphotons,7)

if (not keyword_set(noshow)) then begin
;   print, 'Looking at image... '
   if (keyword_set(zoom)) then look,zm(iphotons,zoom) else look,iphotons
endif

if( keyword_set(mzpshift) ) then begin
	dx=fix(mzpx*mag/24e-6/binx+.5)
	dy=fix(mzpy*mag/24e-6/biny+.5)
	iphotons=shift(iphotons,dx,dy)
	   if not Keyword_Set(noprint) then print, FORMAT= '("Shifted image by: ",i4,i4)',dx,dy
endif

energy=1239.852178/(wavel*1e9)

return,iphotons

finito:
end
