; ; Polarization calibration routine for a single exposure ; time of poise images ; Opal normalization is not yet included ; X matrix is included ; dome = 1 if ( dome ) then begin ;{ window,xsize=678,ysize=678 endif ;} close,1 & close,13 loadct,1 lfilenm = 'list' filenm=' ' close,13 & openr,13,lfilenm a=0. b=0. ; Read in FTS image readf,13,filenm img = readfits(filenm,hdr) xsz = fix(sxpar (hdr , 'NAXIS1')) ysz = fix(sxpar (hdr , 'NAXIS2')) ; Define these arrays now that dimensions are known (right) rimg=fltarr(xsz,ysz) row=fltarr(ysz) ; Fix the byte order byteorder,img ; Correct for integer wraparound img = long(img) img = (img lt 0 ) * ( img + 32768 + 32768) + (img ge 0)*(img) ; Flip top to bottom for IDL img = rotate(img,7) ; Record where intensity is too high for future use a1 = where(img gt 40000) ; Compute an average value for each row and turn that into an image ; Add a fraction of the average row image to eliminate horizontal ; offset in the CCD which is proportional to total signal of row ; Factor of 1/100. is just a guess for jj=0,ysz-2 do begin row(jj)=total(img(*,jj+1)) row(jj)=row(jj)/(1.*xsz) rimg(*,jj)=row(jj) end ; Rather than using a dark, take bias from a bottom row sbias=total(img(4,*))/(1.*ysz) img = img - sbias simg1 = rebin(img,xsz/3,ysz/3 ) tv,simg1 img1=img+rimg/100. ; Here is where opal is divided into data ; ; ; ; Scaling from here on will be changed simg1 = rebin(img1,xsz/3,ysz/3 ) tv,simg1 ; Read in FTS image readf,13,filenm img = readfits(filenm,hdr) ; Fix the byte order byteorder,img ; Correct for integer wraparound img = long(img) img = (img lt 0 ) * ( img + 32768 + 32768) + (img ge 0)*(img) ; Flip top to bottom for IDL img = rotate(img,7) ; Record where intensity is too high for future use a2 = where(img gt 40000) ; Compute an average value for each row and turn that into an image ; Add a fraction of the average row image to eliminate horizontal ; offset in the CCD which is proportional to total signal of row ; Factor of 1/100. is just a guess for jj=0,ysz-2 do begin row(jj)=total(img(*,jj+1)) row(jj)=row(jj)/(1.*xsz) rimg(*,jj)=row(jj) end ; Rather than using a dark, take bias from a bottom row sbias=total(img(4,*))/(1.*ysz) img = img - sbias simg2 = rebin(img,xsz/3,ysz/3 ) tv,simg2 img2=img+rimg/100. ; Here is where opal is divided into data ; Scaling from here on may be changed simg2 = rebin(img2,xsz/3,ysz/3 ) tv,simg2 ; Read in FTS image readf,13,filenm img = readfits(filenm,hdr) ; Fix the byte order byteorder,img ; Correct for integer wraparound img = long(img) img = (img lt 0 ) * ( img + 32768 + 32768) + (img ge 0)*(img) ; Flip top to bottom for IDL img = rotate(img,7) ; Record where intensity is too high for future use a3 = where(img gt 40000) ; Compute an average value for each row and turn that into an image ; Add a fraction of the average row image to eliminate horizontal ; offset in the CCD which is proportional to total signal of row ; Factor of 1/100. is just a guess for jj=0,ysz-2 do begin row(jj)=total(img(*,jj+1)) row(jj)=row(jj)/(1.*xsz) rimg(*,jj)=row(jj) end ; Rather than using a dark, take bias from a bottom row sbias=total(img(4,*))/(1.*ysz) img = img - sbias simg3 = rebin(img,xsz/3,ysz/3 ) tv,simg3 img3=img+rimg/100. ; Here is where opal is divided into data ; Scaling from here on may be changed simg3 = rebin(img3,xsz/3,ysz/3 ) tv,simg3 ; Read in FTS image readf,13,filenm img = readfits(filenm,hdr) ; Fix the byte order byteorder,img ; Correct for integer wraparound img = long(img) img = (img lt 0 ) * ( img + 32768 + 32768) + (img ge 0)*(img) ; Flip top to bottom for IDL img = rotate(img,7) ; Record where intensity is too high for future use a4 = where(img gt 40000) ; Compute an average value for each row and turn that into an image ; Add a fraction of the average row image to eliminate horizontal ; offset in the CCD which is proportional to total signal of row ; Factor of 1/100. is just a guess for jj=0,ysz-2 do begin row(jj)=total(img(*,jj+1)) row(jj)=row(jj)/(1.*xsz) rimg(*,jj)=row(jj) end ; Rather than using a dark, take bias from a bottom row sbias=total(img(4,*))/(1.*ysz) img = img - sbias simg4 = rebin(img,xsz/3,ysz/3 ) tv,simg4 img4=img+rimg/100. ; Here is where opal is divided into data ; Scaling from here on may be changed simg4 = rebin(img4,xsz/3,ysz/3 ) tv,simg4 close,13 intensity=img1+img2+img3+img4 sint = rebin(intensity,xsz/3,ysz/3 ) ; ; Compute matrix form of Q and U, and invert to see how well we did x = fltarr(3,3) x = [[1.,0.04,0.045],[0.,0.627,0.315],[0.,0.32,-0.632]] xinv = invert(x) qmeas = img1-img2-img3+img4 umeas = img1+img2-img3-img4 q = xinv(1,1)*qmeas+xinv(2,1)*umeas u = xinv(1,2)*qmeas+xinv(2,2)*umeas pol = sqrt(q*q+u*u) spol = rebin(pol,xsz/3,ysz/3 ) tv,spol/64. stop print, 'masked pb' ; Read in 1/3 scale mask. See pb.pro for centers and call to make it mask=spol openr,1,'/thisbe/d/elmore/Eclipse/mask' readu,1,mask close,1 tv,spol*mask/256. stop sob=sobel(pol) ssob = rebin(sob,xsz/3,ysz/3 ) tv,ssob/32. stop mad=sob/(pol+250.) smad = rebin(mad,xsz/3,ysz/3 ) tv,smad*128. stop print,"log intensity" tv,bytscl(alog10(sint),2,6) stop print,"log pb" tv,bytscl(alog10(spol),1,5) stop print,"log pb/intensity" ; next good for long exposure ; tv,bytscl(alog10(spol/sint),-1.8,-.1) ; next good for medium tv,bytscl(alog10(spol/sint),-3.4,.4) stop end