pro fit_ang_power_spectrum_meters,$
	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,exam_date=exam_date


If (not keyword_set(Date)) then Date='10'
If (not keyword_set(exam_date)) then exam_date='1_31'
If (not keyword_set(stats_directory)) then stats_directory='stats_'+exam_date+'/'
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=0
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='/nfs/slac/g/ki/ki08/lsst/CPanalysis/2005-05-'+Date+$
          '/ShackHartman/'+stats_directory+'/'
postscript_dir='/nfs/slac/g/ki/ki08/lsst/CPanalysis/'+exam_date+'_figs/'

common func_power_spec, kvals, xpower_fit, ypower_fit

nonfiducial=[1,2,3,4,5,6,7,8,9,13,14,17,18,19,20,21,22,23,24,25,$
          26,27,28,29,30,31,32,38,39,45,46,47,48,49,50,$
          51,52,53,54,55,63,64,71,72,73,74,75,$
          76,77,78,79,88,89,98,98,99,100,$
          101,102,103,113,114,123,124,125,$
          126,127,128,138,139,149,150,$
          151,152,163,164,174,175,$
          176,177,188,189,200,$
          201,210,213,214,225,$
          226,238,239,249,$
          251,261,262,263,264,265,266,$
          276,286,287,288,289,290,291,$
          301,302,303,304,305,306,307,308,309,$
          310,311,312,313,314,315,316,$
          317,318,319,320,321,322,323,324,325,$
          326,327,336,337,338,339,340,341,$
          351,361,362,363,364,365,366,$
          376,388,389,$
          401,413,414,425,$
          426,427,438,439,449,450,$
          451,452,454,463,464,473,474,475,$
          476,477,478,488,489,498,499,500,$
          501,502,503,504,513,514,522,523,524,525,$
          526,527,528,529,530,538,539,546,547,548,549,550,$
          551,552,553,554,555,556,563,564,570,571,572,573,574,575,$
          576,577,578,579,580,581,582,588,589,594,595,596,597,598,599,600,$
          601,602,603,604,605,606,607,608,609,610,613,614,617,618,619,620,621,$
          622,623,624,625]-1

;Extra Index for the scale offset
complement,findgen(626),nonfiducial,fiducial
fitted=fiducial
Ang_factor=float(((24e-6/60e-3)*(3.5e-3)))

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

;CONSTANTS
kmax=12
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

;HUGE DATA CUBE
x_offset_hor=complexarr(num_stats_files,num,num)
y_offset_hor=complexarr(num_stats_files,num,num)
x_offset_ver=complexarr(num_stats_files,num,num)
y_offset_ver=complexarr(num_stats_files,num,num)

;ffts
xfft_array=complexarr(num_stats_files,num,num)
yfft_array=complexarr(num_stats_files,num,num)

xpower=fltarr(num_stats_files,num,num)
ypower=fltarr(num_stats_files,num,num)

;Grid for partvelvec
x_ref_col=12.5*subx+542
y_ref_col=25.0*suby+460
index=findgen(num*num)
gridx=fltarr(num*num)
gridy=fltarr(num*num)
gridx(index)=x_ref_col(index mod num)
gridy(index)=y_ref_col(index/num)

;Indices and such
no_bins=8
a=findgen(13)+12
b=findgen(12)
subn=[a,b]
subm=subn

readcol,offset_files(0),xcoord,ycoord,format='x,f,x,f,x'

;FILL THE DATA CUBES

for file_number=0,num_stats_files-depth-1 do begin

        ;Take the results from the offset.txt files and fill the cube
        readcol,offset_files(file_number),x_offset,y_offset,format='X,X,f,X,f'

        ;Obtain the difference ordered in COLUMN, ROW between centroid and grid         point
	
	x_offset(where(x_offset eq -10))=0
        y_offset(where(y_offset eq -10))=0

        x_offset_hor(file_number,*,*)=Ang_factor*x_offset
        y_offset_hor(file_number,*,*)=Ang_factor*y_offset

        ;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)

 	xfft_array(file_number,*,*)=fft(x_offset_hor(file_number,*,*))
	yfft_array(file_number,*,*)=fft(y_offset_hor(file_number,*,*))
	
	xpower(file_number,*,*)=(abs(xfft_array(file_number,*,*))^2)
 	ypower(file_number,*,*)=(abs(yfft_array(file_number,*,*))^2)
        
endfor

