 pro Puffy_Events_star_region_no_double_count2

;This script will attempt to locate the events on the H2RG-015 that are attributed
;to electron-hole pairs in the CdZnTe 800 um substrate

;The program first finds all pixels that have a count greater than THRESHOLD and flags them.  Next it
;looks at symmetric regions (star-shaped regions of radius r, where the star approaches a circle as r
;goes to infinity) and checks to see if all the pixels in that region are greater than THRESHOLD.  If so
;that event is logged in logfile and all of the pixels are set to zero so that it will not be double counted.


;CONSTANTS******************************************************.
;Force some to be floating point integers

numpixrd=2048.*2048.             ;Size of device
Array_Dimension=2048.          ;Array Dimension
THRESHOLD=10               ;THRESHOLD to flag pixel as Cosmic Ray
Gauss_width=14                 ;Width of Area to use for Gaussian Fitting
min_sig=float(.6)                ;Minimum standard deviation required to flag event
read_num=123                 ;Number of read to use for ramp difference
large_cosm_ind_x=[0]          ;X Pixel value Array to be Appended
large_cosm_ind_y=[0]          ;Y Pixel value Array to be Appended

;FILENAMES******************************************************

Directory='dktest.16May03'

file_number_1=5
file_number_str_1=string(file_number_1,format='(I1)')

file_number_2=2
file_number_str_2=string(file_number_2,format='(I1)')

temp=37
temp_str=string(temp,format='(I2)')

read_num_tot=125
readsstr=string(read_num_tot,format='(I3)')

file_name_1='dark_'+temp_str+'K_'+readsstr+'_min_1_0'+file_number_str_1+'.fits'
file_name_2='dark_'+temp_str+'K_'+readsstr+'_min_1_0'+file_number_str_2+'.fits'

logfile='c:\H2RG_15_EVENTS\star.txt'
;************************************************************
;Open the file that statistics will be written to and make table header

openw,1,logfile,/append

;Enter the directory that holds the dark current files for the device
cd,'\\Rabbit\raid1\Rockwell4\H2RG-015-5.0mu\cold1\'+Directory

;Vectors for successive reads for first image
fits_read,file_name_1,x11,first=long(read_num*numpixrd),last=long((read_num+1)*numpixrd)-1

;Vectors for successive reads for second image
fits_read,file_name_2,x21,first=long(read_num*numpixrd),last=long((read_num+1)*numpixrd)-1

;For last read, generate image to use for overall events
;if read_num eq (read_num_tot-3) then begin

im1=reform(float(x11),2048.,2048.)
im2=reform(float(x21),2048.,2048.)
diff=float(im1-im2)
diff(where(diff lt 0))=0.

filename='C:\H2RG_15_EVENTS\'+Directory+'_final_read_diff_'+'_'+temp_str+'K_'+readsstr+'Reads_'+file_number_str_1+'_'+file_number_str_2+'.fits'
Index_string=Directory+'_'+temp_str+'K_'+readsstr+'_'+file_number_str_1+'_'+file_number_str_2
writefits,filename,diff

