source: sans/SANSReduction/branches/kline_29MAR07/Put in User Procedures/SANS_Reduction_v5.00/Marquee.ipf @ 72

Last change on this file since 72 was 72, checked in by srkline, 16 years ago

more transferring of functions to NCNR_* files.

File size: 15.3 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma version=5.0
3#pragma IgorVersion=4.0
4
5////////////
6// vers 1.21 3 may 06 total transmission incorporated (BSG)
7//
8//**************************
9// Vers 1.2 091901
10//
11// marquee functions are used to:
12//
13// locate the beamcenter
14// set the (x,y) box coordinates for the sum for transmission calculation
15// read out coordinates
16// do a "box sum" of the same box over a range of files
17// do a 2D Gaussian fit over a selected range
18// do a save of the current image (with colorBar) - as a graphics image
19//
20//***************************
21
22//sums the data counts in the box specified by (x1,y1) to (x2,y2)
23//assuming that x1<x2, and y1<y2
24//the x,y values must also be in axis coordinates[0,127] NOT (1,128) detector coords.
25//
26// accepts arbitrary detector coordinates. calling function is responsible for
27// keeping selection in bounds
28Function SumCountsInBox(x1,x2,y1,y2,type)
29        Variable x1,x2,y1,y2
30        String type
31       
32        Variable counts = 0,ii,jj
33       
34        String dest =  "root:"+type
35       
36        //check for logscale data, but don't change the data
37        NVAR gIsLogScale = $(dest + ":gIsLogScale")
38        if (gIsLogScale)
39                wave w=$(dest + ":linear_data")
40        else
41                wave w=$(dest + ":data")
42        endif
43       
44        ii=x1
45        jj=y1
46        do
47                do
48                        counts += w[ii][jj]
49                        jj+=1
50                while(jj<=y2)
51                jj=y1
52                ii+=1
53        while(ii<=x2)
54       
55        Return (counts)
56End
57
58
59//from a marquee selection:
60//calculates the sum of counts in the box, and records this count value in
61//globals for both patch and Trans routines, then records the box coordinates
62//and the count value for the box in the header of the (currently displayed)
63//empty beam file for use later when calculating transmissions of samples
64//these values are written to unused (analysis) fields in the data header
65//4 integers and one real value are written
66//
67//re-written to work with Transmission panel as well as PatchPanel
68//which now bothe work from the same folder (:Patch)
69//
70Function SetXYBoxCoords() :  GraphMarquee
71
72        GetMarquee left,bottom
73        if(V_flag == 0)
74                Abort "There is no Marquee"
75        Endif
76        SVAR dispType=root:myGlobals:gDataDisplayType
77        if(cmpstr(dispType,"RealTime")==0)
78                Print "Can't SetXYBox for a RealTime file"
79                return(1)
80        endif
81                //printf "marquee left in bottom axis terms: %g\r",round(V_left)
82                //printf "marquee right in bottom axis terms: %g\r",round(V_right)
83                //printf "marquee top in left axis terms: %g\r",round(V_top)
84                //printf "marquee bottom in left axis terms: %g\r",round(V_bottom)
85        Variable x1,x2,y1,y2
86        x1 = round(V_left)
87        x2 = round(V_right)
88        y1 = round(V_bottom)
89        y2 = round(V_top)
90       
91        KeepSelectionInBounds(x1,x2,y1,y2)
92       
93        //check to make sure that Patch and Trans data folders exist for writing of global variables
94        If( ! (DataFolderExists("root:myGlobals:Patch"))  )
95                Execute "InitializePatchPanel()"
96        Endif
97        //check to make sure that Patch and Trans data folders exist for writing of global variables
98        If( ! (DataFolderExists("root:myGlobals:TransHeaderInfo"))  )
99                Execute "InitializeTransPanel()"
100        Endif
101       
102        //write string as keyword-packed string, to use IGOR parsing functions
103        String msgStr = "X1="+num2str(x1)+";"
104        msgStr += "X2="+num2str(x2)+";"
105        msgStr += "Y1="+num2str(y1)+";"
106        msgStr += "Y2="+num2str(y2)+";"
107        String/G root:myGlobals:Patch:gPS3 = msgStr
108        String/G root:myGlobals:Patch:gEmpBox = msgStr
109        //changing this global wil update the display variable on the TransPanel
110        String/G root:myGlobals:TransHeaderInfo:gBox = msgStr
111       
112        //sum the counts in the patch - working on the SAM data, to be sure that it's normalized
113        //to the same monitor counts and corrected for detector deadtime
114        String type = "SAM"
115        Variable counts
116        counts = SumCountsInBox(x1,x2,y1,y2,type)
117        Print " marquee counts =",counts
118        //Set the global gTransCts
119        Variable/G root:myGlobals:Patch:gTransCts = counts
120       
121        //now change the extra variables in the empty beam file
122        //get the filename from the SAM folder (there will only be one file)
123        SVAR partialName = root:SAM:FileList
124        //construct valid filename, then prepend path
125        String tempName = FindValidFilename(partialName)
126        Print "in marquee",partialName
127        //Print tempName
128        if(cmpstr(tempName,"")==0)
129                //file not found, get out
130                Abort "file not found, marquee"
131        Endif
132        //name is ok, prepend path to tempName for read routine
133        PathInfo catPathName
134        if (V_flag == 0)
135                //path does not exist - no folder selected
136                Abort "no path selected"
137        else
138                String filename = S_path + tempName
139        endif
140       
141        if(cmpstr(filename,"no file selected")==0)
142                Abort "no file selected"
143        Endif
144       
145        WriteXYBoxToHeader(filename,x1,x2,y1,y2)
146       
147        //ReWriteReal() takes care of open/close on its own
148        Print counts, " counts in XY box"
149        ReWriteReal(filename,counts,494)
150       
151End
152
153//finds the beam center (the centroid) of the selected region
154//and simply prints out the results to the history window
155//values are printed out in detector coordinates, not IGOR coords.
156//
157Function FindBeamCenter() :  GraphMarquee
158
159        //get the current displayed data (so the correct folder is used)
160        SVAR cur_folder=root:myGlobals:gDataDisplayType
161        String dest = "root:" + cur_folder
162       
163        Variable xzsum,yzsum,zsum,xctr,yctr
164        Variable left,right,bottom,top,ii,jj,counts
165       
166        // data wave is hard-wired in as the displayed data
167        NVAR dataIsLog=$(dest + ":gIsLogScale")         //check for log-scaling in current data folder
168        if (dataIsLog)
169                wave data=$(dest + ":linear_data")
170        else
171                wave data=$(dest + ":data")
172        endif
173       
174        GetMarquee left,bottom
175        if(V_flag == 0)
176                Print "There is no Marquee"
177        else
178                left = round(V_left)            //get integer values for selection limits
179                right = round(V_right)
180                top = round(V_top)
181                bottom = round(V_bottom)
182               
183                KeepSelectionInBounds(left,right,bottom,top)
184               
185                // selection valid now, calculate beamcenter
186                xzsum = 0
187                yzsum = 0
188                zsum = 0
189                // count over rectangular selection, doing each row, L-R, bottom to top
190                ii = bottom -1
191                do
192                        ii +=1
193                        jj = left-1
194                        do
195                                jj += 1
196                                counts = data[jj][ii]
197                                xzsum += jj*counts
198                                yzsum += ii*counts
199                                zsum += counts
200                        while(jj<right)
201                while(ii<top)
202               
203                xctr = xzsum/zsum
204                yctr = yzsum/zsum
205               
206                // add 1 to each to get to detector coordinates (1,128)
207                // rather than the data array which is [0,127]
208                xctr+=1
209                yctr+=1
210               
211                Print "X-center (in detector coordinates) = ",xctr
212                Print "Y-center (in detector coordinates) = ",yctr
213        endif
214       
215        //back to root folder (redundant)
216        SetDataFolder root:
217       
218End
219
220//still need to error check - out-of-bounds...waves exist.
221// allows a 2D Gaussian fit to a selected region of data in a SANS_Data window
222//puts up a new graph with the fitted contour
223Function Do_2D_Gaussian_Fit() :  GraphMarquee
224        String topWin=WinName(0,1)              //top *graph* window
225        //exit nicely if not in the Data display window
226        if(cmpstr(topWin,"SANS_Data") != 0)
227                DoAlert 0,"2D Gaussian fitting is only available from the Data Display Window"
228                return(1)
229        Endif
230       
231        GetMarquee/K left,bottom
232        Variable x1,x2,y1,y2,qxlo,qxhi,qylo,qyhi
233        if(V_flag == 0)
234                Print "There is no Marquee"
235        else
236                String junk="",df=""
237               
238                //**hard-wired info about the x-y q-scales
239                qxlo = DimOffset(root:myGlobals:q_x_axis,0)
240                qxhi = DimDelta(root:myGlobals:q_x_axis,0) + qxlo
241//              Print "qxlo,qxhi = ",qxlo,qxhi
242                Wave w=$"root:myGlobals:q_y_axis"
243                qylo=w[0]
244                qyhi=w[1]
245//              print "qylo,qyhi = ",qylo,qyhi
246               
247                junk=ImageInfo("SANS_Data","data",0)
248                df=StringByKey("ZWAVEDF", junk,":",";")
249//              print df
250                Duplicate/O $(df+"data") data,data_err
251                data_err=sqrt(data)             //for weighting
252               
253// comment out the SetScale lines if you want the result in terms of pixels as a way of
254// measuring the beam center. Note that you need to ADD ONE to fitted x0 and y0 to get detector
255// coordinates rather than the zero-indexed array. 2D fitting does have the benefit of
256// reporting error bars on the xy (if you believe that 2D gaussian is correct)         
257                SetScale/I x qxlo,qxhi,"",data
258                SetScale/I y qylo,qyhi,"",data
259               
260                Display /W=(10,50,361,351) /K=1
261                AppendImage data
262                ModifyImage data ctab= {*,*,Grays,0}
263                ModifyGraph width={Plan,1,bottom,left}
264                ModifyGraph mirror=2
265                ModifyGraph lowTrip=1e-04
266                ModifyImage data cindex=$"root:myGlobals:NIHColors"
267                SVAR angst = root:myGlobals:gAngstStr
268                Label bottom "Qx ("+angst+"\\S-1\\M)"
269                Label left "Qy ("+angst+"\\S-1\\M)"
270
271                //keep selection in-bounds
272                x1=V_left
273                x2=V_right
274                y1=V_bottom
275                y2=V_top
276                KeepSelectionInBounds(x1,x2,y1,y2)
277
278                //cross correlation coefficent (K6) must be between 0 and 1, need constraints
279                Make/O/T/N=2 temp_constr
280                temp_constr = {"K6>0","K6<1"}
281               
282                CurveFit/N Gauss2D data[x1,x2][y1,y2] /I=1 /W=data_err /D /R /A=0 /C=temp_constr
283               
284                Killwaves/Z temp_constr
285        endif
286End
287
288// to save the image, simply invoke the IGOR menu item for saving graphics
289//
290Function SaveSANSGraphic() : GraphMarquee
291       
292        NVAR isDemoVersion=root:myGlobals:isDemoVersion
293        if(isDemoVersion==1)
294                //      comment out in DEMO_MODIFIED version, and show the alert
295                DoAlert 0,"This operation is not available in the Demo version of IGOR"
296        else
297                DoIGORMenu "File","Save Graphics"
298        endif
299End
300
301//does a sum over each of the files in the list over the specified range
302// x,y are assumed to already be in-bounds of the data array
303// output is dumped to the command window
304//
305Function DoBoxSum(fileStr,x1,x2,y1,y2)
306        String fileStr
307        Variable x1,x2,y1,y2
308       
309        //parse the list of file numbers
310        String fileList="",item="",pathStr="",fullPath=""
311        Variable ii,num,err,cts
312       
313        PathInfo catPathName
314        If(V_Flag==0)
315                Abort "no path selected"
316        Endif
317        pathStr = S_Path
318       
319        fileList=ParseRunNumberList(fileStr)
320        num=ItemsInList(fileList,",")
321       
322        //loop over the list
323        //add each file to SAM (to normalize to monitor counts)
324        //sum over the box
325        //print the results
326        Make/O/N=(num) FileID,BoxCounts
327        Print "Results are stored in root:FileID and root:BoxCounts waves"
328        for(ii=0;ii<num;ii+=1)
329                item=StringFromList(ii,fileList,",")
330                FileID[ii] = GetRunNumFromFile(item)            //do this here, since the list is now valid
331                fullPath = pathStr+item
332                ReadHeaderAndData(fullPath)
333//              String/G root:myGlobals:gDataDisplayType="RAW"
334//              fRawWindowHook()
335                err = Raw_to_work("SAM")
336                String/G root:myGlobals:gDataDisplayType="SAM" 
337                fRawWindowHook()
338                cts=SumCountsInBox(x1,x2,y1,y2,"SAM")
339                BoxCounts[ii]=cts
340                Print item+" counts = ",cts
341        endfor
342       
343        DoBoxGraph(FileID,BoxCounts)
344       
345        return(0)
346End
347
348Function DoBoxGraph(FileID,BoxCounts)
349        Wave FileID,BoxCounts
350       
351        Sort FileID BoxCounts,FileID            //sort the waves, in case the run numbers were entered out of numerical order
352       
353        Display /W=(5,44,383,306) BoxCounts vs FileID
354        ModifyGraph mode=4
355        ModifyGraph marker=8
356        ModifyGraph grid=2
357        ModifyGraph mirror=2
358        Label left "Counts (per 10^8 monitor counts)"
359        Label bottom "Run Number"
360       
361        return(0)
362End
363// promts the user for a range of file numbers to perform the sum over
364// list must be comma delimited numbers (or dashes) just as in the BuildProtocol panel
365// the (x,y) range is already selected from the marquee
366//
367Function BoxSum() :  GraphMarquee
368        GetMarquee left,bottom
369        if(V_flag == 0)
370                Abort "There is no Marquee"
371        Endif
372        SVAR dispType=root:myGlobals:gDataDisplayType
373        if(cmpstr(dispType,"RealTime")==0)
374                Print "Can't do a BoxSum for a RealTime file"
375                return(1)
376        endif
377        Variable x1,x2,y1,y2
378        x1 = V_left
379        x2 = V_right
380        y1 = V_bottom
381        y2 = V_top
382        KeepSelectionInBounds(x1,x2,y1,y2)
383       
384        String fileStr="",msgStr="Enter a comma-delimited list of run numbers, use dashes for ranges"
385        Prompt fileStr,msgStr
386        DoPrompt "Pick the file range",fileStr
387        Print "fileStr = ",fileStr
388        printf "(x1,x2) (y1,y2) = (%d,%d) (%d,%d)\r",x1,x2,y1,y2
389       
390        DoBoxSum(fileStr,x1,x2,y1,y2)
391       
392        return(0)
393End     
394
395//function that keeps the marquee selection in the range [0,127] inclusive
396// (igor coordinate system)
397// uses pass-by reference!
398//
399// x1 = left
400// x2 = right
401// y1 = bottom
402// y2 = top
403//
404// accepts any detector size
405Function KeepSelectionInBounds(x1,x2,y1,y2)
406        Variable &x1,&x2,&y1,&y2
407       
408        NVAR pixelsX = root:myGlobals:gNPixelsX
409        NVAR pixelsY = root:myGlobals:gNPixelsY
410       
411        //keep selection in-bounds
412        x1 = (round(x1) >= 0) ? round(x1) : 0
413        x2 = (round(x2) <= (pixelsX-1)) ? round(x2) : (pixelsX-1)
414        y1 = (round(y1) >= 0) ? round(y1) : 0
415        y2 = (round(y2) <= (pixelsY-1)) ? round(y2) : (pixelsY-1)
416        return(0)
417End
418
419//testing function, not used
420Function testKeep(x1,x2,y1,y2)
421        Variable x1,x2,y1,y2
422       
423        KeepSelectionInBounds(x1,x2,y1,y2)
424        Print x1,x2,y1,y2
425        return(0)
426End
427
428// generates a histogram of the data as defined by the marquee. The longer dimension of the marquee
429// becomes the x-axis of the histogram (this may need to be changed for some odd case). Pixel range specified
430// by the marquee is inclusive, and is automatically kept in-bounds
431//
432// The counts over the (short) dimension are averaged, and plotted vs. the pixel position.
433// Pixel position is reported as Detector coordinates (1,128). Counts are whatever the current display
434// happens to be.
435//
436Function SANS_Histogram() :  GraphMarquee
437        GetMarquee left,bottom
438        if(V_flag == 0)
439                Abort "There is no Marquee"
440        endif
441        //
442        Variable count,x1,x2,y1,y2,xwidth,ywidth,vsX=1,xx,yy
443        x1 = V_left
444        x2 = V_right
445        y1 = V_bottom
446        y2 = V_top
447        KeepSelectionInBounds(x1,x2,y1,y2)
448        Print "x1,x2,y1,y2 (det) =",x1+1,x2+1,y1+1,y2+1
449        //determine whether to do x vs y or y vs x
450        xwidth=x2-x1
451        ywidth=y2-y1
452        if(xwidth < ywidth)
453                vsX=0           //sum and graph vs Y
454        endif
455        SVAR cur_folder=root:myGlobals:gDataDisplayType
456        WAVE data=$("root:"+cur_folder+":data")         //don't care if it's log or linear scale
457        Make/O/N=(max(xwidth,ywidth)+1) Position,AvgCounts
458        AvgCounts=0
459        //set position wave
460        if(vsX)
461                position=p+x1
462        else
463                position=p+y1
464        endif
465        //convert the position to Detector coordinates
466        position += 1
467       
468        //Compute the histogram (manually)
469        if(vsX)
470                for(xx=x1;xx<=x2;xx+=1)         //outer loop is the "x-axis"
471                        for(yy=y1;yy<=y2;yy+=1)
472                                AvgCounts[xx-x1] += data[xx][yy]
473                        endfor
474                endfor
475                AvgCounts /= (ywidth+1)
476        else
477                for(yy=y1;yy<=y2;yy+=1)
478                        for(xx=x1;xx<=x2;xx+=1)
479                                AvgCounts[yy-y1] += data[xx][yy]
480                        endfor
481                endfor
482                AvgCounts /= (xwidth+1)
483        endif
484        GetMarquee/K            //to keep from drawing the marquee on the new histo graph
485        //draw the graph, or just bring to the front with the new data
486        DoWindow/F SANS_Histo
487        if(V_Flag != 1)
488                Draw_Histo()
489        endif
490       
491        return(0)
492End
493
494//draws the histogram of the 2d data as specified by AvgCounts and Position
495//both wave are assumed to exist in the data folder. The SANS_Histogram() marquee
496//operation is responsible for creating them.
497//
498Function Draw_Histo()
499        Display /W=(197,329,567,461)/K=1 AvgCounts vs Position
500        DoWindow/C SANS_Histo
501        DoWindow/T SANS_Histo,"Histogram"
502        ModifyGraph mode=5,grid=1,mirror=2
503        ModifyGraph rgb=(21845,21845,21845)
504        ModifyGraph standoff=0
505        ModifyGraph hbFill=2
506        ModifyGraph useNegPat=1
507        ModifyGraph usePlusRGB=1
508        ModifyGraph useNegRGB=1
509        ModifyGraph hBarNegFill=2
510        ModifyGraph negRGB=(0,0,65535)
511        SetAxis/A/N=2 left
512        Label left "Counts"
513        Label bottom "Pixel (detector coordinates)"
514End
515
516//function will print marquee coordinates in axis terms, not detector terms
517//since IGOR is [0][127] and detector is (1,128)
518Function PrintMarqueeCoords() :  GraphMarquee
519        GetMarquee left,bottom
520        if(V_flag == 0)
521                Print "There is no Marquee"
522        else
523                Variable count,x1,x2,y1,y2
524                x1 = V_left
525                x2 = V_right
526                y1 = V_bottom
527                y2 = V_top
528                printf "marquee left in bottom axis terms: %g\r",round(V_left)
529                printf "marquee right in bottom axis terms: %g\r",round(V_right)
530                printf "marquee bottom in left axis terms: %g\r",round(V_bottom)
531                printf "marquee top in left axis terms: %g\r",round(V_top)
532                printf "**note that you must add 1 to each axis coordinate to get detector coordinates\r"
533               
534                KeepSelectionInBounds(x1,x2,y1,y2)
535                SVAR cur_folder=root:myGlobals:gDataDisplayType
536                count = SumCountsInBox(x1,x2,y1,y2,cur_folder)
537                Print "counts = "+ num2str(count)
538        endif
539End
Note: See TracBrowser for help on using the repository browser.