;find_center
;
;PURPOSE:
;	find_center will read in the listed coordinates from the files 
;	of type [filename]_mag.txt and try to find the center of the 
; 	coordinate grid     

pro find_center, Directory=Directory, y_center=y_center, x_center=x_center, $
		 letter=letter,file_number=file_number

If (not keyword_set(Directory)) then Directory='/a/pippin01/Volumes/u08/lsst/CPanalysis/2005-05-10/ShackHartman/'
If (not keyword_set(letter)) then letter='a'
If (not keyword_set(file_number)) then file_number='20'
If (not keyword_set(x_center)) then x_center=[542.5]
If (not keyword_set(y_center)) then y_center=[466]

mag_file_name=Directory+'hr5132'+letter+'_'+file_number+'_stats.txt'
num=25

;The x coordinates are symmetric about a point that is not a ShackHartman image
;The y coordinates are symmetric about a ShackHartman image
subx=2*indgen(num)-25
suby=indgen(num)-12

index=findgen(num*num)

fullx=fltarr(num*num)
fully=fltarr(num*num)

;y is constant for 25 entries; x increases for 25 entries and repeats
;i.e., scan is along x followed by y
fullx(index)=subx(index mod num)
fully=rebin(suby,num*num)

readcol,mag_file_name,x_coord,y_coord,format='X,f,f'

point_indices=where(x_coord ne -1)

x_coord=x_coord(point_indices)
y_coord=y_coord(point_indices)
fullx=fullx(point_indices)
fully=fully(point_indices)

x_fit=curvefit(fullx,x_coord,weights,x_center,/noderivative,function_name='linear_function_x',yerror=xerror)
y_fit=curvefit(fully,y_coord,weights,y_center,/noderivative,function_name='linear_function_y',yerror=yerror)

;Errors at each point
y_err=y_coord-y_fit
x_err=x_coord-x_fit

set_plot,'ps'
                                                                                
device,filename=Directory+'hr5132'+letter+file_number+'error.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=1  
!P.POSITION=[0.10,0.10,0.9,0.9]
                                                                              
plot,fully,y_coord,psym=3,title="Reference fit for center of hr5132"$
	+letter+'_'+file_number+'_raw.fits',$
        xtitle="x index",ytitle="x coordinate",background=white,color=black

                                                                       
ploterr,fully,y_coord,y_err,psym=3

device,/close

set_plot,'X'
stop

end