;xy_offsets ; ;PURPOSE: ; This script will return the difference in coordinates between a ; ShackHartman image centroid and the corresponding point on the ; reference grid ; ; Since the center of the grid seems to be drifting pro cn_ind_depth_halfangle_variance_avg_sub,Date=Date,stats_directory=stats_directory,$ x_center_file=x_center_file,y_center_file=y_center_file,$ depth=depth,x_off=x_off,y_off=y_off,num_stats_files=num_stats_files If (not keyword_set(Date)) then Date='10' If (not keyword_set(stats_directory)) then stats_directory='stats' If (not keyword_set(x_center_file)) then x_center_file='x_centers.txt' If (not keyword_set(y_center_file)) then y_center_file='y_centers.txt' If (not keyword_set(depth)) then depth=1 If (not keyword_set(postscript_dir)) then postscript_dir='postscripts' If (not keyword_set(x_off)) then x_off=0 If (not keyword_set(y_off)) then y_off=0 if Date eq '10' then begin letter=['hr5132a','hr5132b','hr5132c','hr5132d','hr5132e','hr5132f',$ 'hr5132g','hr5132h','hr5132i'] endif else if Date eq '11' then begin letter=['h5132g','h5132h','h5132i','hr5132a','hr5132b','hr5132c','hr5132d',$ 'hr5132e','hr5132f','hr5132g','hr5132h'] endif else if Date eq '12' then begin letter=['hr5132a','hr5132b','hr5132c','hr5132d','hr5132e','hr5132f',$ 'hr5132g','hr5132h','hr5132i'] endif stats_directory='/a/pippin01/Volumes/u08/lsst/CPanalysis/2005-05-'+Date+$ '/ShackHartman/'+stats_directory+'/' postscript_dir='/a/pippin01/Volumes/u08/lsst/CPanalysis/2005-05-'+Date+$ '/ShackHartman/'+postscript_dir+'/' Cor_directory='/a/pippin01/Volumes/u08/lsst/CPanalysis/2005-05-'+Date+$ '/ShackHartman/correlations/' for j=0, N_elements(letter)-1 do begin ;ARRAYS OF FILES-Averages have already been subtracted offset_files=file_search(stats_directory+letter(j)+'*offsets.txt',$ count=num_off_files) ;If (not keyword_set(num_stats_files)) then num_stats_files=num_off_files num_stats_files=num_off_files ;CONSTANTS num=25 grid_spots=625 index=findgen(625) ;DATA CUBE CONTAINING ALL ANGLES AND MAGNITUDES phi_arr=fltarr(num_stats_files,grid_spots) magnitude_arr=fltarr(num_stats_files,grid_spots) ;FILL THE DATA CUBES for file_number=0, num_stats_files-depth do begin ;Take the magnitude and angles from the offset.txt files and fill the cube readcol,offset_files(file_number),magnitude,phi,format='X,X,X,X,X,f,f' magnitude_arr(file_number,*)=magnitude phi_arr(file_number,*)=phi endfor ;DATA CUBEs are now in place. Begin the covariance ;CORRELATION ARRAYS Cl_arr=fltarr(num_stats_files,(2*num-1)*(num),depth) Ct_arr=fltarr(num_stats_files,(2*num-1)*(num),depth) mags=fltarr((2*num-1)*(num)) angles=fltarr((2*num-1)*(num)) ;Loop through all of the files in the stack for file_number=0, num_stats_files-depth do begin print,file_number,offset_files(file_number) for time=0, depth-1 do begin spacing_index=0 for delta_x=-(num-1),(num-1) do begin for delta_y=0,(num-1) do begin ;THE ANGLE AND SEPARTATION FOR A GIVEN SET OF DELTAS grid_phi=atan(delta_y,delta_x) grid_sep=sqrt(delta_x^2+delta_y^2) ;Get the grid angles and spacings if file_number eq 0 then begin angles(spacing_index)=grid_phi mags(spacing_index)=grid_sep endif ;Cl=fltarr(grid_spots,depth) Ct=fltarr(grid_spots-abs(delta_x)*abs(delta_y)) Cl=fltarr(grid_spots-abs(delta_x)*abs(delta_y)) ; Cl=fltarr(grid_spots,depth) ;This index represents the number of points going into the sum for ;a given separation c_index=0 for sub_col_index=0,num-1 do begin for sub_row_index=0, num-1 do begin ;Check to see if a valid spot is located at the spot we're interrogating if (((sub_row_index+delta_x) lt num) and $ ((sub_row_index+delta_x) ge 0)) then begin first_spot=num*sub_col_index+sub_row_index second_spot=num*sub_col_index+sub_row_index+num*delta_y+delta_x if ((second_spot ge 0) and (second_spot lt 625) and $ (first_spot ge 0) and (first_spot lt 625)) then begin ;If it is on the grid, check to see if it has a partner if ((phi_arr(file_number,first_spot) ne -10) and $ (phi_arr(file_number,second_spot) ne -10))$ then begin ;Longitudinal will be cosine and transverse will be sine Cl(c_index,time)=magnitude_arr(file_number,first_spot)*$ cos(phi_arr(file_number,first_spot)-grid_phi)*$ magnitude_arr(file_number,second_spot)*$ cos(phi_arr(file_number,second_spot)-grid_phi) Ct(c_index,time)=magnitude_arr(file_number,first_spot)*$ sin(phi_arr(file_number,first_spot)-grid_phi)*$ magnitude_arr(file_number,second_spot)*$ sin(phi_arr(file_number,second_spot)-grid_phi) c_index=c_index+1 endif endif endif endfor endfor ; if grid_phi eq 0 then print,c_index,spacing_index ;Correlations for a given delta_x and delta_y are now in place correlated=where(Cl ne 0) ;If there were any pairs for a covariance, get there mean If correlated(0) ne -1 then begin Cl_arr(file_number,spacing_index,time)=$ mean(Cl(correlated)) Ct_arr(file_number,spacing_index,time)=$ mean(Ct(correlated)) endif else begin ;If there were no pairs, assign a system no value Cl_arr(file_number,spacing_index,time)=!values.f_nan Ct_arr(file_number,spacing_index,time)=!values.f_nan endelse ;Now onto the next set of spacings. Increment delta_y spacing_index=spacing_index+1 endfor endfor endfor endfor Cl_means=fltarr((2*num-1)*(num)) Ct_means=fltarr((2*num-1)*(num)) ;Get the average for each of these guys for n=0, ((2*num-1)*(num))-1 do begin finite_val=where(finite(cl_arr(*,n)) eq 1) if finite_val(0) ne -1 then begin Cl_means(n)=mean(cl_arr(where(finite(cl_arr(*,n)) eq 1),n)) Ct_means(n)=mean(ct_arr(where(finite(ct_arr(*,n)) eq 1),n)) endif else begin Cl_means(n)=!values.f_nan Ct_means(n)=!values.f_nan endelse endfor sort_index=sort(angles) ;Sort things by angle sorted_mags=mags(sort_index) sorted_angles=angles(sort_index) sorted_Cl_means=Cl_means(sort_index) sorted_Ct_means=Ct_means(sort_index) uniq_indices=uniq(sorted_angles) ;Print this to a file get_lun,cor_angles openw,cor_angles,stats_directory+'cor_angles'+letter(j)+'.txt' for i=0, N_elements(sort_index)-1 do begin printf,cor_angles,FORMAT='(%"%f\t%f\t%f\t%f\t%f")',sorted_mags(i),sorted_angles(i),$ sorted_Cl_means(i),sorted_Ct_means(i) endfor close,cor_angles free_lun,cor_angles set_plot,'ps' device,filename=Postscript_dir+letter(j)+'ind_both_angles_half_'+date+'.ps' loadct,39 ;load color table 39 device,/color ;allow color on the postscript device,ysize=8.5,/inches ;Height of plot in y device,xsize=10.0,/inches ;Width of plot in x device,yoffset=1.0,/inches ;Y position of lower left corner white='FFFFFF'x black='000000'x !P.CHARSIZE=.6 !P.THICK=4. color_arr=[0] angle_arr=['0'] sym_arr=[0] ;Set up plot plot,[0,34],[-1.5,1.5],/xstyle,/ystyle,/nodata,xr=[0,26],yr=[-0.5,2],background=white,$ position=[0.10,0.10,0.8,0.9],title='Average Longitudinal Correlations for 5/'+$ Date+'/2005',xtitle='Separation in n',ytitle='Covariance' color_index=0 for ind=0, N_elements(uniq_indices)-2 do begin if ind eq 0 then begin temp_mags=sorted_mags(0:uniq_indices(ind)) temp_angles=sorted_angles(0:uniq_indices(ind)) temp_cl_means=sorted_cl_means(0:uniq_indices(ind)) temp_sorted_indices=sort(temp_mags) oplot,temp_mags,temp_cl_means,color=10*color_index,psym=1 sym_arr=[sym_arr,1] color_arr=[color_arr,(10*color_index)] angle_arr=[angle_arr,string(strtrim(180*temp_angles(0)/!pi,2))] color_index=color_index+1 endif if ((uniq_indices(ind+1)-uniq_indices(ind)) gt 4) then begin temp_mags=sorted_mags(uniq_indices(ind)+1:uniq_indices(ind+1)) temp_angles=sorted_angles(uniq_indices(ind)+1:uniq_indices(ind+1)) temp_cl_means=sorted_cl_means(uniq_indices(ind)+1:uniq_indices(ind+1)) temp_sorted_indices=sort(temp_mags) print,temp_mags(temp_sorted_indices) print,temp_angles(temp_sorted_indices) print,temp_cl_means(temp_sorted_indices) oplot,temp_mags,temp_cl_means,color=10*color_index,psym=1 sym_arr=[sym_arr,1] color_arr=[color_arr,(10*color_index)] angle_arr=[angle_arr,string(strtrim(180*temp_angles(0)/!pi,2))] color_index=color_index+1 endif endfor legend,angle_arr,psym=sym_arr,colors=color_arr,position=[.82,.9],/normal !p.multi=[0,1,1] plot,[0,30],[-1.5,1.5],/xstyle,/ystyle,/nodata,xr=[0,26],yr=[-.5,2],background=white,$ position=[0.10,0.10,0.8,0.9],title='Average Transverse Correlations for 5/'+$ Date+'/2005',xtitle='Separation in n',ytitle='Covariance' color_arr=[0] angle_arr=['0'] sym_arr=[0] color_index=0 for ind=0, N_elements(uniq_indices)-2 do begin if ind eq 0 then begin temp_mags=sorted_mags(0:uniq_indices(ind)) temp_angles=sorted_angles(0:uniq_indices(ind)) temp_ct_means=sorted_ct_means(0:uniq_indices(ind)) temp_sorted_indices=sort(temp_mags) oplot,temp_mags,temp_ct_means,color=10*color_index,psym=1 sym_arr=[sym_arr,1] color_arr=[color_arr,(10*color_index)] angle_arr=[angle_arr,string(strtrim(180*temp_angles(0)/!pi,2))] color_index=color_index+1 endif if ((uniq_indices(ind+1)-uniq_indices(ind)) gt 4) then begin temp_mags=sorted_mags(uniq_indices(ind)+1:uniq_indices(ind+1)) temp_angles=sorted_angles(uniq_indices(ind)+1:uniq_indices(ind+1)) temp_ct_means=sorted_ct_means(uniq_indices(ind)+1:uniq_indices(ind+1)) temp_sorted_indices=sort(temp_mags) ;print,temp_mags(temp_sorted_indices) ;print,temp_angles(temp_sorted_indices) ;print,temp_ct_means(temp_sorted_indices) oplot,temp_mags,temp_ct_means,color=10*color_index,psym=1 sym_arr=[sym_arr,1] color_arr=[color_arr,(10*color_index)] angle_arr=[angle_arr,string(strtrim(180*temp_angles(0)/!pi,2))] color_index=color_index+1 endif endfor legend,angle_arr,psym=sym_arr,colors=color_arr,position=[.82,.9],/normal endfor device,/close set_plot,'X' stop end