set_plot,'ps'
device,filename=Postscript_dir+'Power_Spectrum_Radians_m_plot_5_'+Date+'.ps'                                                                                
loadct,0                      ;load color table 39
device,/color                   ;allow color on the postscript
device,ysize=8.5,/inches        ;Height of plot in y
device,xsize=11.0,/inches        ;Width of plot in x
device,yoffset=1.0,/inches      ;Y position of lower left corner

white='FFFFFF'x
black='000000'x   


xpower_tot=avg(xpower,0)
ypower_tot=avg(ypower,0)

xpower_re=complexarr(num,num)
ypower_re=complexarr(num,num)

        ;Rearrange indices
        for n=0,num-1 do begin
            for m=0,num-1 do begin
                xpower_re(subn(n),subm(m))=xpower_tot(n,m)
                ypower_re(subn(n),subm(m))=ypower_tot(n,m)
            endfor
        endfor

xpower_fit=double(xpower_re(13:24,12))
ypower_fit=double(ypower_re(12,13:24))
kvals=2*!pi*(findgen(12)+1)*1./4.2

rs=findgen(2)
tries=2000

if date eq '11' then begin
rs=amoeba(1.0e-8,scale=[4e-14,.15],p0=[2e-14,-5./3],function_value=fval,$
   function_name='fit_power_spectrum_func',ncalls=tries,$
   nmax=10000)
endif else begin
rs=amoeba(1.0e-18,scale=[1e-13,.05],p0=[4e-14,-5./3],function_value=fval,$
   function_name='fit_power_spectrum_func',ncalls=tries,$
   nmax=1000000)
endelse

items_all=['Data','A*k^-B']
sym_all=[2,3]
plot,xpower_re(12:24,12),xtitle='k',ytitle='Power',$
     title='X Offset Power spectrum for night of '+Date+$
	   ' for kx=0-12 ky=0',psym=2,yrange=[0,max(xpower_re(13:24,12))]

oplot,kvals,rs(0)*kvals^rs(1)
xyouts,.8,.5,'A='+strtrim(rs(0),2),/normal
xyouts,.8,.4,'B='+strtrim(rs(1),2),/normal

legend,items_all,psym=sym_all,position=[.7,.8],/normal

items_all=['Data','.1*k^-11/3']
plot,xpower_re(12,12:24),xtitle='k',ytitle='Power',$
     title='X Offset Power spectrum for night of '+Date+$
           ' for kx=0-12 ky=0',psym=2,yrange=[0,max(xpower_re(12,13:24))]
oplot,.1*(findgen(12)+1)^(-11./3)
legend,items_all,psym=sym_all,position=[.7,.8],/normal

items_all=['Data','A*k^-B']
plot,ypower_re(12,12:24),xtitle='k',ytitle='Power',$
     title='Y Offset Power spectrum for night of '+Date+$
           ' for kx=0, ky=0-12',psym=2,yrange=[0,max(ypower_re(12,13:24))]

oplot,kvals,rs(0)*kvals^rs(1)
xyouts,.8,.5,'A='+strtrim(rs(0),2),/normal
xyouts,.8,.4,'B='+strtrim(rs(1),2),/normal

legend,items_all,psym=sym_all,position=[.7,.8],/normal

items_all=['Data','.1*k^-11/3']
plot,ypower_re(12:24,12),xtitle='k',ytitle='Power',$
     title='Y Offset Power spectrum for night of '+Date+$
           ' for kx=0-12 ky=0',psym=2,yrange=[0,max(ypower_re(13:24,12))]
oplot,.1*(findgen(12)+1)^(-11./3)
legend,items_all,psym=sym_all,position=[.7,.8],/normal

device,/close
erase
 
  set_plot,'z'
  erase
  device, set_font='Courier'
  device,set_resolution=[800,600]
  !p.charsize=.8
  !p.charthick=1.2
  !x.thick=2
  !y.thick=2
  !p.thick=.8
  !p.noerase=0

  items_all=['Data','A*k^-B']
  plot,kvals,ypower_fit,xtitle='k',ytitle='Power',$
     title='Offset Power spectrum for night of '+Date,$
     psym=2,yrange=[0,max(xpower_fit)],$
     color=black,background=white
  oplot,kvals,xpower_fit,psym=1,color=black 
  oplot,kvals,rs(0)*kvals^rs(1),color=black
  xyouts,.8,.5,'A='+strtrim(rs(0),2),/normal,color=black
  xyouts,.8,.4,'B='+strtrim(rs(1),2),/normal,color=black

  jpgimg=tvrd()
  tvlct,reds,greens,blues,/get
  write_png,postscript_dir+'Power_Spectrum_rads_m_5-'+Date+'-2005.png',$
  jpgimg,reds,greens,blues


device,/close
set_plot,'x'

stop

end