; SEQ.PRO
;
; Originally written to show images from Polar Camera initial field test
; at Ryning Farm Observatory
; Modified 1993/1/29 by DPS to handle images from image tree structure
; on /jasper/cnsr3_data1
; Modified 1993/5/11 from SHOW.PRO to display images from both cameras
; in the order they were acquired.
;
; create byte-valued input parameters
cam=0B & filt=cam
day=0 & month=day
mask=BYTARR(256,256)
; load various string-valued arrays needed
initialize,filters,modes,mpaths,firstchars
s=''
!ORDER=1
; set colour table 
loadct,0,/silent
; uncover screen where text is to be printed
IF !D.NAME EQ 'X' THEN WDELETE,0 

; Warn user that individual data files are needed.
PRINT,'This program reads individual raw data files, rather than .ray or .bgs files.'
PRINT,'If you don''t have these files on disk, stop now and load them.'
STOP

; prompt user and get input
PRINT,'Enter the starting year, month, day, UT hour, minute, and second desired'
READ,year,month,day,hr,min,sec
resp=''
Read,'Subtract dark frame before displaying images? ([Y]/N) ',resp
nodk=StrUpCase(resp) EQ 'N'
; ensure a 4-digit year for display; make reasonable assumptions
IF year LT 1900 THEN BEGIN
    IF year LT 50 THEN year=year+2000
    IF year GT 50 THEN year=year+1900
ENDIF

; create title for image plot
;wtitle=STRING(year,months(month-1),day $
;	     ,FORMAT='("Polar Camera",2X,I4,1X,A3,1X,I0.2)')
wtitle="Polar Camera  "+ymd2date(year,month,day,format='y$ n$ d$')
PRINT,wtitle

; define mask for images, with fixed horizon at 150 superpixels 
prad=150
d256=SHIFT(dist(256),128,128) & base=WHERE(d256 GT 150)
wprad=WHERE(d256 LE prad) & mask(wprad)=1

; make lists of available files of desired type
rdilist,0,year,month,day,hr,min,sec,n0,t0,files0
rdilist,1,year,month,day,hr,min,sec,n1,t1,files1

; identify and remove setup images (25 ms exposure time)
wncal=WHERE(STRPOS(files0,'cal') EQ -1,n0)
;wn25ms=WHERE(STRMID(files0,STRLEN(files0)-7,3) NE '000',n0)
t0=t0(wncal)
files0=files0(wncal)
wncal=WHERE(STRPOS(files1,'cal') EQ -1,n1)
;wn25ms=WHERE(STRMID(files1,STRLEN(files1)-7,3) NE '000',n1)
t1=t1(wncal)
files1=files1(wncal)

; zip lists together
IF (n0 EQ n1) THEN BEGIN	; equal numbers in both lists
    IF (MIN(ABS(t0-t1)) LT .5) THEN BEGIN
        t=TRANSPOSE([TRANSPOSE(t0),TRANSPOSE(t1)])
        files=TRANSPOSE([TRANSPOSE(files0),TRANSPOSE(files1)])
    ENDIF
ENDIF ELSE BEGIN		; lists have different lengths
    t=FLTARR(MAX([n0,n1]),2)
    files=STRARR(MAX([n0,n1]),2)
    IF n0 GT n1 THEN BEGIN	; more from Camera 0
        t(*,0)=t0
        files(*,0)=files0
        FOR i=0,n1-1 DO BEGIN	; find best match for each Cam 1 image
            d=MIN(ABS(t0-t1(i)),p)
            IF d LT 0.5 THEN BEGIN
                t(p,1)=t1(i)
                files(p,1)=files1(i)
            ENDIF
        ENDFOR
    ENDIF ELSE BEGIN		; more from Cam 1
        t(*,1)=t1
        files(*,1)=files1
        FOR i=0,n0-1 DO BEGIN   ; find best match for each Cam 0 image
            d=MIN(ABS(t1-t0(i)),p)
            IF d LT 0.5 THEN BEGIN
                t(p,0)=t0(i)
                files(p,0)=files0(i)
            ENDIF
        ENDFOR
    ENDELSE
ENDELSE
ext_list=extension(files)

