; QL.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. ; Modified 1993/6/11 from SEQ.PRO to display large numbers of images ; in compressed form. ; Modified 1993/11/... to show image sequences from Eureka, and to apply ; flat-field corrections and conversion to intensities. ; ; create byte-valued input parameters cam=0B & filt=cam & day=cam & month=cam ; load various string-valued arrays needed initialize,filters,months,days,modes,mpaths,firstchars s='' !ORDER=1 ; uncover screen where text is to be printed IF !D.NAME EQ 'X' THEN WDELETE,0 ; prompt user and get input fac=1. resp='' PRINT,'Do you want to display images as well as save them? ([Y]/N)' READ,resp dispim=(STRUPCASE(resp) NE 'N') IF dispim THEN BEGIN ; set colour table loadct,13,/silent PRINT,'Enter the scale compression factor to be applied (over normal scale; default 1)' READ,resp IF resp NE '' THEN fac=FIX(resp) ENDIF PRINT,'Enter the starting year, month, day, UT hour, minute, and second desired' READ,year,month,day,hr,min,sec IF dispim THEN BEGIN PRINT,'View all images (0; default) or only 1 cycle per hour (1)?' READ,resp short=(STRUPCASE(resp) EQ '1') ENDIF ELSE short=0B PRINT,'Convert DN to intensity? ([Y]/N)' READ,resp rayleighs=(STRUPCASE(resp) NE 'N') IF rayleighs THEN BEGIN PRINT,'Present backgrounds as differential brightnesses (R/nm)? ([Y]/N)' READ,resp RperA=(STRUPCASE(resp) NE 'N') ; For now, disable background subtraction (DPS, 940105) ; PRINT,'Do background subtractions? ([Y]/N)' ; READ,resp resp='N' bgsubtract=(STRUPCASE(resp) NE 'N') IF bgsubtract THEN BEGIN PRINT,'Settle for approximate background subtraction? ([Y]/N)' READ,resp approx=(STRUPCASE(resp) NE 'N') ENDIF ENDIF ELSE BEGIN RperA=0 bgsubtract=0 ENDELSE IF dispim THEN BEGIN PRINT,'Fix the maximum value to be displayed for each camera? (Y/[N])' READ,resp resp=STRUPCASE(resp) ceiling_set=(resp EQ 'Y') IF ceiling_set THEN BEGIN topnum=FLTARR(2) FOR j=0,1 DO BEGIN PRINT,'Enter the ceiling value for camera '+STRTRIM(j,2) READ,top topnum(j)=top ENDFOR ENDIF ENDIF ELSE ceiling_set=0 PRINT,'Display debugging information? (Y/[N])' READ,resp debug=(STRUPCASE(resp) EQ 'Y') IF NOT short THEN BEGIN PRINT,'Save images? ([Y]/N) (This option also creates quick-look data files)' READ,resp imsave=(STRUPCASE(resp) NE 'N') ENDIF ELSE imsave=0 cfsave=BYTARR(10) savunit=BYTARR(10) IF imsave THEN BEGIN PRINT,'Enter the wavelengths to be saved (nm)' savwav='' READ,savwav nwsave=nwrds(savwav) IF nwsave GT 0 THEN BEGIN wsave=FLTARR(nwsave) FOR jsave=0,nwsave-1 DO wsave(jsave)=FLOAT(getwrd(savwav,jsave)) ENDIF ELSE BEGIN MESSAGE,'No wavelengths entered',/INFORMATIONAL ENDELSE init_qlfiles,qlunits ENDIF ELSE nwsave=0 ; ensure a 4-digit year for display; make reasonable assumptions IF year LT 1900 THEN BEGIN yr=year IF year LT 50 THEN year=year+2000 IF year GT 50 THEN year=year+1900 ENDIF ; create title for image plot IF dispim THEN BEGIN wtitle=STRING(year,months(month-1),day $ ,FORMAT='("Polar Camera",2X,I4,1X,A3,1X,I0.2)') ; define mask for images, with fixed horizon at 150 superpixels mask0=BYTARR(256,256) prad=150 d256=SHIFT(dist(256),128,128) & base=WHERE(d256 GT 150) sig0=WHERE(d256 LE prad) & mask0(sig0)=1 sig1=sig0 & mask1=mask0 ENDIF ; make lists of available files of desired type rdilist,2,year,month,day,hr,min,sec,n,t,files ; Reverse the indices of the time and file name arrays as returned ; by rdilist t=TRANSPOSE(t) files=TRANSPOSE(files) ext_list=extension(files) ; Identify and delete calibration files from time and file name arrays, ; so no impossibly long cycles appear match=STRING(day,FORMAT='("/c",I2.2)') calfiles=WHERE(STRPOS(files,match) NE -1) one2two,calfiles,files,numrow,camcol keep=setminus(INDGEN(n),numrow) n=N_ELEMENTS(keep) t=t(keep,*) files=files(keep,*) PRINT,n,FORMAT='(I0," image pairs excluding calibration images")' ; Identify the first images or image pairs of repeated sequences. ; Typical exposure times are 1, 10, and 60 s, with correponding ; inter-image elapsed times of 5, 14, and 64 s; any consistently- ; occurring elapsed times much longer than these probably indicate ; the end of a data-taking sequence and the start of another one. dt=(t(*,0)-SHIFT(t(*,0),1))>0 cycle=[WHERE(dt GT 75),n] ; extra element simplifies finding ; last image in last cycle cyclen=(cycle-SHIFT(cycle,1))>0 cycle=cycle(0:N_ELEMENTS(cycle)-2) cyclen=cyclen(1:*) ; If a cycle is too long, investigate. IF MAX(cyclen) GT 5 THEN BEGIN PRINT,'Maximum cycle length is '+STRTRIM(MAX(cyclen),2)+'; check it out!' STOP ENDIF ; If only hourly cycles are to be shown, reduce cycle and cyclen ; appropriately, to display only dark cycle and the following ; open-shutter cycle each hour. IF dispim THEN BEGIN IF short THEN BEGIN wdkcyc=WHERE(cyclen EQ 2) wnextcyc=wdkcyc+1 dkcyc=cycle(wdkcyc) nextcyc=setminus(cycle(wnextcyc),dkcyc) scycle=[dkcyc,nextcyc] scyclen=[cyclen(wdkcyc),cyclen(wnextcyc)] scyclen=scyclen(0:N_ELEMENTS(scycle)-1) ssc=SORT(scycle) cycle=scycle(ssc) cyclen=scyclen(ssc) ENDIF ENDIF maxcyclen=MAX(cyclen) cycleparms=STRARR(5,2*maxcyclen); to store image annotation strings IF dispim AND (maxcyclen GT 5*fac) THEN BEGIN ; avoid display problems longcyc=STRTRIM(maxcyclen,2) ; with long cycles maxninrow=STRTRIM(FIX(5*fac),2) minnewfac=STRTRIM(CEIL(maxcyclen/5.),2) PRINT,'The longest cycle is '+longcyc+' images long; no more than' PRINT,maxninrow+' images will fit side by side on the screen.' PRINT,'You must change the image scale compression factor.' PRINT,'Enter a new factor by which to shrink the images (min. '+minnewfac+')' READ,fac ENDIF ; create and label image display window imgsize=256/fac IF dispim THEN BEGIN X_charheight=10 ; try to keep scalings that apply in X windows X_charwidth=8 barwidth=16 ninrow=maxcyclen wwidth=ninrow*imgsize dheight=2*imgsize+barwidth ; run colour bar along bottom display=BYTARR(wwidth,dheight) ; array to hold one cycle of images wheight=dheight+4*X_charheight ENDIF ; Load the pseudo-flat-field images generated at Eureka, for later ; application to the images before display. ffimg=BYTARR(5*256,2*256) ffim=BYTARR(256,256) FOR j=0,1 DO BEGIN FOR i=0,4 DO BEGIN ffname=STRING(j,i,FORMAT='("~cnsr3/flat/ffimg_",2I1,".dat")') OPENR,ffunit,ffname,/GET_LUN READU,ffunit,ffim CLOSE,ffunit FREE_LUN,ffunit ffimg(i*256,j*256)=ffim ENDFOR ENDFOR ; Load the Rayleigh-second/DN values for converting to intensities IF rayleighs THEN BEGIN getsens,year,month,day,sens0,sens1 sens=[sens0,sens1] ; Load the center wavelengths and rectangular passbands for background subtraction unitstring='R' cwl=FLTARR(10) rpb=cwl OPENR,funit,'/jasper/cnsr3_data1/passband/cwl.rpb',/GET_LUN READF,funit,cwl,rpb CLOSE,funit FREE_LUN,funit ENDIF ELSE unitstring='DN' ; If images are to be saved, open the associated files to hold them IF nwsave GT 0 THEN BEGIN seqnum=INTARR(10) FOR jsave=0,nwsave-1 DO BEGIN tmp=WHERE(ABS(wsave(jsave)-cwl) LT rpb/2.) tmp=tmp(0) cfsave(tmp)=1 fname=STRING(yr,month,day,tmp/5,tmp MOD 5 $ ,FORMAT='("/jasper/cnsr3_data1/save/",3I2.2,2I1,".")') CASE 1 OF (NOT rayleighs): fname=fname+'cDN' (rayleighs AND (NOT bgsubtract)): fname=fname+'ray' bgsubtract: fname=fname+'bgs' ENDCASE OPENW,tt,fname,/GET_LUN savunit(tmp)=tt CASE tmp OF 0: save00=ASSOC(savunit(0),BYTARR(131584)) 1: save01=ASSOC(savunit(1),BYTARR(131584)) 2: save02=ASSOC(savunit(2),BYTARR(131584)) 3: save03=ASSOC(savunit(3),BYTARR(131584)) 4: save04=ASSOC(savunit(4),BYTARR(131584)) 5: save10=ASSOC(savunit(5),BYTARR(131584)) 6: save11=ASSOC(savunit(6),BYTARR(131584)) 7: save12=ASSOC(savunit(7),BYTARR(131584)) 8: save13=ASSOC(savunit(8),BYTARR(131584)) 9: save14=ASSOC(savunit(9),BYTARR(131584)) ENDCASE ENDFOR ENDIF ; Create the format string for describing background subtractions IF debug THEN $ bgform='("Image ",I1,", camera ",I1,";"/" subtract image ",I1,", camera ",I1)' IF dispim THEN BEGIN IF !D.NAME EQ 'X' 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/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 SPAWN,'whoami',username 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 SPAWN,'date -u',datetime XYOUTS,scl(res=120.,wwidth-X_charwidth),scl(res=120.,wheight-2.5*X_charheight) $ ,STRMID(datetime,11,9)+' UT',ALIGNMENT=1.0,/DEVICE,CHARSIZE=csz XYOUTS,scl(res=120.,wwidth-X_charwidth),scl(res=120.,wheight-1.1*X_charheight) $ ,STRMID(datetime,4,7)+STRMID(datetime,24,4),ALIGNMENT=1.0,/DEVICE,CHARSIZE=csz ; create and display color bar bar=REBIN(BINDGEN(256,barwidth),ninrow*imgsize,barwidth) display(0,dheight-barwidth)=bar ; mtv,bar,0,0,res=120. ; show image orientation, one time only XYOUTS,scl(3*X_charwidth,res=120.) $ ,scl(barwidth+2*imgsize+3*X_charheight,res=120.) $ ,'N',/DEVICE,ALIGNMENT=0.5,CHARSIZE=csz XYOUTS,scl(res=120.,X_charwidth) $ ,scl(res=120.,barwidth+2*imgsize+1.5*X_charheight),'W' $ ,/DEVICE,ALIGNMENT=0.5,CHARSIZE=csz XYOUTS,scl(res=120.,5*X_charwidth) $ ,scl(res=120.,barwidth+2*imgsize+1.5*X_charheight),'E' $ ,/DEVICE,ALIGNMENT=0.5,CHARSIZE=csz XYOUTS,scl(res=120.,3*X_charwidth) $ ,scl(res=120.,barwidth+2*imgsize),'S',/DEVICE $ ,ALIGNMENT=0.5,CHARSIZE=csz ENDIF ; Calculate the azimuth and zenith angle images for each camera col=BINDGEN(256)#REPLICATE(1B,1,256) row=TRANSPOSE(col) ; Decide where Polar Camera was located CASE 1 OF (yr EQ 93) AND ((month EQ 1) OR (month EQ 2)): PCsite=2 ; RL ((yr EQ 93) AND (month GE 10)) OR (yr GE 94): PCsite=3 ; Eureka ELSE: PCsite=1 ; U of C ENDCASE cr2az,0,0,col,row,az0,za0,site=PCsite cr2az,1,0,col,row,az1,za1,site=PCsite ; Identify image pixels within 2 degrees of the zenith zen0=WHERE(za0 LT 2.) zen1=WHERE(za1 LT 2.) IF bgsubtract THEN BEGIN IF dispim THEN BEGIN WINDOW,3,XSIZE=256,YSIZE=512,XPOS=0,YPOS=0 WSET,0 ENDIF ; Restrict future attention to image points where data are ; significant - use flat field correction images as guides ; sig0.EU is LONARR(57649) ; sig1.EU is LONARR(58383) IF debug THEN PRINT,'Identifying image interior...' OPENR,s0unit,'/jasper/cnsr3_data1/sig0.EU',/GET_LUN s0=0L READU,s0unit,s0 sig0=LONARR(s0) READU,s0unit,sig0 CLOSE,s0unit FREE_LUN,s0unit OPENR,s1unit,'/jasper/cnsr3_data1/sig1.EU',/GET_LUN s1=0L READU,s1unit,s1 sig1=LONARR(s1) READU,s1unit,sig1 CLOSE,s1unit FREE_LUN,s1unit ; Triangulate over significant parts of images ; triangles0.EU is LONARR(3,115213) ; triangles1.EU is LONARR(3,116598) IF debug THEN PRINT,'Triangulating az-za maps...' OPENR,t0unit,'/jasper/cnsr3_data1/triangles0.EU',/GET_LUN triangles0=ASSOC(t0unit,LONARR(3,115213)) OPENR,t1unit,'/jasper/cnsr3_data1/triangles1.EU',/GET_LUN triangles1=ASSOC(t1unit,LONARR(3,116598)) gridspacing=[1.0,1.0] ; 1 deg spacing in az and za boundaries=[0.0,0.0,360.0,72.0] ; minaz,minza,maxaz,maxza ; Generate coefficients for background subtraction coeffs=FLTARR(10,3) ; 10 filters, 3 coeffs per filter IF approx THEN BEGIN coeffs(*,0)=[0.,0.,1.,0.,1.,0.,0.,0.,1.,0.] coeffs(*,1)=[0.,0.,0.,1.,0.,0.,0.,0.,0.,0.] coeffs(*,2)=[1.,0.,0.,0.,0.,1.,0.,1.,0.,0.] same=[0,0,1,1,1,1,0,1,0,0] ENDIF ELSE BEGIN ; Compute Lagrange polynomial coefficients for extrapolating ; backgrounds from measured wavelengths IF debug THEN PRINT,'Generating Lagrange interpolation coefficients...' denom0=(cwl(1)-cwl(9))*(cwl(1)-cwl(6)) denom1=(cwl(9)-cwl(1))*(cwl(9)-cwl(6)) denom2=(cwl(6)-cwl(1))*(cwl(6)-cwl(9)) FOR lw=0,9 DO BEGIN IF (lw NE 1) AND (lw NE 6) AND (lw NE 9) THEN BEGIN coeffs(lw,0)=((cwl(lw)-cwl(9))*(cwl(lw)-cwl(6)))/denom0 coeffs(lw,1)=((cwl(lw)-cwl(1))*(cwl(lw)-cwl(6)))/denom1 coeffs(lw,2)=((cwl(lw)-cwl(1))*(cwl(lw)-cwl(9)))/denom2 ENDIF ENDFOR ENDELSE pp=FINDGEN(100)*0.01+0.005 ENDIF ; read in each cycle of images, decode the headers, rebin to ; (imgsize)x(imgsize) if necessary, scale to 8 bits and display blank=BYTARR(imgsize,imgsize) itime=BYTARR(4) ; load the dark array with representative dark frames for now, ; and supersede them with actual measured ones as they are available darks=INTARR(512,512) ; dummy variables for now rdmeandk,0,10,dk,/nointeractive darks(0,0)=dk rdmeandk,0,60,dk,/nointeractive darks(256,0)=dk rdmeandk,1,10,dk,/nointeractive darks(0,256)=dk rdmeandk,1,60,dk,/nointeractive darks(256,256)=dk dkexp=[10,60] FOR k=0,N_ELEMENTS(cycle)-1 DO BEGIN ; create (or refresh) arrays for temporary storage keep=FLTARR(cyclen(k)*256,512) ; processed images before background subtraction keepexst=BYTARR(cyclen(k),2); does image exist? keeptype=keepexst ; image type code, 1-6 keepfilt=keepexst ; filters images taken through keeptime=STRARR(cyclen(k),2); times of images keepexp=keepexst ; exposure lengths keephb=BYTARR(512,2*cyclen(k)) IF NOT dispim THEN PRINT,k,cycle(k),cyclen(k) $ ,FORMAT='("Cycle ",I3," starts at pair "' $ +',I3,", and includes ",I1," pairs")' ip=-1 FOR i=cycle(k),cycle(k)+cyclen(k)-1 DO BEGIN FOR j=0,1 DO BEGIN ip=ip+1 ; update index number for screen position of image p=ip/2 ; position in row ; Is the file path string a valid file name? image_exists=(STRPOS(files(i,j),'___') EQ -1) IF image_exists THEN BEGIN rdkimg,files(i,j),h,im IF debug THEN PRINT,files(i,j)+' read in' KIh=gethd(h) FOR l=0,3 DO itime(l)=KIh.misc.tm.(l) exp=KIh.misc.exp_scale*KIh.exp.exposure/1000. exptype=KIh.misc.exposure_type ; If this image is a dark frame, prepare to save it and all others ; in this cycle, otherwise do the appropriate dark subtraction IF exptype EQ 5 THEN BEGIN ; dark frames - save 'em! first_time=(NOT KEYWORD_SET(width)) width=512/KIh.exp.sbin ; set up the dark frames length=512/KIh.exp.pbin ; and exp. times arrays IF NOT first_time THEN BEGIN n_darks=N_ELEMENTS(WHERE(dkexp)) old_darks=darks(0:n_darks*width-1,*) old_dkexp=dkexp(0:n_darks-1) ENDIF darks=INTARR(cyclen(k)*width,2*length) dkexp=FLTARR(cyclen(k)) IF NOT first_time THEN BEGIN darks(0,0)=old_darks dkexp(0)=old_dkexp ENDIF IF debug THEN PRINT,'Dark frames being loaded!' sim=MEDIAN(im,5) ; lose hot spots darks(p*width,j*length)=sim ; store the dark frame IF j EQ 0 THEN dkexp(p)=exp ; and the exposure time ENDIF ELSE BEGIN ; otherwise subtract 'em pd=WHERE(dkexp EQ exp) ; find right exposure time IF pd(0) NE -1 THEN BEGIN pd=pd(0) IF KEYWORD_SET(width) THEN BEGIN dwidth=width dlength=length ENDIF ELSE BEGIN dwidth=256 dlength=256 ENDELSE dkim=darks(pd*dwidth:(pd+1)*dwidth-1,j*dlength:(j+1)*dlength-1) im=(im-dkim)>0 ENDIF ENDELSE ; if image isn't 256x256, rebin it IF (KIh.exp.sbin NE 2) OR (KIh.exp.pbin NE 2) THEN BEGIN 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 ENDIF ; Apply pseudo-flat-field correction to dark-subtracted image IF (exptype NE 3) AND (exptype NE 5) THEN BEGIN filt=KIh.misc.fw im=im*(ffimg(filt*256:(filt+1)*256-1,j*256:(j+1)*256-1)/100.) ; Convert to Rayleighs if desired IF rayleighs THEN BEGIN calfac=sens(5*j+filt)/(exp+0.045) im=im*calfac IF debug THEN PRINT,calfac $ ,FORMAT='("*",F5.2," to get Rayleighs")' IF RperA THEN BEGIN IF (j EQ 0) AND (filt EQ 1) THEN im=im/rpb(1) IF (j EQ 1) AND (filt EQ 1) THEN im=im/rpb(6) IF (j EQ 1) AND (filt EQ 4) THEN im=im/rpb(9) ENDIF ENDIF ENDIF ; Store processed images for subsequent background subtraction or other use keepexst(p,j)=1B keep(p*256:(p+1)*256-1,j*256:(j+1)*256-1)=im keeptype(p,j)=BYTE(exptype) keepfilt(p,j)=filt keeptime(p,j)=STRING(itime(0:2),FORMAT='(3(I2.2,:,":"))') keepexp(p,j)=exp keephb(0,ip)=h ; Store zenith intensity in appropriate quick-look data file IF ((exptype NE 3) AND (exptype NE 5)) AND imsave THEN BEGIN IF j EQ 0 THEN zdata=im(zen0) ELSE zdata=im(zen1) save_qldata,j,filt,qlunits,itime,exp,zdata ENDIF ; Store background images separately for later use IF bgsubtract THEN BEGIN IF (j EQ 0) AND (filt EQ 1) THEN BEGIN bg608=im IF NOT RperA THEN bg608=bg608/rpb(1) ENDIF IF (j EQ 1) AND (filt EQ 1) THEN BEGIN bg820=im IF NOT RperA THEN bg820=bg820/rpb(6) ENDIF IF (j EQ 1) AND (filt EQ 4) THEN BEGIN bg714=im IF NOT RperA THEN bg714=bg714/rpb(9) ENDIF ENDIF ENDIF ENDFOR ENDFOR ; For now we only do approximate background subtraction darkcyc=(TOTAL(keeptype EQ 5) EQ cyclen(k)*2) IF bgsubtract AND (NOT darkcyc) THEN BEGIN ; Create uniform grids of background intensities vs az and za gbg608=TRIGRID(az0(sig0),za0(sig0),bg608(sig0),triangles0(0) $ ,gridspacing,boundaries) gbg714=TRIGRID(az1(sig1),za1(sig1),bg714(sig1),triangles1(0) $ ,gridspacing,boundaries) gbg820=TRIGRID(az1(sig1),za1(sig1),bg820(sig1),triangles1(0) $ ,gridspacing,boundaries) bgmade=BYTARR(10) OPENU,bgunit,'/jasper/cnsr3_data1/bgsave.dat',/GET_LUN bgs=ASSOC(bgunit,INTARR(256,256)) FOR i=0,cyclen(k)-1 DO BEGIN FOR j=0,1 DO BEGIN IF keepexst(i,j) THEN BEGIN fptr=5*j+keepfilt(i,j) IF NOT ((fptr EQ 1) OR (fptr EQ 6) OR (fptr EQ 9)) THEN BEGIN IF bgmade(fptr) THEN bg=bgs(fptr) ELSE BEGIN ; attempt at general method IF j EQ 0 THEN BEGIN sig=sig0 az=az0(sig0) za=za0(sig0) ENDIF ELSE BEGIN sig=sig1 az=az1(sig1) za=za1(sig1) ENDELSE im=keep(i*256:(i+1)*256-1,j*256:(j+1)*256-1) ; interpolate background from recorded images IF approx AND same(fptr) THEN BEGIN bg=coeffs(fptr,0)*bg608 $ +coeffs(fptr,1)*bg714 $ +coeffs(fptr,2)*bg820 bg=bg*rpb(fptr) ENDIF ELSE BEGIN bga=INTERPOLATE(coeffs(fptr,0)*gbg608 $ +coeffs(fptr,1)*gbg714 $ +coeffs(fptr,2)*gbg820 $ ,az,za)*rpb(fptr) az2cr,j,keepfilt(i,j),az,za,cl,rw,/small,site=PCsite ; store into a 256x256 array bg=FLTARR(256,256) bg(cl,rw)=bga ; find out least value a real data point can have dp=percentile(bg(sig),pp) thr=MIN(dp(WHERE(dp GT 0))) IF debug THEN PRINT,j,keepfilt(i,j),thr ; truncate all background values below that value bg(WHERE(bg LE thr))=0 ; fill holes in background image FOR bgr=0,255 DO bg(*,bgr)=fixgaps(bg(*,bgr)) bg=bg>0 ENDELSE IF dispim THEN BEGIN WSET,3 imax=percentile(im,0.99) TV,BYTSCL(im,MIN=0,MAX=imax,TOP=!D.N_COLORS-1),0,0 TV,BYTSCL(bg,MIN=0,MAX=imax,TOP=!D.N_COLORS-1),0,256 WSET,0 ENDIF bgs(fptr)=bg bgmade(fptr)=1 ENDELSE STOP,'check out im,bg; .con to continue' keep(i*256:(i+1)*256-1,0:255)=(im-bg)>0 ENDIF ENDIF ; ; OH 829 nm channel ; IF (j EQ 0) AND (keepfilt(i,j) EQ 0) THEN BEGIN ; im=keep(i*256:(i+1)*256-1,0:255) ; bgi=WHERE(keepfilt(*,1) EQ 1) ; bgi=bgi(0) ; bg=keep(bgi*256:(bgi+1)*256-1,256:511) ; IF debug THEN PRINT,i,j,bgi,1,FORMAT=bgform ; IF NOT RperA THEN bg=bg/rpb(6) ; bg=bg*rpb(0) ; fptr=5*j+keepfilt(i,j) ; STOP ; bg=INTERPOLATE(coeffs(fptr,0)*gbg608 $ ; +coeffs(fptr,1)*gbg714 $ ; +coeffs(fptr,2)*gbg820 $ ; ,az0(sig0),za0(sig0))*rpb(fptr) ; WSET,3 ; imax=MAX(im) ; TV,BYTSCL(im,MIN=0,MAX=imax,TOP=!D.N_COLORS-1),0,0 ; TV,BYTSCL(bg,MIN=0,MAX=imax,TOP=!D.N_COLORS-1),256,0 ; WSET,0 ; keep(i*256:(i+1)*256-1,0:255)=(im-bg)>0 ; ENDIF ENDFOR ENDFOR CLOSE,bgunit FREE_LUN,bgunit ENDIF ip=-1 FOR i=cycle(k),cycle(k)+cyclen(k)-1 DO BEGIN ic=i-cycle(k) FOR j=0,1 DO BEGIN ip=ip+1 IF keepexst(ic,j) THEN BEGIN ; identify useful min and max values for byte-scaling im=keep(ic*256:(ic+1)*256-1,j*256:(j+1)*256-1) IF dispim THEN BEGIN IF NOT ceiling_set THEN BEGIN IF MAX(im) GT 1 THEN BEGIN IF j EQ 0 $ THEN lmaxi=percentile(im(sig0),0.999) $ ELSE lmaxi=percentile(im(sig1),0.999) ENDIF ELSE lmaxi=1 ; hbin=(10^CEIL(ALOG10((MAX(im)/32767)>1)))>1 ; ih=HISTOGRAM(im(wprad),MIN=0,BIN=hbin); 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=hbin*MAX(WHERE(cih LE .999)); top 0.1% of pixels saturate ENDIF ELSE lmaxi=topnum(j) ; if viewing sky image, take dark 'skirt' as min value exptype=keeptype(ic,j) shut=(exptype EQ 3) OR (exptype EQ 5) ; bias or dark frame IF NOT shut THEN BEGIN IF lmaxi NE 1 THEN BEGIN lmaxbh=percentile(im(base),0.5)>0 ; bh=HISTOGRAM(im(base),MIN=0,BIN=hbin) ; maxbh=MAX(bh,lmaxbh) ; bottom = most prob. skirt value ; lmaxbh=hbin*lmaxbh ENDIF ELSE lmaxbh=0 ENDIF ELSE BEGIN ; otherwise take lowest level in dark frame ; lmini=MIN(WHERE(cih GE .001)) IF j EQ 0 $ THEN lmaxbh=percentile(im(sig0),0.001) $ ELSE lmaxbh=percentile(im(sig1),0.001) ENDELSE ; rebin the image to the scale requested by the user IF debug THEN PRINT,lmaxbh,lmaxi,unitstring $ ,FORMAT='("Useful values span ",I5,"-",I5,1X,A)' IF j EQ 0 $ THEN im=REBIN(im*mask0,imgsize,imgsize,/SAMPLE) $ ELSE im=REBIN(im*mask1,imgsize,imgsize,/SAMPLE) display(ic*imgsize,j*imgsize)= $ ROTATE(BYTSCL(im,MIN=lmaxbh,MAX=lmaxi),5) ENDIF ; save desired images in associated files fptr=5*j+keepfilt(ic,j) exptype=keeptype(ic,j) IF ((exptype NE 3) AND (exptype NE 5)) AND cfsave(fptr) THEN BEGIN hb=keephb(*,ip) bbin=2^CEIL(ALOG(MAX(im)/32767.)/ALOG(2.))>1 hb(290)=bbin outrec=[hb,BYTE(FIX(ROUND(im/FLOAT(bbin))),0,131072)] seq=seqnum(fptr) CASE fptr OF 0: save00(seq)=outrec 1: save01(seq)=outrec 2: save02(seq)=outrec 3: save03(seq)=outrec 4: save04(seq)=outrec 5: save10(seq)=outrec 6: save11(seq)=outrec 7: save12(seq)=outrec 8: save13(seq)=outrec 9: save14(seq)=outrec ENDCASE seqnum(fptr)=seqnum(fptr)+1 ENDIF ; create and save image annotation strings IF dispim THEN BEGIN camfilt=STRMID(files(i,j),STRLEN(files(i,j))-9,2) cycleparms(0,ip)=camfilt+':'+ext_list(i,j) cycleparms(1,ip)=keeptime(ic,j) cycleparms(2,ip)=STRING(keepexp(ic,j),FORMAT='(F5.1,1X,"s")') cycleparms(3,ip)=STRING(ROUND(lmaxi),unitstring $ ,FORMAT='("< ",I0.0,A)') cycleparms(4,ip)=STRING(ROUND(lmaxbh),unitstring $ ,FORMAT='("> ",I0.0,A)') ENDIF ENDIF ELSE BEGIN ; No image found; load blank image and use dummy annotation strings IF dispim THEN BEGIN display(ic*imgsize,j*imgsize)=blank cycleparms(0,ip)='??:-1' cycleparms(1,ip)='??:??:??' cycleparms(2,ip)='??.? s' cycleparms(3,ip)='< 1R' cycleparms(4,ip)='> 0R' ENDIF ENDELSE ENDFOR ; loop through cameras ENDFOR ; loop through images in this cycle ; Blank the rest of the display, if necessary IF dispim THEN BEGIN IF ip LT 2*ninrow-1 THEN BEGIN FOR l=ip+1,2*ninrow-1 DO BEGIN p=l/2 q=l MOD 2 display(p*imgsize,q*imgsize)=blank ENDFOR ENDIF ENDIF ; Display the images all at once IF dispim THEN BEGIN cyclestart=cycleparms(1,0) IF STRPOS(cyclestart,'?') NE -1 THEN cyclestart=cycleparms(1,1) PRINT,cyclestart,FORMAT='(/"Ready at ",A,":",$)' IF GET_KBRD(1) EQ 'q' THEN GOTO,DONE mtv,display,res=120. ; Now annotate the images FOR l=0,ip DO BEGIN p=l/2 q=l MOD 2 ; Image camera, filter, and sequence number XYOUTS,scl(res=120.,p*imgsize+1) $ ,scl(res=120.,barwidth+(2-q)*imgsize-X_charheight/fac-1) $ ,cycleparms(0,l) $ ,/DEVICE,ALIGNMENT=0.0,CHARSIZE=csz/fac,COLOR=0.8*!D.N_COLORS ; and image time XYOUTS,scl(res=120.,p*imgsize+1) $ ,scl(res=120.,barwidth+(1-q)*imgsize+X_charheight/fac+3) $ ,cycleparms(1,l) $ ,/DEVICE,ALIGNMENT=0.0,CHARSIZE=csz/fac,COLOR=0.8*!D.N_COLORS ; and exposure time XYOUTS,scl(res=120.,p*imgsize+1) $ ,scl(res=120.,barwidth+(1-q)*imgsize+1) $ ,cycleparms(2,l) $ ,/DEVICE,ALIGNMENT=0.0,CHARSIZE=csz/fac,COLOR=0.8*!D.N_COLORS ; and maximum displayed value XYOUTS,scl(res=120.,(p+1)*imgsize-1) $ ,scl(res=120.,barwidth+(2-q)*imgsize-X_charheight/fac-1) $ ,cycleparms(3,l) $ ,/DEVICE,ALIGNMENT=1.0,CHARSIZE=csz/fac,COLOR=0.8*!D.N_COLORS ; and minimum displayed value XYOUTS,scl(res=120.,(p+1)*imgsize-1) $ ,scl(res=120.,barwidth+(1-q)*imgsize+1) $ ,cycleparms(4,l) $ ,/DEVICE,ALIGNMENT=1.0,CHARSIZE=csz/fac,COLOR=0.8*!D.N_COLORS ENDFOR ENDIF ENDFOR ; Close associated files IF imsave THEN BEGIN close_qlfiles,year,month,day,qlunits FOR jsave=0,9 DO $ IF savunit(jsave) NE 0 THEN BEGIN CLOSE,savunit(jsave) FREE_LUN,savunit(jsave) ENDIF ENDIF DONE: END