star_vector_5=where(((shift(diff,-5) gt THRESHOLD) and (shift(diff,-4) gt THRESHOLD) and (shift(diff,-3) gt THRESHOLD) and (shift(diff,-2) gt THRESHOLD) and (shift(diff,-1) gt THRESHOLD) $
                 and (diff gt THRESHOLD) and (shift(diff,1) gt THRESHOLD) and (shift(diff,2) gt THRESHOLD) and (shift(diff,3) gt THRESHOLD) and (shift(diff,4) gt THRESHOLD) and (shift(diff,5) gt THRESHOLD) $
                 and (shift(diff,-2044) gt THRESHOLD) and (shift(diff,-2045) gt THRESHOLD) and (shift(diff,-2046) gt THRESHOLD) $
                 and (shift(diff,-2047) gt THRESHOLD) and (shift(diff,-2048) gt THRESHOLD) and (shift(diff,-2049) gt THRESHOLD) $
                 and (shift(diff,-2050) gt THRESHOLD) and (shift(diff,-2051) gt THRESHOLD) and (shift(diff,-2052) gt THRESHOLD) $
                 and (shift(diff,2044) gt THRESHOLD) and (shift(diff,2045) gt THRESHOLD) and (shift(diff,2046) gt THRESHOLD) $
                 and (shift(diff,2047) gt THRESHOLD) and (shift(diff,2048) gt THRESHOLD) $
                 and (shift(diff,2049) gt THRESHOLD) and (shift(diff,2050) gt THRESHOLD) and (shift(diff,2051) gt THRESHOLD) and (shift(diff,2052) gt THRESHOLD) $
                 and (shift(diff,-4084) gt THRESHOLD) and (shift(diff,-4085) gt THRESHOLD) and (shift(diff,-4086) gt THRESHOLD) $
                 and (shift(diff,-4087) gt THRESHOLD) and (shift(diff,-4088) gt THRESHOLD)  $
                 and (shift(diff,4084) gt THRESHOLD) and (shift(diff,4085) gt THRESHOLD) and (shift(diff,4086) gt THRESHOLD) $
                 and (shift(diff,4087) gt THRESHOLD) and (shift(diff,4088) gt THRESHOLD) $
                 and (shift(diff,-6143) gt THRESHOLD) and (shift(diff,-6144) gt THRESHOLD) and (shift(diff,-6145) gt THRESHOLD) $
                 and (shift(diff,6143) gt THRESHOLD) and (shift(diff,6144) gt THRESHOLD) and (shift(diff,6145) gt THRESHOLD) $
                 and (shift(diff,-8191) gt THRESHOLD) and (shift(diff,-8192) gt THRESHOLD) and (shift(diff,-8193) gt THRESHOLD) $
                 and (shift(diff,8191) gt THRESHOLD) and (shift(diff,8192) gt THRESHOLD)) and (shift(diff,8193) gt THRESHOLD)  $
                 and (shift(diff,-10240) gt THRESHOLD) $
                 and (shift(diff,10240) gt THRESHOLD) $
                 and (shift(diff,-5) lt diff) and (shift(diff,-4) lt diff) and (shift(diff,-3) lt diff) and (shift(diff,-2) lt diff) and (shift(diff,-1) lt diff) $
                 and (shift(diff,1) lt diff) $
                 and (shift(diff,2) lt diff) and (shift(diff,3) lt diff) and (shift(diff,4) lt diff) and (shift(diff,5) lt diff) $
                 and (shift(diff,-2044) lt diff) and (shift(diff,-2045) lt diff) and (shift(diff,-2046) lt diff) $
                 and (shift(diff,-2047) lt diff) and (shift(diff,-2048) lt diff) and (shift(diff,-2049) lt diff) $
                 and (shift(diff,-2050) lt diff) and (shift(diff,-2051) lt diff) and (shift(diff,-2052) lt diff) $
                 and (shift(diff,2044) lt diff) and (shift(diff,2045) lt diff) and (shift(diff,2046) lt diff) $
                 and (shift(diff,2047) lt diff) and (shift(diff,2048) lt diff) $
                 and (shift(diff,2049) lt diff) and (shift(diff,2050) lt diff) and (shift(diff,2051) lt diff) $
                 and (shift(diff,2052) lt diff) $
                 and (shift(diff,-4084) lt diff) and (shift(diff,-4085) lt diff) and (shift(diff,-4086) lt diff) $
                 and (shift(diff,-4087) lt diff) and (shift(diff,-4088) lt diff)  $
                 and (shift(diff,4084) lt diff) and (shift(diff,4085) lt diff) and (shift(diff,4086) lt diff) $
                 and (shift(diff,4087) lt diff) and (shift(diff,4088) lt diff) $
                 and (shift(diff,-6143) lt diff) and (shift(diff,-6144) lt diff) and (shift(diff,-6145) lt diff) $
                 and (shift(diff,6143) lt diff) and (shift(diff,6144) lt diff) and (shift(diff,6145) lt diff) $
                 and (shift(diff,-8191) lt diff) and (shift(diff,-8192) lt diff) and (shift(diff,-8193) lt diff) $
                 and (shift(diff,8191) lt diff) and (shift(diff,8192) lt diff)) and (shift(diff,8193) lt diff)  $
                 and (shift(diff,-10240) lt diff) $
                 and (shift(diff,10240) lt diff))$
                 ,length_5)