; read in or generate required dark frames
dkt=0
dkim=INTARR(1,512)
dkpair=INTARR(256,512)
FOR i=0,n0-1 DO BEGIN
    IF STRPOS(files0(i),'_No') EQ -1 THEN BEGIN
        exp0=FIX(STRMID(files0(i),STRLEN(files0(i))-7,3))
        we=WHERE(exp0 EQ dkt,nwe)
        IF nwe EQ 0 THEN BEGIN	; haven't got this exposure time yet
            rdmeandk,0,exp0,dk0
            rdmeandk,1,exp0,dk1
            dkpair(*,0:255)=dk0
            dkpair(*,256:511)=dk1
            dkt=[dkt,exp0]
            dkim=[dkim,dkpair]
        ENDIF
    ENDIF
ENDFOR
FOR i=0,n1-1 DO BEGIN
    IF STRPOS(files1(i),'_No') EQ -1 THEN BEGIN
        exp1=FIX(STRMID(files1(i),STRLEN(files1(i))-7,3))
        we=WHERE(exp1 EQ dkt,nwe)
        IF nwe EQ 0 THEN BEGIN	; haven't got this exposure time yet
            rdmeandk,0,exp1,dk0
            rdmeandk,1,exp1,dk1
            dkpair(*,0:255)=dk0
            dkpair(*,256:511)=dk1
            dkt=[dkt,exp1]
            dkim=[dkim,dkpair]
        ENDIF
    ENDIF
ENDFOR
dkt=dkt(1:*)		; drop dummy entry at start of array
dkim=dkim(1:*,*)	; drop dummy column ...

; create and label image display window
X_charheight=10		; try to keep scalings that apply in X windows
X_charwidth=8
barwidth=16
ninrow=4
nincol=3
ninwin=ninrow*nincol
wwidth=ninrow*256+barwidth
wheight=nincol*256+4*X_charheight
IF (!D.NAME EQ 'X') OR (!D.NAME EQ 'WIN') THEN $
    WINDOW,0,XSIZE=wwidth,YSIZE=wheight $
    ELSE DEVICE,/COLOR,BITS_PER_PIXEL=8	; PostScript
csz=((1.*X_charheight)/wheight)/((1.*!D.Y_CH_SIZE)/!D.Y_VSIZE)	; keep 'X' char size
XYOUTS,scl(res=120.,(wwidth-barwidth)/2) $
      ,scl(res=120.,wheight-3*X_charheight),wtitle $
      ,/DEVICE,ALIGNMENT=0.5,CHARSIZE=2.0*csz,WIDTH=twidthn
twidthd=CONVERT_COORD([twidthn,0.],/NORMAL,/TO_DEVICE)
twidthd=twidthd(0)			; we only need the x-component
; get user's name and display it
IF (STRLOWCASE(STRMID(!VERSION.OS,0,3)) NE 'win') $
    THEN SPAWN,'whoami',username $
    ELSE username='CNSR3'
XYOUTS,scl(res=120.,wwidth-X_charwidth),scl(res=120.,wheight-3.9*X_charheight) $
      ,'User '+STRUPCASE(username),ALIGNMENT=1.0,/DEVICE,CHARSIZE=csz
; get date/time and display them
dt_tm_brk,systime(),date,time
XYOUTS,scl(res=120.,wwidth-X_charwidth),scl(res=120.,wheight-2.5*X_charheight) $
      ,time+' UT',ALIGNMENT=1.0,/DEVICE,CHARSIZE=csz
XYOUTS,scl(res=120.,wwidth-X_charwidth),scl(res=120.,wheight-1.1*X_charheight) $
      ,date,ALIGNMENT=1.0,/DEVICE,CHARSIZE=csz
; create and display color bar
bar=REBIN(BINDGEN(256,barwidth),nincol*256,barwidth)
bar=ROTATE(bar,3)
mtv,BYTSCL(bar,TOP=!D.N_COLORS-1),ninrow*256,0,res=120.

