;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 mpeg_phase_screen,Date=Date,stats_directory=stats_directory,$ x_center_file=x_center_file,y_center_file=y_center_file,$ low=low,lcut=lcut,phase_dim=phase_dim,$ xwindspeed=xwindspeed,ywindspeed=ywindspeed If (not keyword_set(Date)) then Date='10' If (not keyword_set(stats_directory)) then stats_directory='stats/' If (not keyword_set(letter)) then letter=['a','b','c','d','e','f','g','h','i'] If (not keyword_set(xwindspeed)) then xwindspeed=25 If (not keyword_set(ywindspeed)) then ywindspeed=25 If (not keyword_set(phase_dim)) then phase_dim=long(1024) 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(low)) then low=100.5 If (not keyword_set(lcut)) then lcut=10000.0 ;WINDSPEED THAT PHASE SCREEN IS DRAGGED xwindspeed_m=xwindspeed*.17 ywindspeed_m=ywindspeed*.17 windspeeds=[xwindspeed,ywindspeed] max_windspeed=max(windspeeds) ;DIMENSIONS OF PHASE SCREEN AND SUB_SECTION OF SCREEN sub_sec=512 offset=sub_sec/2-12 grid_spots=long(625) num=25 num_of_cor=(phase_dim-sub_sec)/max_windspeed ;INNER AND OUTER SCALES lowstr=strtrim(string(low,format='(%"%10.1f")'),2) lcutstr=strtrim(string(lcut,format='(%"%10.1f")'),2) lowinmeters=phase_dim*.17/low highinmeters=phase_dim*.17/lcut lowinmeters=strtrim(string(lowinmeters,format='(%"%10.2f")'),2) highinmeters=strtrim(string(highinmeters,format='(%"%10.4f")'),2) ;FIDUCIAL 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 ;DIRECTORIES Shack_Directory='/nfs/slac/g/ki/ki08/lsst/CPanalysis/2005-05-'+Date+$ '/ShackHartman/' Stats_Directory=Shack_Directory+stats_directory screen_directory='/nfs/slac/g/ki/ki08/lsst/CPanalysis/lsst_sim/atmosphere/' phase_files=file_search(screen_directory+'phase_'+strtrim(lowstr,2)+'_'+$ strtrim(lcutstr,2)+'*.out',$ count=num_phase_files) total_num_of_cor=num_of_cor*num_phase_files ;FILES stats_files=file_search(stats_directory+'*stats.txt',count=num_files) offset_files=file_search(Stats_directory+'*5132*newest_offsets.txt',count=num_o_files) ;GRID readcol,offset_files(0),grid_index,x_coord,x_off,y_coord,y_off,$ format='i,f,f,f,f' x_coord=[x_coord,0] y_coord=[y_coord,0] x_off=[x_off,4] y_off=[y_off,4] ;MAKE SUB ARRAY FOR SLOPES ACROSS SUBAPERTURE x_slopes_sub_arr=fltarr(total_num_of_cor,25,25) y_slopes_sub_arr=fltarr(total_num_of_cor,25,25) magnitude=fltarr(total_num_of_cor,25,25) phi=fltarr(total_num_of_cor,25,25) ;FFTs xfft_array=complexarr(total_num_of_cor,num,num) yfft_array=complexarr(total_num_of_cor,num,num) fiducial=1 for screen_num=0,num_phase_files-1 do begin ;Open the phase screen binary file and read floating point (data_type=4) openr,phase,phase_files(screen_num),/get_lun screen=read_binary(phase,data_dims=[phase_dim,phase_dim],data_type=4) close,phase free_lun,phase ;Use Central Differencing to get slopes and normalize to 1 xslopes=fltarr(phase_dim,phase_dim) yslopes=fltarr(phase_dim,phase_dim) yslopes(1:1022,1:1022)=(screen(2:1023,1:1022)-screen(0:1021,1:1022))/2 xslopes(1:1022,1:1022)=(screen(1:1022,2:1023)-screen(1:1022,0:1021))/2 yslopes=yslopes/sqrt(mean(yslopes^2)) xslopes=xslopes/sqrt(mean(xslopes^2)) ;WINDSPEED TO MOVE PHASE SCREEN ACROSS SUBAPERTURE ;FILL THE ARRAYS FOR LATER PROCESSING for increment=0,num_of_cor-1 do begin tot_increment=num_of_cor*screen_num+increment ;CORRELATIONS x_segment=increment*xwindspeed y_segment=increment*ywindspeed ;SLOPES x_slopes_sub_arr(tot_increment,*,*)=xslopes[$ phase_dim-offset-x_segment:phase_dim-x_segment-1-offset+num,$ phase_dim-offset-y_segment:phase_dim-y_segment-1-offset+num] y_slopes_sub_arr(tot_increment,*,*)=yslopes[$ phase_dim-offset-x_segment:phase_dim-x_segment-1-offset+num,$ phase_dim-offset-y_segment:phase_dim-y_segment-1-offset+num] ;25x25 grid, make one-dimensional and apply fiducial if fiducial eq 1 then begin xtemp=fltarr(num*num) ytemp=fltarr(num*num) xtemp=x_slopes_sub_arr(tot_increment,*,*) ytemp=y_slopes_sub_arr(tot_increment,*,*) xtemp(nonfiducial)=0 ytemp(nonfiducial)=0 x_slopes_sub_arr(tot_increment,*,*)=xtemp y_slopes_sub_arr(tot_increment,*,*)=ytemp endif ;GET FFT OF SLOPES xfft_array(tot_increment,*,*)=fft(x_slopes_sub_arr(tot_increment,*,*)) yfft_array(tot_increment,*,*)=fft(y_slopes_sub_arr(tot_increment,*,*)) endfor endfor readcol,Stats_Directory+x_center_file,x_fit_center,format='X,X,f' readcol,Stats_Directory+y_center_file,y_fit_center,format='X,X,f' ;Time will be coded as 2005*1000+Day_Number+hr/24+min/60+second/3600 year_offset=double(2005000) time_gap=0.0001 num_of_telescope_shifts=0 ;MOVIE mpegname='/nfs/slac/g/ki/ki08/lsst/CPanalysis/mpegs/Offsets_outer_'+$ strtrim(lowinmeters,2)+'_inner_'+strtrim(highinmeters,2)+'_stats.mpeg' movie=mpeg_open([1000,768]) frame_no=0 frame_gap=25 ;CONSTANTS num=25 ;ARRAY REFERENCE 25x25=625 FOR GRID subx=2*indgen(num)-25 ;The x-reference numbers suby=indgen(num)-12 ;The y-reference numbers ;GET THE OFFSETS FROM THE PHASE SCREENS ;SET PLOT TO MEMORY BUFFER IN ORDER TO BE ABLE TO DO OTHER THINGS set_plot,'Z' device,set_resolution=[1000,768] for file_number=0, 20 do begin ;num_of_files-1 do begin x_offset=fltarr(num*num) y_offset=fltarr(num*num) x_offset(*)=x_slopes_sub_arr(file_number,*,*) y_offset(*)=y_slopes_sub_arr(file_number,*,*) x_offset=[x_offset,0] y_offset=[y_offset,0] erase !p.noerase=1 !p.charsize=1 partvelvec,x_offset(fitted),y_offset(fitted),x_coord(fitted),$ y_coord(fitted),yrange=[100,820],xrange=[180,880],xstyle=17,$ ystyle=17,position=[.05,.05,.7,.9] vec_field=tvrd() for i=0,frame_gap do begin mpeg_put,movie,IMAGE=vec_field,frame=(frame_no*frame_gap+i),$ /order endfor frame_no=frame_no+1 endfor mpeg_save,movie,FILENAME=mpegname mpeg_close,movie end