Recipe for the coincidence loss correction from Poole et al. (2008)

coincidence_loss_corr.pro

;compute the size of the images (you can also do this manually rather than calling these keywords from the header)
nxpix_uw1_ext1 = sxpar(hima_sk_uw1_ext1,'NAXIS1')
nypix_uw1_ext1 = sxpar(hima_sk_uw1_ext1,'NAXIS2')

;make a new image file to save the correction factors
coicorr_uw1_ext1 = fltarr(nxpix_uw1_ext1,nypix_uw1_ext1)
;make a new image file to save the corrected image
ima_sk_coicorr_uw1_ext1 = fltarr(nxpix_uw1_ext1,nypix_uw1_ext1)

;define the variables from Poole et al. 2008
alpha = 0.9842
ft = 0.0110329 ;0.0110322?
a1 = 0.0658568
a2 = -0.0907142
a3 = 0.0285951
a4 = 0.0308063

for i=0,nxpix_uw1_ext1-1 do begin
   for j=0,nypix_uw1_ext1-1 do begin
      if ((i LE 4) || (i GE (nxpix_uw1_ext1-4)) || (j LE 4) || (j GE (nypix_uw1_ext1-4))) then begin
         ;UVW1
         coicorr_uw1_ext1[i,j] = 0
         ima_sk_coicorr_uw1_ext1[i,j] = 0
      endif else begin
         xpixmin = i-4
         xpixmax = i+4
         ypixmin = j-4
         ypixmax = j+4
         ;UVW1
         ima_UVW1sum = total(ima_norm_uw1_ext1[xpixmin:xpixmax,ypixmin:ypixmax])

         xvec_UVW1 = ft*ima_UVW1sum

 
         fxvec_UVW1 = 1 + (a1*xvec_UVW1) + (a2*xvec_UVW1*xvec_UVW1) + (a3*xvec_UVW1*xvec_UVW1*xvec_UVW1) + (a4*xvec_UVW1*xvec_UVW1*xvec_UVW1*xvec_UVW1)

 
         Ctheory_UVW1 = - alog(1-(alpha*ima_UVW1sum*ft))/(alpha*ft)
        
 
         coicorr_uw1_ext1[i,j] = Ctheory_UVW1*fxvec_UVW1/ima_UVW1sum
        
 
         ima_sk_coicorr_uw1_ext1[i,j] = coicorr_uw1_ext1[i,j]*ima_norm_uw1_ext1[i,j]
      endelse
   endfor
endfor