if length_5 ne 0 then begin

    radius=5
    star_indices_5=array_indices(diff,star_vector_5)
    print,'5 guys'
    print,star_indices_5
    print,length_5

    for k=double(0),(length_5-1) do begin

        printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_5(0,k),star_indices_5(1,k),radius,Index_string

    endfor

    for i=-5,5  do begin
       for j=-5,5 do begin

         diff(star_vector_5+i*Array_Dimension+j)=0

       endfor
    endfor

endif

star_vector_4=where(((shift(diff,-4) gt THRESHOLD) and (shift(diff,-3) gt THRESHOLD) and (shift(diff,-2) gt THRESHOLD) and (shift(diff,-1) gt THRESHOLD) and (diff gt THRESHOLD) and (shift(diff,1) gt THRESHOLD) and (shift(diff,2) gt THRESHOLD)  $
                 and (shift(diff,3) gt THRESHOLD) and (shift(diff,4) gt THRESHOLD) $
                 and (shift(diff,-2045) gt THRESHOLD) and (shift(diff,-2046) gt THRESHOLD) and (shift(diff,-2047) gt THRESHOLD) and (shift(diff,-2048) gt THRESHOLD) and (shift(diff,-2049) gt THRESHOLD) and (shift(diff,-2050) gt THRESHOLD) $
                 and (shift(diff,-2051) gt THRESHOLD) $
                 and (shift(diff,2045) gt THRESHOLD) and (shift(diff,2046) gt THRESHOLD) and (shift(diff,2047) gt THRESHOLD) and (shift(diff,2048) gt THRESHOLD) and (shift(diff,2049) gt THRESHOLD) $
                 and (shift(diff,2051) gt THRESHOLD) $
                 and (shift(diff,-4085) gt THRESHOLD) and (shift(diff,-4086) gt THRESHOLD) and (shift(diff,-4087) gt THRESHOLD)  $
                 and (shift(diff,4085) gt THRESHOLD) and (shift(diff,4086) gt THRESHOLD) and (shift(diff,4087) gt THRESHOLD) $
                 and (shift(diff,-6143) gt THRESHOLD) and (shift(diff,-6144) gt THRESHOLD) and (shift(diff,-6145) gt THRESHOLD) $
                 and (shift(diff,6143) gt THRESHOLD) and (shift(diff,6144) gt THRESHOLD) and (shift(diff,6145) gt THRESHOLD) $
                 and (shift(diff,-8192) gt THRESHOLD) $
                 and (shift(diff,8192) gt THRESHOLD)), $
                 length_4)

if length_4 ne 0 then begin

    radius=4
    star_indices_4=array_indices(diff,star_vector_4)
    print,'4 guys'
    print,length_4
    print,star_indices_4
    printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_4(0),star_indices_4(1),radius,Index_string

    for k=double(0),length_4 do begin

        printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_4(0,k),star_indices_4(1,k),radius,Index_string

    endfor

    for i=-4,4 do begin
       for j=-4,4 do begin

         diff(star_vector_4+i*Array_Dimension+j)=0

       endfor
    endfor

endif

