source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_Marquee_Operations.ipf @ 1141

Last change on this file since 1141 was 1141, checked in by srkline, 3 years ago

updates to the beam center corrections that take into account the lateral variation of zero points of the tubes. This is an important correction to the reference beam center, especially when the graphite monochrmoator is used.

various bug fixes included as part of reduction , largely to catch missing or incorrect information in the file header.

New "patch" functions have been added for number of guides, beam stops, source aperture, etc. which were missing in the file and caused the resolution calculation to be totally incorrect.

File size: 14.1 KB
Line 
1#pragma TextEncoding = "MacRoman"
2#pragma rtGlobals=3             // Use modern global access method and strict wave access.
3
4
5// TODO:
6// x- get BoxSum Operational, for rudimentary Transmission calculations and ABS scale
7// x- need to be able to identify the detector panel from the Marquee selection
8// x- then do the count in (unscaled) coordinates of the actual data array
9//
10//
11// -- still many operations in (SANS)Marquee.ipf that will be useful to implement:
12//
13// -- writing the box coordinates and the counts (and error) to the data file
14// x- determining the beam center (centroid) of the selection
15//  -- writing beam center (centroid) to the file?
16//  -- a box sum over a range of files (with a plot)
17// -- box sum over annular regions
18// -- box sum over arcs
19// --  saving Graphics image
20// -- histogram of the counts
21//
22//
23//
24
25
26
27
28//
29//function will print marquee coordinates in axis terms, not detector terms
30//
31// will also calculate the sum in the box, and the error, and print it out
32//
33Function V_PrintMarqueeCoords() :  GraphMarquee
34        GetMarquee left,bottom
35        if(V_flag == 0)
36                Print "There is no Marquee"
37        else
38                Variable count,x1,x2,y1,y2,ct_err
39                x1 = V_left
40                x2 = V_right
41                y1 = V_bottom
42                y2 = V_top
43                printf "marquee left in bottom axis terms: %g\r",round(V_left)
44                printf "marquee right in bottom axis terms: %g\r",round(V_right)
45                printf "marquee bottom in left axis terms: %g\r",round(V_bottom)
46                printf "marquee top in left axis terms: %g\r",round(V_top)
47//              printf "**note that you must add 1 to each axis coordinate to get detector coordinates\r"
48               
49                // NOTE:
50                // this function MODIFIES x and y values on return, converting them to panel coordinates
51                // detector panel is identified from the (left,top) coordinate (x1,y2)
52                String detStr = V_FindDetStrFromLoc(x1,x2,y1,y2)               
53//              Printf "Detector = %s\r",detStr
54
55//
56                SVAR gCurDispType = root:Packages:NIST:VSANS:Globals:gCurDispType       
57                // this function will modify the x and y values (passed by reference) as needed to keep on the panel
58                V_KeepSelectionInBounds(x1,x2,y1,y2,detStr,gCurDispType)
59                Printf "%d;%d;%d;%d;\r",x1,x2,y1,y2
60
61
62                count = V_SumCountsInBox(x1,x2,y1,y2,ct_err,gCurDispType,detStr)
63
64               
65                Print "total counts = ",count
66                Print "err/counts = ",ct_err/count
67                print "average counts per pixel = ",count/(x2-x1)/(y2-y1)
68        endif
69End
70
71
72// NOTE:
73// this function MODIFIES x and y values on return, converting them to panel coordinates
74// detector panel is identified from the (left,top) coordinate (x1,y2)
75Function/S V_FindDetStrFromLoc(x1,x2,y1,y2)
76        Variable &x1,&x2,&y1,&y2
77       
78        // which images are here?
79        String detStr="",imStr,carriageStr
80        String currentImageRef,activeSubwindow,imageList
81        Variable ii,nIm,testX,testY,tab
82       
83        // which tab is selected? -this is the main graph panel (subwindow may not be the active one!)
84        ControlInfo/W=VSANS_Data tab0
85        tab = V_Value
86        if(tab == 0)
87                activeSubwindow = "VSANS_Data#det_panelsF"
88        elseif (tab == 1)
89                activeSubwindow = "VSANS_Data#det_panelsM"
90        else
91                activeSubwindow = "VSANS_Data#det_panelsB"
92        endif
93       
94        imageList = ImageNameList(activeSubwindow,";") 
95
96        nIm = ItemsInList(imageList,";")
97        if(nIm==0)
98                return("")              //problem, get out
99        endif
100
101        // images were added in the order TBLR, so look back through in the order RLBT, checking each to see if
102        // the xy value is found on that (scaled) array
103                               
104        // loop backwards through the list of panels (may only be one if on the back)
105        for(ii=nIm-1;ii>=0;ii-=1)
106                Wave w = ImageNameToWaveRef(activeSubwindow,StringFromList(ii, imageList,";"))
107               
108                // which, if any image is the mouse xy location on?
109                // use a multidemensional equivalent to x2pnt: (ScaledDimPos - DimOffset(waveName, dim))/DimDelta(waveName,dim)
110
111               
112                testX = ScaleToIndex(w,x1,0)
113                testY = ScaleToIndex(w,y2,1)
114               
115                if( (testX >= 0 && testX < DimSize(w,0)) && (testY >= 0 && testY < DimSize(w,1)) )
116                        // we're in-bounds on this wave
117                       
118                        // deduce the detector panel
119                        currentImageRef = StringFromList(ii, imageList,";")     //the image instance ##
120                        // string is "data", or "data#2" etc. - so this returns "", "1", "2", or "3"
121                        imStr = StringFromList(1, currentImageRef,"#")         
122                        carriageStr = activeSubWindow[strlen(activeSubWindow)-1]
123                       
124                        if(cmpstr(carriageStr,"B")==0)
125                                detStr = carriageStr
126                        else
127                                if(strlen(imStr)==0)
128                                        imStr = "9"                     // a dummy value so I can replace it later
129                                endif
130                                detStr = carriageStr+imStr              // "F2" or something similar
131                                detStr = ReplaceString("9", detStr, "T")        // ASSUMPTION :::: instances 0123 correspond to TBLR
132                                detStr = ReplaceString("1", detStr, "B")        // ASSUMPTION :::: this is the order that the panels
133                                detStr = ReplaceString("2", detStr, "L")        // ASSUMPTION :::: are ALWAYS added to the graph
134                                detStr = ReplaceString("3", detStr, "R")        // ASSUMPTION ::::
135                        endif
136                       
137                        Printf "Detector panel %s=(%d,%d)\r",detStr,testX,testY
138                       
139                        x1 = ScaleToIndex(w,x1,0)               // get all four marquee values to pass back to the calling function
140                        x2 = ScaleToIndex(w,x2,0)               // converted into detector coordinates
141                        y1 = ScaleToIndex(w,y1,1)
142                        y2 = ScaleToIndex(w,y2,1)
143                       
144                       
145                        ii = -1         //look no further, set ii to bad value to exit the for loop
146
147                endif
148               
149        endfor
150       
151        return(detStr)
152       
153End
154
155//testing function only, not called anywhere
156Function V_testKeepInBounds(x1,x2,y1,y2,detStr,folderStr)
157        Variable x1,x2,y1,y2
158        String detStr,folderStr
159       
160        V_KeepSelectionInBounds(x1,x2,y1,y2,detStr,folderStr)
161        Print x1,x2,y1,y2
162        return(0)
163End
164
165// for the given detector
166Function V_KeepSelectionInBounds(x1,x2,y1,y2,detStr,folderStr)
167        Variable &x1,&x2,&y1,&y2
168        String detStr,folderStr
169       
170        Variable pixelsX = V_getDet_pixel_num_x(folderStr,detStr)
171        Variable pixelsY = V_getDet_pixel_num_y(folderStr,detStr)
172
173        //keep selection in-bounds
174        x1 = (round(x1) >= 0) ? round(x1) : 0
175        x2 = (round(x2) <= (pixelsX-1)) ? round(x2) : (pixelsX-1)
176        y1 = (round(y1) >= 0) ? round(y1) : 0
177        y2 = (round(y2) <= (pixelsY-1)) ? round(y2) : (pixelsY-1)
178       
179        return(0)
180End
181
182
183
184//sums the data counts in the box specified by (x1,y1) to (x2,y2)
185//assuming that x1<x2, and y1<y2
186//the x,y values must also be in array coordinates[0] NOT scaled detector coords.
187//
188// accepts arbitrary detector coordinates. calling function is responsible for
189// keeping selection in bounds
190//
191Function V_SumCountsInBox(x1,x2,y1,y2,ct_err,type,detStr)
192        Variable x1,x2,y1,y2,&ct_err
193        String type,detStr
194       
195        Variable counts = 0,ii,jj,err2_sum
196       
197// get the waves of the data and the data_err
198        Wave w = V_getDetectorDataW(type,detStr)
199        Wave data_err = V_getDetectorDataErrW(type,detStr)
200
201       
202        err2_sum = 0            // running total of the squared error
203        ii=x1
204        jj=y1
205        do
206                do
207                        counts += w[ii][jj]
208                        err2_sum += data_err[ii][jj]*data_err[ii][jj]
209                        jj+=1
210                while(jj<=y2)
211                jj=y1
212                ii+=1
213        while(ii<=x2)
214       
215        err2_sum = sqrt(err2_sum)
216        ct_err = err2_sum
217       
218//      Print "error = ",ct_err
219//      Print "error/counts = ",ct_err/counts
220       
221        Return (counts)
222End
223
224
225//sums the data counts in the box specified by (x1,y1) to (x2,y2)
226//assuming that x1<x2, and y1<y2
227//the x,y values must also be in array coordinates[0] NOT scaled detector coords.
228//
229// accepts arbitrary detector coordinates. calling function is responsible for
230// keeping selection in bounds
231//
232// basically the same as V_SumCountsInBox, except the PBR value has been removed so that the
233// function can be used from the command line
234//
235Function V_SumCountsInBox_Cmd(x1,x2,y1,y2,type,detStr)
236        Variable x1,x2,y1,y2
237        String type,detStr
238       
239        Variable counts = 0,ii,jj,err2_sum,ct_err
240       
241// get the waves of the data and the data_err
242        Wave w = V_getDetectorDataW(type,detStr)
243        Wave data_err = V_getDetectorDataErrW(type,detStr)
244
245                       
246        err2_sum = 0            // running total of the squared error
247        ii=x1
248        jj=y1
249        do
250                do
251                        counts += w[ii][jj]
252                        err2_sum += data_err[ii][jj]*data_err[ii][jj]
253                        jj+=1
254                while(jj<=y2)
255                jj=y1
256                ii+=1
257        while(ii<=x2)
258       
259        err2_sum = sqrt(err2_sum)
260        ct_err = err2_sum
261
262        Print "sum of counts = ",counts
263        Print "error = ",ct_err
264        Print "error/counts = ",ct_err/counts
265       
266        Return (counts)
267End
268
269Proc pV_SumCountsInBox_Cmd(x1,x2,y1,y2,type,detStr)
270        Variable x1=280,x2=430,y1=350,y2=1020
271        String type="RAW",detStr="B"
272       
273        V_SumCountsInBox_Cmd(x1,x2,y1,y2,type,detStr)
274End
275
276Function V_FindCentroid() :  GraphMarquee
277
278//      //get the current displayed data (so the correct folder is used)
279//      SVAR cur_folder=root:myGlobals:gDataDisplayType
280//      String dest = "root:Packages:NIST:" + cur_folder
281       
282        Variable xzsum,yzsum,zsum,xctr,yctr
283        Variable left,right,bottom,top,ii,jj,counts
284        Variable x_mm_sum,y_mm_sum,x_mm,y_mm
285
286       
287        GetMarquee left,bottom
288        if(V_flag == 0)
289                Print "There is no Marquee"
290        else
291                left = round(V_left)            //get integer values for selection limits
292                right = round(V_right)
293                top = round(V_top)
294                bottom = round(V_bottom)
295
296                // detector panel is identified from the (left,top) coordinate (x1,y2)
297                String detStr = V_FindDetStrFromLoc(left,right,bottom,top)             
298        //      Printf "Detector = %s\r",detStr
299       
300        //
301                SVAR gCurDispType = root:Packages:NIST:VSANS:Globals:gCurDispType       
302                // this function will modify the x and y values (passed by reference) as needed to keep on the panel
303                V_KeepSelectionInBounds(left,right,bottom,top,detStr,gCurDispType)
304                Print left,right,bottom,top
305                       
306               
307                // selection valid now, calculate beamcenter
308                // get the waves of the data and the data_err
309                Wave data = V_getDetectorDataW(gCurDispType,detStr)
310                Wave data_err = V_getDetectorDataErrW(gCurDispType,detStr)
311               
312                // get the real-space information
313                String destPath = "root:Packages:NIST:VSANS:"+gCurDispType
314                Wave data_realDistX = $(destPath + ":entry:instrument:detector_"+detStr+":data_realDistX")
315                Wave data_realDistY = $(destPath + ":entry:instrument:detector_"+detStr+":data_realDistY")
316       
317                xzsum = 0
318                yzsum = 0
319                zsum = 0
320                x_mm_sum = 0
321                y_mm_sum = 0
322               
323                // count over rectangular selection, doing each row, L-R, bottom to top
324                ii = bottom -1
325                do
326                        ii +=1
327                        jj = left-1
328                        do
329                                jj += 1
330                                counts = data[jj][ii]
331                                xzsum += jj*counts
332                                yzsum += ii*counts
333                                zsum += counts
334                               
335                                x_mm_sum += data_realDistX[jj][ii]*counts
336                                y_mm_sum += data_realDistY[jj][ii]*counts
337                        while(jj<right)
338                while(ii<top)
339               
340                xctr = xzsum/zsum
341                yctr = yzsum/zsum
342               
343                x_mm = x_mm_sum/zsum
344                y_mm = y_mm_sum/zsum
345                // add 1 to each to get to detector coordinates (1,128)
346                // rather than the data array which is [0,127]
347//              xctr+=1
348//              yctr+=1
349
350// correct for the zero position (y-position) on the L/R panels not being exactly equal
351// the lateral scan data from Dec 2018 is used to correct this. The span of zero points
352// is relatively small (+- 0.5 pixel) but is significant for data using graphite monochromator
353//     
354        // check that the correction waves exist, if not, generate them V_TubeZeroPointTables()
355                Wave/Z tube_num = $("root:Packages:NIST:VSANS:Globals:tube_" + detStr)
356                Wave/Z yCtr_tube = $("root:Packages:NIST:VSANS:Globals:yCtr_" + detStr)
357                if(!WaveExists(tube_num))
358                        Execute "V_TubeZeroPointTables()"
359                        Wave/Z tube_num = $("root:Packages:NIST:VSANS:Globals:tube_" + detStr)
360                        Wave/Z yCtr_tube = $("root:Packages:NIST:VSANS:Globals:yCtr_" + detStr)
361                endif
362               
363                Variable yCorrection = interp(xCtr,tube_num,yCtr_tube)
364                Variable yPixSize = V_getDet_y_pixel_size(gCurDispType,detStr)
365                yPixSize /= 10          // convert mm to cm
366                // offsets were determined in Dec 2018 using:
367                // FR tube # 7 = 61.70 pix
368                // MR tube # 10 = 61.94 pix
369               
370                Print "X-center (in array coordinates 0->n-1 ) = ",xctr
371                Print "Y-center (in array coordinates 0->n-1 ) = ",yctr
372               
373                Print "X-center (cm) = ",x_mm/10
374                Print "Y-center (cm) = ",y_mm/10
375
376                if(cmpstr(detStr,"FR") == 0)
377                        Print "Reference Y-Center is corrected for tube #7 zero position"               
378
379                        yCorrection = 61.70 - yCorrection
380                        Print "yCorrection (pix) = ",yCorrection
381                        Print "yCorrection (cm) = ",yCorrection*yPixSize
382                        Print "FRONT Reference X-center (cm) = ",x_mm/10
383                        Print "FRONT Reference Y-center (cm) = ",y_mm/10 + yCorrection*yPixSize
384                endif
385
386                if(cmpstr(detStr,"MR") == 0)
387                        Print "Reference Y-Center is corrected for tube #10 zero position"             
388
389                        yCorrection = 61.94 - yCorrection
390                        Print "yCorrection (pix) = ",yCorrection
391                        Print "yCorrection (cm) = ",yCorrection*yPixSize
392                        Print "MIDDLE Reference X-center (cm) = ",x_mm/10
393                        Print "MIDDLE Reference Y-center (cm) = ",y_mm/10 + yCorrection*yPixSize
394                endif
395               
396// if measured on the LEFT panel, convert to the RIGHT coordinates for the reference value     
397// these corrections are exactly the opposite (subtract, not add) of what is done in V_fDeriveBeamCenters(xFR,yFR,xMR,yMR)
398                if(cmpstr(detStr,"FL") == 0)
399                        Print "FRONT Reference X-center (cm) = ",x_mm/10 - kBCtrOffset_FL_x     // NEW Dec 2018 values
400                        Print "FRONT Reference Y-center (cm) = ",y_mm/10 - kBCtrOffset_FL_y
401                endif
402               
403                if(cmpstr(detStr,"ML") == 0)
404                        Print "MIDDLE Reference X-center (cm) = ",x_mm/10 - kBCtrOffset_ML_x
405                        Print "MIDDLE Reference Y-center (cm) = ",y_mm/10 - kBCtrOffset_ML_y
406                endif
407        endif
408       
409        //back to root folder (redundant)
410        SetDataFolder root:
411       
412End
413
414
415
416//
417//function writes new box coordinates to the data file
418//
419// also writes the panel where the coordinates were set (non-nice field in /reduction)
420//
421Function V_UpdateBoxCoords() :  GraphMarquee
422        GetMarquee left,bottom
423        if(V_flag == 0)
424                Print "There is no Marquee"
425        else
426                Variable count,x1,x2,y1,y2,ct_err
427                x1 = V_left
428                x2 = V_right
429                y1 = V_bottom
430                y2 = V_top
431
432//              Print x1,x2,y1,y2
433
434                // NOTE:
435                // this function MODIFIES x and y values on return, converting them to panel coordinates
436                // detector panel is identified from the (left,top) coordinate (x1,y2)
437                String detStr = V_FindDetStrFromLoc(x1,x2,y1,y2)               
438//
439                SVAR gCurDispType = root:Packages:NIST:VSANS:Globals:gCurDispType
440                string boxStr
441                // this function will modify the x and y values (passed by reference) as needed to keep on the panel
442                V_KeepSelectionInBounds(x1,x2,y1,y2,detStr,gCurDispType)
443                sprintf boxStr,"%d;%d;%d;%d;",x1,x2,y1,y2
444
445//              Print x1,x2,y1,y2
446               
447                SVAR gCurrentFile = root:Packages:NIST:VSANS:Globals:gLastLoadedFile            //for the status of the display
448
449                V_writeBoxCoordinates(gCurrentFile,V_List2NumWave(boxStr,";","inW"))
450       
451                V_writeReduction_BoxPanel(gCurrentFile,detStr)
452
453//              count = V_SumCountsInBox(x1,x2,y1,y2,ct_err,gCurDispType,detStr)               
454//              Print "counts = ",count
455//              Print "err/counts = ",ct_err/count
456
457                // kill the file from RawVSANS so that the updated box coordinates will be re-read in
458                //
459                V_KillNamedDataFolder(gCurrentFile)
460        endif
461End
Note: See TracBrowser for help on using the repository browser.