PRO grid,year,month,day,cam

PRINT,'Loading image area'
sname='/jasper/cnsr3_data1/sig'+STRTRIM(cam,2)+'.EU'
OPENR,sunit,sname,/GET_LUN
ns=0L
READU,sunit,ns
sig=LONARR(ns)
READU,sunit,sig
CLOSE,sunit
FREE_LUN,sunit

PRINT,'Finding image boundary'
mk,mask
mask(sig)=1
bmask=ROTATE(REBIN(mask,512,512,/SAMPLE),2)
bedge=(SOBEL(bmask) GT 0)
bedge(0,134:419)=1
bedge(1,134:419)=1
bedge(142:351,0)=1
bedge(142:351,1)=1

PRINT,'Setting up plot coordinates'
IF !D.NAME EQ 'PS' THEN BEGIN
    psres=512/(!D.Y_SIZE/2540.)
    len=scl(512,res=psres)
    ENDIF ELSE BEGIN
    psres=90.7
    len=512
    ERASE
    ENDELSE
xoff=(!D.X_SIZE-len)/2.
yoff=0
xp0=xoff
xp1=xp0+len
yp0=0
yp1=len

PRINT,'Plotting image boundary'
mtv,BYTSCL(bedge,TOP=!D.N_COLORS-1),xoff,yoff,res=psres

PRINT,'Getting azimuth and zenith angle coordinates'
col=BINDGEN(256)#REPLICATE(1B,1,256)
row=TRANSPOSE(col)
cr2az,year,month,day,cam,0,col,row,az,za,site=3,/degrees

PRINT,'Setting up plot'
PLOT,INDGEN(512),INDGEN(512),/DEVICE,/NODATA,/NOERASE $
    ,POSITION=[xp0,yp0,xp1,yp1],XMARGIN=[0,0],YMARGIN=[0,0] $
    ,XRANGE=[0,511],YRANGE=[0,511],XSTYLE=5,YSTYLE=5

PRINT,'Plotting zenith angle contours'
CONTOUR,ROTATE(REBIN(za*mask,512,512),2),LEVELS=FINDGEN(8)*10 $
       ,C_LABELS=REPLICATE(1B,8),/OVERPLOT,/DEVICE $
       ,POSITION=[xp0,yp0,xp1,yp1],XRANGE=[0,511],YRANGE=[0,511] $
       ,XSTYLE=5,YSTYLE=5,COLOR=!D.N_COLORS-1

PRINT,'Plotting azimuth contours'
CONTOUR,ROTATE(REBIN(az*mask,512,512),2),LEVELS=FINDGEN(12)*30 $
       ,C_LABELS=REPLICATE(1B,12),/OVERPLOT,/DEVICE $
       ,POSITION=[xp0,yp0,xp1,yp1],XRANGE=[0,511],YRANGE=[0,511] $
       ,XSTYLE=5,YSTYLE=5,COLOR=!D.N_COLORS-1

PRINT,'Labelling image'
xcs=!D.X_CH_SIZE
ycs=!D.Y_CH_SIZE
STOP
XYOUTS,xp0-3*xcs,yp1-ycs,'SEQUENCE #',ALIGNMENT=0,/DEVICE
XYOUTS,xp0-3*xcs,yp0,'EXPOSURE TIME',ALIGNMENT=0,/DEVICE
XYOUTS,xp0-3*xcs,yp0+1.25*ycs,'EXPOSURE START',ALIGNMENT=0,/DEVICE
XYOUTS,xp1+3*xcs,yp0,'MIN. BRIGHTNESS',ALIGNMENT=1,/DEVICE
XYOUTS,xp1+3*xcs,yp1-ycs,'MAX. BRIGHTNESS',ALIGNMENT=1,/DEVICE

RETURN
END