star_vector_3=where(((shift(diff,-3) gt THRESHOLD) and (shift(diff,-2) gt THRESHOLD) and (shift(diff,-1) gt THRESHOLD) and (diff gt THRESHOLD) and (shift(diff,1) gt THRESHOLD) and (shift(diff,2) gt THRESHOLD) and (shift(diff,3) gt THRESHOLD) $
                 and (shift(diff,-2046) gt THRESHOLD) and (shift(diff,-2047) gt THRESHOLD) and (shift(diff,-2048) gt THRESHOLD) and (shift(diff,-2049) gt THRESHOLD) and (shift(diff,-2050) gt THRESHOLD) $
                 and (shift(diff,2046) gt THRESHOLD) and (shift(diff,2047) gt THRESHOLD) and (shift(diff,2048) gt THRESHOLD) and (shift(diff,2049) gt THRESHOLD) and (shift(diff,2050) gt THRESHOLD) $
                 and (shift(diff,-4085) gt THRESHOLD) and (shift(diff,-4086) gt THRESHOLD) and (shift(diff,-4087) gt THRESHOLD)  $
                 and (shift(diff,4085) gt THRESHOLD) and (shift(diff,4086) gt THRESHOLD) and (shift(diff,4087) gt THRESHOLD) $
                 and (shift(diff,6144) gt THRESHOLD) $
                 and (shift(diff,-6144) gt THRESHOLD)), $
                 length_3)

if length_3 ne 0 then begin

    radius=3
    star_indices_3=array_indices(diff,star_vector_3)
    print,'3 guys'
    print,length_3
    print,star_indices_3

    for k=double(0),(length_3-1) do begin

        printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_3(0,k),star_indices_3(1,k),radius,Index_string

    endfor

    for i=-3,3 do begin
       for j=-3,3 do begin

         diff(star_vector_3+i*Array_Dimension+j)=0

       endfor
    endfor
endif


;       -
;      ---
;     -----   Star
;      ---
;       -

star_vector_2=where(((shift(diff,-2) gt THRESHOLD) and (shift(diff,-1) gt THRESHOLD) and ((diff gt THRESHOLD) and (shift(diff,1) gt THRESHOLD) and (shift(diff,2) gt THRESHOLD) $
                 and (shift(diff,-2047) gt THRESHOLD) and (shift(diff,-2048) gt THRESHOLD) and (shift(diff,-2049) gt THRESHOLD)) $
                 and (shift(diff,2047) gt THRESHOLD) and (shift(diff,2048) gt THRESHOLD) and (shift(diff,2049) gt THRESHOLD) $
                 and (shift(diff,4086) gt THRESHOLD) $
                 and (shift(diff,-4086) gt THRESHOLD)), $
                 length_2)

if length_2 ne 0 then begin

    radius=2
    star_indices_2=array_indices(diff,star_vector_2)
    print,'2 guys'
    print,length_2
    ;print,star_indices_2

    for k=double(0),(length_2-1) do begin

        printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_2(0,k),star_indices_2(1,k),radius,Index_string

    endfor

    for i=-2,2 do begin
       for j=-2,2 do begin

         diff(star_vector_2+i*Array_Dimension+j)=0

       endfor
    endfor

endif

star_vector_1=where(((diff gt THRESHOLD) and (shift(diff,1) gt THRESHOLD) and (shift(diff,-1) gt THRESHOLD) $
                 and (shift(diff,-2048) gt THRESHOLD) and (shift(diff,2048) gt THRESHOLD)),length_1)

if length_1 ne 0 then begin

    radius=1
    star_indices_1=array_indices(diff,star_vector_1)
    print,'1 guys'
    print,length_1
    printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_1(0),star_indices_1(1),radius,Index_string

    for k=double(0),(length_1-1) do begin

        printf,1,FORMAT='(%"%d\t\t%d\t\t%d\t\t%s")',star_indices_1(0,k),star_indices_1(1,k),radius,Index_string

    endfor

    for i=-1,1 do begin
       for j=-1,1 do begin

         diff(star_vector_1+i*Array_Dimension+j)=0

       endfor
    endfor
endif


close,1

end


