;xy_offsets
;
;PURPOSE:
;	This script will calculate the covariances of the ShackHartman images 
;	using the distorted grid obtained by averaging all of the positions.

pro fft_offsets,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


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='newest_x_centers.txt'
If (not keyword_set(y_center_file)) then y_center_file='newest_y_centers.txt'
If (not keyword_set(depth)) then depth=1
If (not keyword_set(postscript_dir)) then postscript_dir='distorted_postscripts'
If (not keyword_set(x_off)) then x_off=0
If (not keyword_set(y_off)) then y_off=0

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+'/'

;ARRAYS OF FILES
offset_files=file_search(stats_directory+'*5132*newest_offsets.txt',count=num_stats_files)

;CONSTANTS
num=25
grid_spots=625
index=findgen(625)

;ARRAY REFERENCE 25x25=625 FOR GRID
subx=2*indgen(num)-25           ;The x-reference numbers
suby=indgen(num)-12             ;The y-reference numbers

total_images=num_stats_files

;HUGE DATA CUBE
x_offset_hor=fltarr(total_images,num,num)
y_offset_hor=fltarr(total_images,num,num)
x_offset_ver=fltarr(total_images,num,num)
y_offset_ver=fltarr(total_images,num,num)
xcoord_hor=fltarr(total_images,num,num)
ycoord_hor=fltarr(total_images,num,num)
xcoord_ver=fltarr(total_images,num,num)
ycoord_ver=fltarr(total_images,num,num)

;Correlations
Cxt_full=fltarr(total_images,num,depth)
Cxl_full=fltarr(total_images,num,depth)
Cyt_full=fltarr(total_images,num,depth)
Cyl_full=fltarr(total_images,num,depth)


;FILL THE DATA CUBES
for file_number=0, total_images-1 do begin

	;Take the results from the offset.txt files and fill the cube
	readcol,offset_files(file_number),x_coord,x_offset,y_coord,y_offset,$
	format='X,f,f,f,f'
	
	;Obtain the difference ordered in COLUMN, ROW between centroid and grid 	point
	zeros=where(x_coord eq -10)
        x_coord(zeros)=0
        y_coord(zeros)=0
	x_offset(zeros)=0
	y_offset(zeros)=0

	x_offset_hor(file_number,*,*)=x_offset
	y_offset_hor(file_number,*,*)=y_offset
	xcoord_hor(file_number,*,*)=x_coord
	ycoord_hor(file_number,*,*)=y_coord

	;Obtain the difference ordered in ROW, COLUMN
	x_offset_ver(file_number,*,*)=x_offset((index mod 25)*25+index/25)
	y_offset_ver(file_number,*,*)=y_offset((index mod 25)*25+index/25)
	xcoord_ver(file_number,*,*)=x_coord((index mod 25)*25+index/25)
        ycoord_ver(file_number,*,*)=y_coord((index mod 25)*25+index/25)
	
endfor

	
;Fourier Transform info
;
;x_array(i,j,k)
;i=time		j=column	k=row
;fft(array,dimension=dimension)
;
;DIMENSION=3 (across rows    (0,0,*))-----AVG(DIM=2) avg columns
;DIMENSION=2 (across columns (0,*,0))-----AVG(DIM=1) avg rows
;DIMENSION=1 (across depth   (*,0,0))-----AVG(DIM=0)
;
fxoffv=fft(x_offset_hor,dimension=3)
fyoffv=fft(y_offset_hor,dimension=3)

fxoffh=fft(x_offset_hor,dimension=2)
fyoffh=fft(y_offset_hor,dimension=2)

;Avg the transform coefficents for the rows
fxoffha1=avg(fxoffh,2)
fxoffha2=avg(fxoffha1,0)
fyoffha1=avg(fyoffh,2)
fyoffha2=avg(fyoffha1,0)

fyoffva1=avg(fyoffv,1)
fyoffva2=avg(fyoffva1,0)
fxoffva1=avg(fxoffv,1)
fxoffva2=avg(fxoffva1,0)

set_plot,'ps'
device,filename=Postscript_dir+'FFT_offsets.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=.7
!P.THICK=4.
                                                                                    
items_time=['t=0']
sym_time=[1]
items_k_tran=['Data Points','n!e-7/32!n','n!e-13/64!n','n!e-27/128!n']
sym_k_tran=[1,0,0,0]
line_k_tran=[0,1,2,3]
items_k_lon=['Data Points','n!e-1/3!n']
sym_k_long=[0,1]
colors_time=[0]

!p.multi=[0,1,2]

ytop=max([max(fxoffha2),max(fyoffva2)])

plot,abs(fxoffha2),title='Fourier Transform of Spot Deflections',$
	xtitle='K',ytitle='Magnitude',yrange=[0,ytop+0.002]
oplot,abs(fyoffva2)

ytop=max([max(fyoffha2),max(fxoffva2)])
plot,abs(fyoffha2),$
         xtitle='K',ytitle='Magnitude',yrange=[0,ytop+0.002]
oplot,abs(fxoffva2)

!p.multi=[0,1,2]
plot,atan(fxoffha2,/phase),title='Fourier Transform of Spot Deflections',$
         xtitle='K',ytitle='Phase',yrange=[-3.5,3.5]
oplot,atan(fyoffva2,/phase)

plot,atan(fxoffva2,/phase),$
	 xtitle='K',ytitle='Phase',yrange=[-3.5,3.5]

oplot,atan(fyoffha2,/phase)

device,/close
set_plot,'x'

stop

end


