PRO reducecals,yr,mo,day,cam,R_s_per_DN,verbose=vb,notv=notv ; find calibration images and dark images for given day ; read calibration images and dark images ; perform dark subtraction ; read flat-field images and do flat-field corrections ; determine center and approximate radius of PLBS images ; average camera signal over central 90% of images ; read in filter throughput values ; calculate Rayleigh-seconds/DN for various channels, using all ; available images for each channel, and adding 0.045 s to nominal ; exposure times ; return array of Rayleigh-seconds/DN values dos=(!VERSION.OS EQ 'windows') ; Initialize OS-specific items @isitdos ; find calibration images and dark images for given day rdilist,cam,yr,mo,day,0,0,0,n,t,f IF n EQ 0 THEN RETURN cal=WHERE(STRPOS(f,calname) NE -1,ncal) dk=WHERE(STRPOS(f,dkname) NE -1,ndk) IF ncal EQ 0 THEN RETURN ; first get header bytes and identify the real calibration images seq=BYTARR(ncal) FOR j=0,ncal-1 DO BEGIN rdkihd,f(cal(j)),hb h=gethd(hb) seq(j)=(h.misc.exp_scale*h.exp.exposure GT 25) ENDFOR IF TOTAL(seq) EQ 0 THEN RETURN cal=cal(WHERE(seq)) ; calibration images take longer than 25 ms ncal=N_ELEMENTS(cal) ; display calibration image file names PRINT,'Camera '+STRTRIM(cam,2)+' calibration files: FOR j=0,ncal-1 DO PRINT,f(cal(j)) ; not enough files for a real calibration sequence IF ncal LT 8 THEN RETURN ; read calibration images and dark images ; Now load the dark frames nearest in time to the calibration images firstcal=MIN(t(cal),MAX=lastcal) ; times of first and last cal images dkfr=INTARR(512,768) ; array to hold dark frames ccdt=FLTARR(2,3) ; array to hold CCD temps ; find last 1-s dark frame taken before first calibration image bef1s=MAX(WHERE(t(dk) LT firstcal) AND (STRPOS(f(dk),'001.0') NE -1)) ; read in the header bytes and image rdkimg,f(dk(bef1s)),hb,dkim ; store the image for later use dkfr(0,0)=MEDIAN(dkim,3) ; decode the header bytes h=gethd(hb) ; and store the CCD temperature for later reference ccdt(0,0)=h.misc.temp.tcf.mean/10. ; now repeat the process for the first 1-s dark frame after calibration aft1s=MIN(WHERE(t(dk) GT lastcal) AND (STRPOS(f(dk),'001.0') NE -1)) rdkimg,f(dk(aft1s)),hb,dkim dkfr(256,0)=MEDIAN(dkim,3) h=gethd(hb) ccdt(1,0)=h.misc.temp.tcf.mean/10. ; now 10-s dark frames bef10s=MAX(WHERE(t(dk) LT firstcal) AND (STRPOS(f(dk),'010.0') NE -1)) rdkimg,f(dk(bef10s)),hb,dkim dkfr(0,256)=MEDIAN(dkim,3) h=gethd(hb) ccdt(0,1)=h.misc.temp.tcf.mean/10. aft10s=MIN(WHERE(t(dk) GT lastcal) AND (STRPOS(f(dk),'010.0') NE -1)) rdkimg,f(dk(aft10s)),hb,dkim dkfr(256,256)=MEDIAN(dkim,3) h=gethd(hb) ccdt(1,1)=h.misc.temp.tcf.mean/10. ; and finally 60-s dark frames bef60s=MAX(WHERE(t(dk) LT firstcal) AND (STRPOS(f(dk),'060.0') NE -1)) rdkimg,f(dk(bef60s)),hb,dkim dkfr(0,512)=MEDIAN(dkim,3) h=gethd(hb) ccdt(0,2)=h.misc.temp.tcf.mean/10. aft60s=MIN(WHERE(t(dk) GT lastcal) AND (STRPOS(f(dk),'060.0') NE -1)) rdkimg,f(dk(aft60s)),hb,dkim dkfr(256,512)=MEDIAN(dkim,3) h=gethd(hb) ccdt(1,2)=h.misc.temp.tcf.mean/10. ; Before actually processing, load the path names for the throughput ; data tppath=STRARR(2,5) tppath(0,0)=tproot+dd+'829c' tppath(0,1)=tproot+dd+'608c' tppath(0,2)=tproot+dd+'630c' tppath(0,3)=tproot+dd+'668c' tppath(0,4)=tproot+dd+'520c' tppath(1,0)=tproot+dd+'835c' tppath(1,1)=tproot+dd+'820c' tppath(1,2)=tproot+dd+'785c' tppath(1,3)=tproot+dd+'558c' tppath(1,4)=tproot+dd+'714c' ; Now go through the real calibration data files one at a time wtsum=FLTARR(5) sumwt=wtsum R_s_per_DN=wtsum FOR j=0,ncal-1 DO BEGIN rdkimg,f(cal(j)),hb,im h=gethd(hb) nsat=TOTAL(im EQ 16383) IF KEYWORD_SET(vb) THEN PRINT,nsat,FORMAT='(/I0.0," saturated pixels")' IF nsat LT 100 THEN BEGIN texp=h.misc.exp_scale*h.exp.exposure/1000 CASE texp OF 1: j2=0 10: j2=1 60: j2=2 ENDCASE tcf=h.misc.temp.tcf.mean/10. IF ABS(ccdt(0,j2)-tcf) LT ABS(ccdt(1,j2)-tcf) THEN j1=0 ELSE j1=1 ; subtract appropriate dark frame im=(im-dkfr(256*j1:256*(j1+1)-1,256*j2:256*(j2+1)-1))>0 ; load and apply appropriate flat-field correction image cam=h.misc.cam filt=h.misc.fw flatname=STRING(flatroot,dd,cam,filt $ ,FORMAT='(A,A1,"ffimg_",2I1,".dat")') OPENR,unit,flatname,/GET_LUN ffim=BYTARR(256,256) READU,unit,ffim CLOSE,unit FREE_LUN,unit im=FIX(ROUND(im*FLOAT(ffim)/100.)) ; determine position and size of PLBS image meancent,im,clbs,rlbs,npix,radius=radlbs,notv=notv ; IF the PLBS radius is reasonable, carry on IF radlbs GT 10 THEN BEGIN ; restrict averaging to central 90% of PLBS image rad90=SQRT(0.9)*radlbs one2two,LINDGEN(65536L),im,cmask,rmask mk,mask mask(WHERE(SQRT((cmask-clbs)^2+(rmask-rlbs)^2) LE rad90))=1 avgdn=TOTAL(mask*im)/TOTAL(mask) ; read in throughput data for this filter tpfile=getfile(tppath(cam,filt)+dd+'thruput') ndata=N_ELEMENTS(tpfile)-1 tptemp=FLTARR(ndata) tpkr=tptemp FOR k=1,ndata DO BEGIN tptemp(k-1)=FLOAT(getwrd(tpfile(k),0)) tpkr(k-1)=FLOAT(getwrd(tpfile(k),1)) ENDFOR tp=spline(tptemp,tpkr,h.misc.temp.tfw.mean/10.) ; Now calculate sensitivity sens=1000.*tp*(texp+0.045)/avgdn wtsum(filt)=wtsum(filt)+(texp+0.045)*sens sumwt(filt)=sumwt(filt)+(texp+0.045) IF KEYWORD_SET(vb) THEN PRINT,cam,filt,texp,tp,avgdn $ ,FORMAT='("cam ",I1," filt ",I1,F6.1," s ",F7.2," kR ",F7.1," DN")' IF KEYWORD_SET(vb) THEN PRINT,sens,FORMAT='(F6.1," Rayleigh-s/DN")' ENDIF ENDIF ENDFOR FOR j=0,4 DO IF sumwt(j) GT 0. THEN R_s_per_DN(j)=wtsum(j)/sumwt(j) ; Save results to a file if they mean anything IF N_ELEMENTS(WHERE(R_s_per_DN)) GE 2 THEN BEGIN ; First see whether results are reasonable, i.e., within 25% ; of other results nzR=WHERE(R_s_per_DN) IF cam EQ 0 THEN frac=R_s_per_DN/[128.8,61.6,71.3,61.9,90.0] $ ELSE frac=R_s_per_DN/[112.1,102.8,86.7,68.4,53.9] ok=(TOTAL(ABS(frac-1.) LT 0.25) EQ N_ELEMENTS(nzR)) IF ok THEN BEGIN savename=STRING(calsroot,dd,yr,mo,day,cam,cal_suffix $ ,FORMAT='(A,A1,3I2.2,".",I1,A)') OPENW,unit,savename,/GET_LUN PRINTF,unit,R_s_per_DN,FORMAT='(5F6.1)' CLOSE,unit FREE_LUN,unit ENDIF ELSE BEGIN PRINT,yr,mo,day,cam $ ,FORMAT='("Wonky calibration data on ",3I2.2," for camera ",I1)' ENDELSE ENDIF RETURN END