; read in each image, decode the header, subtract dark frame, rebin to 
; 256x256 if necessary, scale to 8 bits and display
ip=-1
mid=-1
itime=BYTARR(4)
FOR i=0,MIN([n0,n1])-1 DO BEGIN
    ti=t(i,0)
    p=WHERE(ABS(t(*,1)-ti),nm)
    IF nm GT -1 THEN p=[i,p(0)] ELSE p=[i,i]
    FOR j=0,1 DO BEGIN
        fnlen=STRLEN(files(i,j))
        IF fnlen GT 0 THEN BEGIN
            IF StrPos(files(i,j),'__No') EQ -1 THEN BEGIN
                rdkimg,files(i,j),h,im
                KIh=gethd(h)
                FOR k=0,3 DO itime(k)=KIh.misc.tm.(k)
                ENDIF ELSE BEGIN
                mk,im           ; make a 256x256 array of 0B
                ENDELSE
; if image isn't 256x256, rebin it
            im=REBIN(im,256,256,/SAMPLE)
; ensure counts are actually summed if a 512x512 image is compressed
            cfac=512./(256*KIh.exp.sbin)
            rfac=512./(256*KIh.exp.pbin)
            im=im*cfac*rfac
; subtract appropriate dark frame
            expt=FIX(KIh.misc.exp_scale*KIh.exp.exposure/1000.)
            dp=WHERE(dkt EQ expt)
; if this is a dark frame, store it in the dkim array
	    IF STRPOS(files(i,j),'dk') NE -1 THEN $
		dkim(256*dp(0),256*j)=MEDIAN(im,3)
            IF NOT nodk THEN im=FIX((im-dkim(256*dp(0):256*dp(0)+255,256*j:256*j+255))>0)
; identify useful min and max values for byte-scaling
            ih=HISTOGRAM(im(wprad),MIN=0,BIN=1)	; histogram of exposed image
            cih=ih	; next line computes cumulative pixel value pdf
            FOR jj=1,N_ELEMENTS(cih)-1 DO cih(jj)=cih(jj)+cih(jj-1)
            cih=cih/FLOAT(MAX(cih))		; normalize to unity
            lmaxi=MAX(WHERE(cih LE .999))	; top 0.1% of pixels saturate
; if viewing sky image, take dark 'skirt' as min value
            exptype=STRMID(files(i,j),fnlen-12,1)
            shut=(exptype EQ 'b') OR (exptype EQ 'd')	; bias or dark frame
            IF NOT shut THEN BEGIN
                bh=HISTOGRAM(im(base),MIN=0,BIN=1)
                maxbh=MAX(bh,lmaxbh)		; bottom = most prob. skirt value
            ENDIF ELSE BEGIN
; otherwise take lowest level in dark frame
                lmini=MIN(WHERE(cih GE .001))
                lmaxbh=lmini	; rename this parameter for compatibility below
            ENDELSE
        ENDIF ELSE BEGIN
; if no file name, generate an identifiably fake image
            im=BYTARR(256,256)
            im(RANDOMU(s,100))=1
            itime=BYTARR(4)
            expt=-1
            lmaxi=1
            lmaxbh=0
        ENDELSE

        ip=ip+1		; index number for screen position of image
; show image orientation, one time only
        IF ip EQ 0 THEN BEGIN
            nflag=isnumber(KIh.misc.az,azimuth)
            IF nflag NE 0 THEN BEGIN 	; if the camera azimuth is known,
                azrad=!DTOR*azimuth		; plot arrows showing the azimuth
                leftend=scl(res=120.,wwidth/2.)-twidthd/2.
                maxlen=(leftend-scl(res=120.,1.5*X_charwidth)) $
                       /scl(res=120.,(COS(azrad)+SIN(azrad))*X_charheight)
                arlen=MIN([FIX(maxlen),4])*X_charheight
                tailx=REPLICATE(scl(res=120.,0.5*X_charwidth+SIN(azrad)*arlen),2)
                taily=REPLICATE(scl(res=120.,nincol*256),2)
                headx=[tailx(0)-scl(res=120.,arlen*SIN(azrad)), $
                       tailx(1)+scl(res=120.,arlen*COS(azrad))]
                heady=[scl(res=120.,nincol*256+arlen*COS(azrad)), $
                       scl(res=120.,nincol*256+arlen*SIN(azrad))]
                arrow,tailx,taily,headx,heady,/device,hsize=-0.2
                XYOUTS,tailx(0),heady(0)-scl(res=120.,X_charheight),'N' $
                       ,ALIGNMENT=0,/DEVICE,CHARSIZE=csz
                XYOUTS,headx(1)+scl(res=120.,X_charwidth),heady(1),'E' $
                       ,ALIGNMENT=0.5,/DEVICE,CHARSIZE=csz
            ENDIF ELSE BEGIN		; otherwise show cardinal directions
                XYOUTS,scl(3*X_charwidth,res=120.) $
                      ,scl(nincol*256+3*X_charheight,res=120.) $
                      ,'N',/DEVICE,ALIGNMENT=0.5,CHARSIZE=csz
                XYOUTS,scl(res=120.,X_charwidth) $
                      ,scl(res=120.,nincol*256+1.5*X_charheight),'W' $
                      ,/DEVICE,ALIGNMENT=0.5,CHARSIZE=csz
                XYOUTS,scl(res=120.,5*X_charwidth) $
                      ,scl(res=120.,nincol*256+1.5*X_charheight),'E' $
                      ,/DEVICE,ALIGNMENT=0.5,CHARSIZE=csz
                XYOUTS,scl(res=120.,3*X_charwidth) $
                      ,scl(res=120.,nincol*256),'S',/DEVICE $
                      ,ALIGNMENT=0.5,CHARSIZE=csz
            ENDELSE
        ENDIF
