pro center_drift, Directory=Directory,stats_directory=stats_directory,$
		  center_file=center_file,x_centers=x_centers,$
		  y_centers=y_centers

If (not keyword_set(Directory)) then Directory='/a/pippin01/Volumes/u08/lsst/CPanalysis/2005-05-10/ShackHartman/'
If (not keyword_set(stats_directory)) then stats_directory='stats_test'
If (not keyword_set(center_file)) then center_file=Directory+stats_directory+'/centers.txt'
If (not keyword_set(x_centers_file)) then x_centers_file=Directory+stats_directory+'/x_centers.txt'
If (not keyword_set(y_centers_file)) then y_centers_file=Directory+stats_directory+'/y_centers.txt'


;FILES AND HOLDING
readcol,center_file,x_centers,y_centers,format='f,f'
num_files=N_elements(x_centers)
get_lun,xfile
get_lun,yfile

x_slopes=fltarr(num_files)
y_slopes=fltarr(num_files)
x_start_index=0
x_stop_index=0
y_start_index=0
y_stop_index=0
x_fit_guess=fltarr(2)
y_fit_guess=fltarr(2)
file_index_arr=findgen(num_files)

openw,xfile,x_centers_file
openw,yfile,y_centers_file

printf,xfile,FORMAT='(%"%s\t%s\t\t%s\t%s\t%s")',$
      'X center','X Fitted Center','m','b','diff'
printf,yfile,FORMAT='(%"%s\t%s\t\t%s\t%s\t%s")',$
      'X center','X Fitted Center','m','b','diff'

set_plot,'ps'

device,filename=Directory+'postscripts/y_centers.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=1.2
!P.THICK=4.
;!P.NOERASE=0
;!P.POSITION=[0.10,0.10,0.9,0.9]
                                                                                

for file_index=0,num_files-2 do begin
	
	;Build array of x and y coordinates and corresponding slopes to feed 
	;to the linear fitting functions

	x_slopes(file_index)=x_centers(file_index+1)-x_centers(file_index)
	y_slopes(file_index)=y_centers(file_index+1)-y_centers(file_index)

	If ((y_slopes(file_index) gt 3) and (file_index-y_start_index gt 1))$
	   or (file_index eq num_files-2) or (((file_index+1) mod 99) eq 0) then begin
                y_fit_guess(0)=float(y_centers(y_start_index))
                y_fit_guess(1)=float(mean(y_slopes[y_start_index:file_index]))
                
		y_fit=curvefit(findgen(file_index-y_start_index+1),$
                y_centers[y_start_index:file_index],weights,y_fit_guess,$
                /noderivative,function_name='linear_function_2_y',$
                yerror=yerror)
			
		plot,y_centers[y_start_index:file_index],$
       		 psym=1,title="Fitted Y coordinates for "+string(y_start_index)+' through '+$
		 strtrim(string(file_index)),xtitle="File Sub Number",ytitle="Y Coordinate",$
		 background=white,color=black,yrange=[460,472],POSITION=[0.10,0.10,0.9,0.9]                                                                           
		oplot,y_fit

		for y_sub_index=0,n_elements(y_fit)-1 do begin   
			y_t_center=y_centers(y_start_index+y_sub_index)
			y_f_center=y_fit(y_sub_index)
			y_diff=y_f_center-y_t_center                                                             
                	printf,yfile,FORMAT='(%"%d\t%f\t%f\t%f\t%f\t%f")',$
			y_start_index+y_sub_index,y_t_center,$
			y_f_center,y_fit_guess(0),$
			y_fit_guess(1),y_diff

		endfor
			
			y_start_index=file_index+1
                	
		                                                                
        endif


	If ((x_slopes(file_index) gt 1.5) and (file_index-x_start_index gt 2)) $
		or (file_index eq num_files-2) then begin
		x_fit_guess(0)=x_centers(x_start_index)
		x_fit_guess(1)=float(mean(x_slopes[x_start_index:file_index]))
		
		x_fit=curvefit(findgen(file_index-x_start_index+1),$
		x_centers[x_start_index:file_index],weights,x_fit_guess,$
		/noderivative,function_name='linear_function_2_x',$
		yerror=xerror)

		for x_sub_index=0,n_elements(x_fit)-1 do begin
                        x_t_center=x_centers(x_start_index+x_sub_index)
                        x_f_center=x_fit(x_sub_index)
                        x_diff=x_f_center-x_t_center                                                                                 
                        printf,xfile,FORMAT='(%"%f\t%f\t%f\t%f\t%f")',$
                        x_t_center,x_f_center,x_fit_guess(0),$
                        x_fit_guess(1),x_diff
                                                                                
                endfor

		x_start_index=file_index+1
	endif


	print,file_index                                                                                


endfor	

close,xfile
close,yfile


free_lun,xfile
free_lun,yfile

device,/close

stop
end