PRO dropcals,n,t,files,filtlist,typelist,texplist ; This little routine removes calibration images from the image ; list sequence to simplify cycle identification. imgcateg,files,filtlist,typelist,texplist calfiles=WHERE(typelist EQ 2) 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")' imgcateg,files,filtlist,typelist,texplist RETURN END ;--------------------------------------------------------------------- PRO getcycles,n,t,filtlist,typelist,texplist,cycle,cyclen,debug=debug ; 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 ; If the above technique doesn't work, we probably have continuous ; image acquisition without dropping back to DOS. Try another tack. IF (((N_ELEMENTS(cycle) EQ 2) AND (MAX(cycle-[-1,n]) EQ 0)) OR $ (N_ELEMENTS(cycle) LT n/10)) THEN BEGIN flag=0 cycle=tryshift(filtlist,typelist,texplist,flag,debug=debug) ENDIF IF cycle(0) NE 0 $ ; first cycle SHOULD start at image 0 THEN cycle=[0,cycle] IF N_ELEMENTS(flag) EQ 0 THEN flag=0 IF flag AND NOT KEYWORD_SET(debug) THEN $ PRINT,'Had a spot of trouble sorting out the cycles, you know' cyclen=(cycle-SHIFT(cycle,1))>0 cycle=cycle(0:N_ELEMENTS(cycle)-2) cyclen=cyclen(1:*) RETURN END ;--------------------------------------------------------------------- PRO opnsaveu,yr,month,day,nwsave,wsave,cwl,rpb,cfsave,rayleighs $ ,bgsubtract,savunit,save00,save01,save02,save03,save04 $ ,save10,save11,save12,save13,save14 @isitdos FOR jsave=0,nwsave-1 DO BEGIN tmp=WHERE(ABS(wsave(jsave)-cwl) LT rpb/2.) tmp=tmp(0) cfsave(tmp)=1 fname=saveroot+dd+STRING(yr,month,day,tmp/5,tmp MOD 5 $ ,FORMAT='(3I2.2,2I1,".")') CASE 1 OF (NOT rayleighs): fname=fname+'cDN' (rayleighs AND (NOT bgsubtract)): fname=fname+'ray' bgsubtract: fname=fname+'bgs' ENDCASE ; (1) OPENW,tt,fname,/GET_LUN savunit(tmp)=tt CASE tmp OF 0: save00=ASSOC(savunit(0),BYTARR(66048)) 1: save01=ASSOC(savunit(1),BYTARR(66048)) 2: save02=ASSOC(savunit(2),BYTARR(66048)) 3: save03=ASSOC(savunit(3),BYTARR(66048)) 4: save04=ASSOC(savunit(4),BYTARR(66048)) 5: save10=ASSOC(savunit(5),BYTARR(66048)) 6: save11=ASSOC(savunit(6),BYTARR(66048)) 7: save12=ASSOC(savunit(7),BYTARR(66048)) 8: save13=ASSOC(savunit(8),BYTARR(66048)) 9: save14=ASSOC(savunit(9),BYTARR(66048)) ENDCASE ; (tmp) ENDFOR ; (jsave) RETURN END ;--------------------------------------------------------------------- ; NEWQL.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. ; Modified 1994/03/29 to be compatible with filename formats under DOS. ; PRO mnewql,year,month,day,hr,min,sec $ ,batch=batch,intens=intens,diffbg=diffbg,debug=debug $ ,save=save,wlsave=wlsave,qlimageplot=qlimageplot COMMON KI_image,h,im COMMON ffstuff,ffimg IF N_PARAMS() LT 6 THEN BEGIN MESSAGE,'Must enter year,month,day,hr,min,sec at least!' RETURN ENDIF setstart=SYSTIME(1) ; Set up OS-specific items @isitdos ; create byte-valued input parameters cam=0B & filt=cam s='' !ORDER=1 ; uncover screen where text is to be printed IF !D.NAME EQ 'X' THEN WDELETE,0 ; If interactive operation desired, prompt user and get input IF NOT KEYWORD_SET(batch) THEN batch=0B fac=1. resp='' IF NOT batch THEN BEGIN PRINT,'Do you want to display images as well as save them? ([Y]/N)' READ,resp dispim=(STRUPCASE(resp) NE 'N') ENDIF ELSE dispim=0B ; (NOT batch) 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 ; (dispim) ; 1-cycle-per-hour option only works for Eureka data! PCsite=PoCasite(year,month,day) IF dispim AND (PCsite EQ 3) 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 ; (dispim) IF NOT batch THEN BEGIN PRINT,'Convert DN to intensity? ([Y]/N)' READ,resp rayleighs=(STRUPCASE(resp) NE 'N') ENDIF ELSE BEGIN IF KEYWORD_SET(intens) THEN rayleighs=intens ELSE rayleighs=1B ENDELSE ; (NOT batch) IF rayleighs THEN BEGIN IF NOT batch THEN BEGIN PRINT,'Present backgrounds as differential brightnesses (R/nm)? ([Y]/N)' READ,resp RperA=(STRUPCASE(resp) NE 'N') ENDIF ELSE BEGIN IF KEYWORD_SET(diffbg) THEN RperA=diffbg ELSE RperA=1B ENDELSE ; (NOT batch) ; 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 ; (bgsubtract) ENDIF ELSE BEGIN RperA=0 bgsubtract=0 ENDELSE ; (rayleighs) 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 ; (j) ENDIF ; (ceiling_set) ENDIF ELSE ceiling_set=0 ; (dispim) IF NOT batch THEN BEGIN PRINT,'Display debugging information? (Y/[N])' READ,resp debug=(STRUPCASE(resp) EQ 'Y') ENDIF ; (NOT batch) IF NOT batch THEN BEGIN 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 ; (NOT short) ENDIF ELSE IF KEYWORD_SET(save) THEN imsave=save ELSE imsave=0B cfsave=BYTARR(10) savunit=BYTARR(10) IF NOT batch THEN BEGIN 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 ; (nwsave GT 0) ENDIF ELSE nwsave=0 ; (imsave) ENDIF ELSE BEGIN nwsave=N_ELEMENTS(wlsave) IF nwsave GT 0 THEN BEGIN wsave=wlsave nwsave=N_ELEMENTS(wsave) ENDIF ; (nwsave GT 0) ENDELSE ; (NOT batch) ; make lists of available files of desired type rdilist,2,year,month,day,hr,min,sec,n,t,files ; For production purposes, only run NEWQL if all data from the ; given day are available. IF (hr*60L+min)*60+sec EQ 0 THEN BEGIN mint=MIN(t/3600.,MAX=maxt) IF (mint GT 0.5) OR (maxt LT 23.5) THEN BEGIN MESSAGE,"Full day's data not available; not worth running" $ ,/INFORMATIONAL SETENV,'ACTION=skip' RETURN ENDIF ENDIF IF imsave THEN iniqlfil,qlunits IF KEYWORD_SET(qlimageplot) THEN initqlif,year,month,day ; ensure a 4-digit year for display; make reasonable assumptions yr=year IF year LT 1900 THEN BEGIN IF year LT 50 THEN year=year+2000 IF year GT 50 THEN year=year+1900 ENDIF ELSE yr=year MOD 100 ; (year LT 1900) ; create title for image plot IF dispim $ THEN wtitle='Polar Camera '+ymd2date(year,month,day,format='d$ n$ y$') ; 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 IF n GT 0 THEN BEGIN ; Reverse the indices of the time and file name arrays as returned ; by rdilist t=TRANSPOSE(t) files=TRANSPOSE(files) ext_list=extensio(files) ; Identify and delete calibration files from time and file name arrays, ; so no impossibly long cycles appear dropcals,n,t,files,filtlist,typelist,texplist ENDIF ; (n GT 0) IF n GT 0 THEN BEGIN getcycles,n,t,filtlist,typelist,texplist,cycle,cyclen,debug=debug ; If a cycle is too long, investigate. IF NOT batch THEN BEGIN IF MAX(cyclen)/fac GT 5 THEN BEGIN PRINT,'Maximum cycle length is '+STRTRIM(MAX(cyclen),2)+'; check it out!' STOP ENDIF ; (MAX(cyclen) GT 5) ENDIF ; (NOT batch) ; 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. ; ; NOTA BENE: This won't work for continuous data acquisition, such ; as that done at Rabbit lake and Spy Hill. Only good for Eureka! 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 ; (short) ENDIF ; (dispim) maxcyclen=MAX(cyclen) cycleparms=STRARR(5,2*maxcyclen); to store image annotation strings fac=1 IF maxcyclen GT 5 THEN fac=2 ; IF dos THEN scrwid=4. ELSE scrwid=5. ; IF dispim AND (maxcyclen GT scrwid*fac) THEN BEGIN ; avoid display problems ; longcyc=STRTRIM(maxcyclen,2) ; with long cycles ; maxninrow=STRTRIM(FIX(scrwid*fac),2) ; minnewfac=STRTRIM(CEIL(maxcyclen/scrwid),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 ; (dispim AND (maxcyclen GT scrwid*fac)) ; 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 ; (dispim) ; Load the pseudo-flat-field images generated at Eureka, for later ; application to the images before display. getffimg,year,month,day ; Load the Rayleigh-second/DN values for converting to intensities IF rayleighs THEN BEGIN getsens,year,month,day,sens0,sens1,nointeractive=batch sens=[sens0,sens1] ; Load the center wavelengths and rectangular passbands for background subtraction rdcwlrpb,cwl,rpb,unitstring ENDIF ELSE unitstring='DN' ; (rayleighs) ; If images are to be saved, open the associated files to hold them IF nwsave GT 0 THEN BEGIN opnsaveu,yr,month,day,nwsave,wsave,cwl,rpb,cfsave,rayleighs $ ,bgsubtract,savunit,save00,save01,save02,save03,save04 $ ,save10,save11,save12,save13,save14 seqnum=INTARR(10) ENDIF ; (nwsave GT 0) ; Create the format string for describing background subtractions IF NOT KEYWORD_SET(debug) THEN debug=0B IF debug THEN $ bgform='("Image ",I1,", camera ",I1,";"/" subtract image ",I1,", camera ",I1)' IF dispim THEN BEGIN 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/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 dos THEN username='CNSR3' ELSE 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 IF dos THEN BEGIN datetime=SYSTIME() date=STRMID(datetime,4,7)+STRMID(datetime,20,4) ENDIF ELSE BEGIN SPAWN,'date -u',datetime date=STRMID(datetime,4,7)+STRMID(datetime,24,4) ENDELSE time=STRMID(datetime,11,9)+' UT' XYOUTS,scl(res=120.,wwidth-X_charwidth) $ ,scl(res=120.,wheight-2.5*X_charheight) $ ,time,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=congrid(BINDGEN(256,barwidth),ninrow*imgsize,barwidth,/minus_one) 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 ; (dispim) ; Decide where Polar Camera was located PCsite=PoCasite(year,month,day) ; Calculate the azimuth and zenith angle images for each camera col=BINDGEN(256)#REPLICATE(1B,1,256) row=TRANSPOSE(col) cr2az,year,month,day,0,0,col,row,az0,za0,/deg,site=PCsite cr2az,year,month,day,1,0,col,row,az1,za1,/deg,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 ; (dispim) ; 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,droot+dd+'sig0.EU',/GET_LUN s0=0L READU,s0unit,s0 sig0=LONARR(s0) READU,s0unit,sig0 CLOSE,s0unit FREE_LUN,s0unit OPENR,s1unit,droot+dd+'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,droot+dd+'triangl0.EU',/GET_LUN triangles0=ASSOC(t0unit,LONARR(3,115213)) OPENR,t1unit,droot+dd+'triangl1.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 ; ((lw NE 1) AND (lw NE 6) AND (lw NE 9)) ENDFOR ; (lw) ENDELSE ; (approx) pp=FINDGEN(100)*0.01+0.005 ENDIF ; (bgsubtract) ; 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] procstart=SYSTIME(1) PRINT,procstart-setstart,FORMAT='("Setup took ",F5.1," seconds.")' FOR k=0,N_ELEMENTS(cycle)-1 DO BEGIN ; create (or refresh) arrays for temporary storage IF dispim OR (nwsave GT 0) THEN BEGIN keep=INTARR(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)) keepaz=FLTARR(cyclen(k),2) ; solar azimuths ENDIF ; (dispim OR (nwsave GT 0)) IF NOT dispim $ THEN PRINT,k,cycle(k),t(cycle(k)),cyclen(k) $ ,FORMAT='("Cycle ",I3," starts at pair "' $ +',I4," (",F6.0," s), and includes "' $ +',I0," 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) kihsize=SIZE(KIh) good_header=(kihsize(N_ELEMENTS(kihsize)-2) EQ 8) ENDIF ELSE good_header=0B ; If 'gethd' returned a valid header, carry on. IF good_header THEN BEGIN FOR l=0,3 DO itime(l)=KIh.misc.tm.(l) exp=FIX(ROUND(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=(N_ELEMENTS(width) EQ 0) 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 ; (NOT first_time) ;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 ; (NOT first_time) IF debug THEN PRINT,'Dark frames being loaded!' dkp=WHERE(dkexp EQ exp,ndkp) ; find exposure IF ndkp EQ 0 $ THEN MESSAGE,'Weird dark exposure!' $ ,/INFORMATIONAL $ ELSE BEGIN ; store the dark frame dkp=dkp(0) darks(dkp*width,j*length)=MEDIAN(im,5) ENDELSE ;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 N_ELEMENTS(width) GT 0 THEN BEGIN dwidth=width dlength=length ENDIF ELSE BEGIN dwidth=256 dlength=256 ENDELSE ; (N_ELEMENTS(width) GT 0) dkim=darks(pd*dwidth:(pd+1)*dwidth-1,j*dlength:(j+1)*dlength-1) im=(im-dkim)>0 ENDIF ; (pd(0) NE -1) ENDELSE ; (exptype EQ 5) ; 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 ; ((KIh.exp.sbin NE 2) OR (KIh.exp.pbin NE 2)) ; Apply pseudo-flat-field correction to dark-subtracted image IF (exptype NE 3) AND (exptype NE 5) THEN BEGIN filt=KIh.misc.fw im=FIX(ROUND(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=FIX(ROUND(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=FIX(ROUND(im/rpb(1))) IF (j EQ 1) AND (filt EQ 1) THEN im=FIX(ROUND(im/rpb(6))) IF (j EQ 1) AND (filt EQ 4) THEN im=FIX(ROUND(im/rpb(9))) ENDIF ; (RperA) ENDIF ; (rayleighs) ENDIF ; ((exptype NE 3) AND (exptype NE 5)) ; Store processed images for subsequent background subtraction or other use IF dispim OR (nwsave GT 0) THEN BEGIN keepexst(p,j)=1B ; If we're making processed files, despike the images before saving them IF nwsave GT 0 THEN im=cleanmod(im,3,2,3,/quiet) keep(p*256,j*256)=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 sunazel,KIh.misc.dt.year,KIh.misc.dt.month $ ,KIh.misc.dt.day,KIh.misc.tm.hour $ ,KIh.misc.tm.minute,KIh.misc.tm.second $ ,solaz,solel keepaz(p,j)=solaz keephb(0,ip)=h ENDIF ; (dispim OR (nwsave GT 0)) ; 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) savqldat,j,filt,qlunits,itime,exp,zdata ENDIF ; (((exptype NE 3) AND (exptype NE 5)) AND imsave) ; Save image in quick-look image file for later plotting, if ; appropriate IF KEYWORD_SET(qlimageplot) THEN sto2qlif,unitstring ; 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 ; ((j EQ 0) AND (filt EQ 1)) IF (j EQ 1) AND (filt EQ 1) THEN BEGIN bg820=im IF NOT RperA THEN bg820=bg820/rpb(6) ENDIF ; ((j EQ 1) AND (filt EQ 1)) IF (j EQ 1) AND (filt EQ 4) THEN BEGIN bg714=im IF NOT RperA THEN bg714=bg714/rpb(9) ENDIF ; ((j EQ 1) AND (filt EQ 4)) ENDIF ; (bgsubtract) ENDIF ; (valid_header) ENDFOR ; (j) ENDFOR ; (i) ; For now we only do approximate background subtraction IF bgsubtract $ THEN darkcyc=(TOTAL(keeptype EQ 5) EQ cyclen(k)*2) $ ELSE darkcyc=0B ; 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,droot+dd+'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 ; (j EQ 0) ; 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,year,month,day,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=pcentile(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 ; (approx AND same(fptr)) ; IF dispim THEN BEGIN ; WSET,3 ; imax=pcentile(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 ; (dispim) ; bgs(fptr)=bg ; bgmade(fptr)=1 ; ENDELSE ; (bgmade(fptr)) ; STOP,'check out im,bg; .con to continue' ; keep(i*256:(i+1)*256-1,0:255)=(im-bg)>0 ; ENDIF ; (NOT ((fptr EQ 1) OR ; ; (fptr EQ 6) OR (fptr EQ 9)) ; ENDIF ; (keepexst(i,j)) ; ; 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 ; (j) ; ENDFOR ; (i) ; CLOSE,bgunit ; FREE_LUN,bgunit ; ENDIF ; (bgsubtract AND (NOT darkcyc)) IF dispim OR (nwsave GT 0) THEN BEGIN 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 imsig=im(sig0) $ ELSE imsig=im(sig1) lmaxi=pcentile(imsig,0.995) hlp,j,imsig ENDIF ELSE lmaxi=1 ; (MAX(im) GT 1) ENDIF ELSE lmaxi=topnum(j) ; (NOT ceiling_set) ; 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=pcentile(im(base),0.5)>0 hlp,im(base) ENDIF ELSE lmaxbh=0 ENDIF ELSE BEGIN ; otherwise take lowest level in dark frame IF j EQ 0 $ THEN imsig=im(sig0) $ ELSE imsig=im(sig1) lmaxbh=pcentile(imsig,0.001) hlp,j,imsig ENDELSE ; (NOT shut) ; 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 mask=mask0 ELSE mask=mask1 IF (((KIh.exp.sbin EQ 2) AND (KIh.exp.pbin EQ 2)) $ AND (fac LE 1.)) $ THEN im=im*mask $ ELSE im=REBIN(im*mask,imgsize,imgsize,/SAMPLE) display(ic*imgsize,j*imgsize)= $ ROTATE(BYTSCL(im,MIN=lmaxbh,MAX=lmaxi),5) ENDIF ; (dispim) ; 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) ; only save values up to the 99.5th percentile p995=pcentile(im,0.995) hlp,im ; bbin is the number such that 2**(bbin-1) LE p995 LE 2**bbin bbin=BYTE((CEIL(ALOG((MAX((im>0)0) hb(290)=bbin hb(291)=1B ; save im/(2.**bbin) in a 1-byte space outrec=[hb,REFORM(BYTE(ROUND(im/FLOAT(2^bbin))),65536,1)] 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 ; (fptr) seqnum(fptr)=seqnum(fptr)+1 ENDIF ; (((exptype NE 3) AND ; (exptype NE 5)) AND cfsave(fptr) ; 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 ; (dispim) 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 ; (dispim) ENDELSE ; (keepexst(ic,j)) ENDFOR ; (j) loop through cameras ENDFOR ; (i) 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 ; (l) ENDIF ; (ip LT 2*ninrow-1) ENDIF ; (dispim) ; 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 ; and direction to the Sun solaz=keepaz(p,q) IF solaz EQ 0. THEN sunptr='' ELSE BEGIN IF solaz LT 180. THEN BEGIN sunptr='Sun>' ENDIF ELSE BEGIN sunptr='