; display the image
        mtv,ROTATE(im*mask,5),low=lmaxbh,high=lmaxi,res=120. $
           ,(ip MOD ninrow)*256,(nincol-1-(ip MOD ninwin)/ninrow)*256
; display image camera, filter, and sequence numbers
        IF fnlen GT 0 $
            THEN camfilt=STRMID(files(i,j),fnlen-9,2) $
            ELSE camfilt=STRMID(files(i,1-j),fnlen-9,2)
        XYOUTS,scl(res=120.,(ip MOD ninrow)*256+1) $
              ,scl(res=120.,(nincol-(ip MOD ninwin)/ninrow)*256-1-X_charheight) $
              ,camfilt+':'+ext_list(i,j),/DEVICE,ALIGNMENT=0.0,CHARSIZE=csz $
              ,COLOR=!D.N_COLORS
; and image time
        XYOUTS,scl(res=120.,(ip MOD ninrow)*256+1) $
              ,scl(res=120.,(nincol-1-(ip MOD ninwin)/ninrow)*256+3+X_charheight) $
              ,/DEVICE,ALIGNMENT=0.0,CHARSIZE=csz,COLOR=!D.N_COLORS $
              ,STRING(itime(0:2),FORMAT='(I2.2,":",I2.2,":",I2.2)')
; and exposure time
        XYOUTS,scl(res=120.,(ip MOD ninrow)*256+1) $
              ,scl(res=120.,(nincol-1-(ip MOD ninwin)/ninrow)*256+1) $
              ,STRING(expt,FORMAT='(F5.1,1X,"s")'),/DEVICE,ALIGNMENT=0.0 $
              ,CHARSIZE=csz,COLOR=!D.N_COLORS
; and maximum displayed DN
        XYOUTS,scl(res=120.,((ip MOD ninrow)+1)*256-1) $
              ,scl(res=120.,(nincol-(ip MOD ninwin)/ninrow)*256-1-X_charheight) $
              ,STRING(FIX(lmaxi),FORMAT='("< ",I0.0," DN")') $
              ,/DEVICE,ALIGNMENT=1.0,CHARSIZE=csz,COLOR=!D.N_COLORS
; and minimum displayed DN
        XYOUTS,scl(res=120.,((ip MOD ninrow)+1)*256-1) $
              ,scl(res=120.,(nincol-1-(ip MOD ninwin)/ninrow)*256+1) $
              ,STRING(FIX(lmaxbh),FORMAT='("> ",I0.0," DN")'),/DEVICE $
              ,ALIGNMENT=1.0,CHARSIZE=csz,COLOR=!D.N_COLORS
        PRINT,'.',FORMAT='(A,$)'
        IF (ip+1) MOD ninwin EQ 0 THEN PRINT,'|',FORMAT='(A,$)'
    ENDFOR
    IF GET_KBRD(1) EQ 'q' THEN GOTO,DONE
ENDFOR
DONE:
PRINT,''
END