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

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

more file adjustments to clean out NCNR bits

File size: 15.8 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 needs some changes after 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        //go find the file, open it and write 4 integers to the file
146        //in the positions for analysis.rows(2), .cols(2) = 4 unused 4byte integers
147        Variable refnum
148        Open/A/T="????TEXT" refnum as filename
149        FSetPos refnum,478
150        FBinWrite/F=3/B=3 refNum, x1
151        FBinWrite/F=3/B=3 refNum, x2
152        FBinWrite/F=3/B=3 refNum, y1
153        FBinWrite/F=3/B=3 refNum, y2
154        //move to the end of the file before closing
155        FStatus refnum
156        FSetPos refnum,V_logEOF
157        Close refnum
158        //ReWriteReal() takes care of open/close on its own
159        Print counts, " counts in XY box"
160        ReWriteReal(filename,counts,494)
161       
162End
163
164//finds the beam center (the centroid) of the selected region
165//and simply prints out the results to the history window
166//values are printed out in detector coordinates, not IGOR coords.
167//
168Function FindBeamCenter() :  GraphMarquee
169
170        //get the current displayed data (so the correct folder is used)
171        SVAR cur_folder=root:myGlobals:gDataDisplayType
172        String dest = "root:" + cur_folder
173       
174        Variable xzsum,yzsum,zsum,xctr,yctr
175        Variable left,right,bottom,top,ii,jj,counts
176       
177        // data wave is hard-wired in as the displayed data
178        NVAR dataIsLog=$(dest + ":gIsLogScale")         //check for log-scaling in current data folder
179        if (dataIsLog)
180                wave data=$(dest + ":linear_data")
181        else
182                wave data=$(dest + ":data")
183        endif
184       
185        GetMarquee left,bottom
186        if(V_flag == 0)
187                Print "There is no Marquee"
188        else
189                left = round(V_left)            //get integer values for selection limits
190                right = round(V_right)
191                top = round(V_top)
192                bottom = round(V_bottom)
193               
194                KeepSelectionInBounds(left,right,bottom,top)
195               
196                // selection valid now, calculate beamcenter
197                xzsum = 0
198                yzsum = 0
199                zsum = 0
200                // count over rectangular selection, doing each row, L-R, bottom to top
201                ii = bottom -1
202                do
203                        ii +=1
204                        jj = left-1
205                        do
206                                jj += 1
207                                counts = data[jj][ii]
208                                xzsum += jj*counts
209                                yzsum += ii*counts
210                                zsum += counts
211                        while(jj<right)
212                while(ii<top)
213               
214                xctr = xzsum/zsum
215                yctr = yzsum/zsum
216               
217                // add 1 to each to get to detector coordinates (1,128)
218                // rather than the data array which is [0,127]
219                xctr+=1
220                yctr+=1
221               
222                Print "X-center (in detector coordinates) = ",xctr
223                Print "Y-center (in detector coordinates) = ",yctr
224        endif
225       
226        //back to root folder (redundant)
227        SetDataFolder root:
228       
229End
230
231//still need to error check - out-of-bounds...waves exist.
232// allows a 2D Gaussian fit to a selected region of data in a SANS_Data window
233//puts up a new graph with the fitted contour
234Function Do_2D_Gaussian_Fit() :  GraphMarquee
235        String topWin=WinName(0,1)              //top *graph* window
236        //exit nicely if not in the Data display window
237        if(cmpstr(topWin,"SANS_Data") != 0)
238                DoAlert 0,"2D Gaussian fitting is only available from the Data Display Window"
239                return(1)
240        Endif
241       
242        GetMarquee/K left,bottom
243        Variable x1,x2,y1,y2,qxlo,qxhi,qylo,qyhi
244        if(V_flag == 0)
245                Print "There is no Marquee"
246        else
247                String junk="",df=""
248               
249                //**hard-wired info about the x-y q-scales
250                qxlo = DimOffset(root:myGlobals:q_x_axis,0)
251                qxhi = DimDelta(root:myGlobals:q_x_axis,0) + qxlo
252//              Print "qxlo,qxhi = ",qxlo,qxhi
253                Wave w=$"root:myGlobals:q_y_axis"
254                qylo=w[0]
255                qyhi=w[1]
256//              print "qylo,qyhi = ",qylo,qyhi
257               
258                junk=ImageInfo("SANS_Data","data",0)
259                df=StringByKey("ZWAVEDF", junk,":",";")
260//              print df
261                Duplicate/O $(df+"data") data,data_err
262                data_err=sqrt(data)             //for weighting
263               
264// comment out the SetScale lines if you want the result in terms of pixels as a way of
265// measuring the beam center. Note that you need to ADD ONE to fitted x0 and y0 to get detector
266// coordinates rather than the zero-indexed array. 2D fitting does have the benefit of
267// reporting error bars on the xy (if you believe that 2D gaussian is correct)         
268                SetScale/I x qxlo,qxhi,"",data
269                SetScale/I y qylo,qyhi,"",data
270               
271                Display /W=(10,50,361,351) /K=1
272                AppendImage data
273                ModifyImage data ctab= {*,*,Grays,0}
274                ModifyGraph width={Plan,1,bottom,left}
275                ModifyGraph mirror=2
276                ModifyGraph lowTrip=1e-04
277                ModifyImage data cindex=$"root:myGlobals:NIHColors"
278                SVAR angst = root:myGlobals:gAngstStr
279                Label bottom "Qx ("+angst+"\\S-1\\M)"
280                Label left "Qy ("+angst+"\\S-1\\M)"
281
282                //keep selection in-bounds
283                x1=V_left
284                x2=V_right
285                y1=V_bottom
286                y2=V_top
287                KeepSelectionInBounds(x1,x2,y1,y2)
288
289                //cross correlation coefficent (K6) must be between 0 and 1, need constraints
290                Make/O/T/N=2 temp_constr
291                temp_constr = {"K6>0","K6<1"}
292               
293                CurveFit/N Gauss2D data[x1,x2][y1,y2] /I=1 /W=data_err /D /R /A=0 /C=temp_constr
294               
295                Killwaves/Z temp_constr
296        endif
297End
298
299// to save the image, simply invoke the IGOR menu item for saving graphics
300//
301Function SaveSANSGraphic() : GraphMarquee
302       
303        NVAR isDemoVersion=root:myGlobals:isDemoVersion
304        if(isDemoVersion==1)
305                //      comment out in DEMO_MODIFIED version, and show the alert
306                DoAlert 0,"This operation is not available in the Demo version of IGOR"
307        else
308                DoIGORMenu "File","Save Graphics"
309        endif
310End
311
312//does a sum over each of the files in the list over the specified range
313// x,y are assumed to already be in-bounds of the data array
314// output is dumped to the command window
315//
316Function DoBoxSum(fileStr,x1,x2,y1,y2)
317        String fileStr
318        Variable x1,x2,y1,y2
319       
320        //parse the list of file numbers
321        String fileList="",item="",pathStr="",fullPath=""
322        Variable ii,num,err,cts
323       
324        PathInfo catPathName
325        If(V_Flag==0)
326                Abort "no path selected"
327        Endif
328        pathStr = S_Path
329       
330        fileList=ParseRunNumberList(fileStr)
331        num=ItemsInList(fileList,",")
332       
333        //loop over the list
334        //add each file to SAM (to normalize to monitor counts)
335        //sum over the box
336        //print the results
337        Make/O/N=(num) FileID,BoxCounts
338        Print "Results are stored in root:FileID and root:BoxCounts waves"
339        for(ii=0;ii<num;ii+=1)
340                item=StringFromList(ii,fileList,",")
341                FileID[ii] = GetRunNumFromFile(item)            //do this here, since the list is now valid
342                fullPath = pathStr+item
343                ReadHeaderAndData(fullPath)
344//              String/G root:myGlobals:gDataDisplayType="RAW"
345//              fRawWindowHook()
346                err = Raw_to_work("SAM")
347                String/G root:myGlobals:gDataDisplayType="SAM" 
348                fRawWindowHook()
349                cts=SumCountsInBox(x1,x2,y1,y2,"SAM")
350                BoxCounts[ii]=cts
351                Print item+" counts = ",cts
352        endfor
353       
354        DoBoxGraph(FileID,BoxCounts)
355       
356        return(0)
357End
358
359Function DoBoxGraph(FileID,BoxCounts)
360        Wave FileID,BoxCounts
361       
362        Sort FileID BoxCounts,FileID            //sort the waves, in case the run numbers were entered out of numerical order
363       
364        Display /W=(5,44,383,306) BoxCounts vs FileID
365        ModifyGraph mode=4
366        ModifyGraph marker=8
367        ModifyGraph grid=2
368        ModifyGraph mirror=2
369        Label left "Counts (per 10^8 monitor counts)"
370        Label bottom "Run Number"
371       
372        return(0)
373End
374// promts the user for a range of file numbers to perform the sum over
375// list must be comma delimited numbers (or dashes) just as in the BuildProtocol panel
376// the (x,y) range is already selected from the marquee
377//
378Function BoxSum() :  GraphMarquee
379        GetMarquee left,bottom
380        if(V_flag == 0)
381                Abort "There is no Marquee"
382        Endif
383        SVAR dispType=root:myGlobals:gDataDisplayType
384        if(cmpstr(dispType,"RealTime")==0)
385                Print "Can't do a BoxSum for a RealTime file"
386                return(1)
387        endif
388        Variable x1,x2,y1,y2
389        x1 = V_left
390        x2 = V_right
391        y1 = V_bottom
392        y2 = V_top
393        KeepSelectionInBounds(x1,x2,y1,y2)
394       
395        String fileStr="",msgStr="Enter a comma-delimited list of run numbers, use dashes for ranges"
396        Prompt fileStr,msgStr
397        DoPrompt "Pick the file range",fileStr
398        Print "fileStr = ",fileStr
399        printf "(x1,x2) (y1,y2) = (%d,%d) (%d,%d)\r",x1,x2,y1,y2
400       
401        DoBoxSum(fileStr,x1,x2,y1,y2)
402       
403        return(0)
404End     
405
406//function that keeps the marquee selection in the range [0,127] inclusive
407// (igor coordinate system)
408// uses pass-by reference!
409//
410// x1 = left
411// x2 = right
412// y1 = bottom
413// y2 = top
414//
415// accepts any detector size
416Function KeepSelectionInBounds(x1,x2,y1,y2)
417        Variable &x1,&x2,&y1,&y2
418       
419        NVAR pixelsX = root:myGlobals:gNPixelsX
420        NVAR pixelsY = root:myGlobals:gNPixelsY
421       
422        //keep selection in-bounds
423        x1 = (round(x1) >= 0) ? round(x1) : 0
424        x2 = (round(x2) <= (pixelsX-1)) ? round(x2) : (pixelsX-1)
425        y1 = (round(y1) >= 0) ? round(y1) : 0
426        y2 = (round(y2) <= (pixelsY-1)) ? round(y2) : (pixelsY-1)
427        return(0)
428End
429
430//testing function, not used
431Function testKeep(x1,x2,y1,y2)
432        Variable x1,x2,y1,y2
433       
434        KeepSelectionInBounds(x1,x2,y1,y2)
435        Print x1,x2,y1,y2
436        return(0)
437End
438
439// generates a histogram of the data as defined by the marquee. The longer dimension of the marquee
440// becomes the x-axis of the histogram (this may need to be changed for some odd case). Pixel range specified
441// by the marquee is inclusive, and is automatically kept in-bounds
442//
443// The counts over the (short) dimension are averaged, and plotted vs. the pixel position.
444// Pixel position is reported as Detector coordinates (1,128). Counts are whatever the current display
445// happens to be.
446//
447Function SANS_Histogram() :  GraphMarquee
448        GetMarquee left,bottom
449        if(V_flag == 0)
450                Abort "There is no Marquee"
451        endif
452        //
453        Variable count,x1,x2,y1,y2,xwidth,ywidth,vsX=1,xx,yy
454        x1 = V_left
455        x2 = V_right
456        y1 = V_bottom
457        y2 = V_top
458        KeepSelectionInBounds(x1,x2,y1,y2)
459        Print "x1,x2,y1,y2 (det) =",x1+1,x2+1,y1+1,y2+1
460        //determine whether to do x vs y or y vs x
461        xwidth=x2-x1
462        ywidth=y2-y1
463        if(xwidth < ywidth)
464                vsX=0           //sum and graph vs Y
465        endif
466        SVAR cur_folder=root:myGlobals:gDataDisplayType
467        WAVE data=$("root:"+cur_folder+":data")         //don't care if it's log or linear scale
468        Make/O/N=(max(xwidth,ywidth)+1) Position,AvgCounts
469        AvgCounts=0
470        //set position wave
471        if(vsX)
472                position=p+x1
473        else
474                position=p+y1
475        endif
476        //convert the position to Detector coordinates
477        position += 1
478       
479        //Compute the histogram (manually)
480        if(vsX)
481                for(xx=x1;xx<=x2;xx+=1)         //outer loop is the "x-axis"
482                        for(yy=y1;yy<=y2;yy+=1)
483                                AvgCounts[xx-x1] += data[xx][yy]
484                        endfor
485                endfor
486                AvgCounts /= (ywidth+1)
487        else
488                for(yy=y1;yy<=y2;yy+=1)
489                        for(xx=x1;xx<=x2;xx+=1)
490                                AvgCounts[yy-y1] += data[xx][yy]
491                        endfor
492                endfor
493                AvgCounts /= (xwidth+1)
494        endif
495        GetMarquee/K            //to keep from drawing the marquee on the new histo graph
496        //draw the graph, or just bring to the front with the new data
497        DoWindow/F SANS_Histo
498        if(V_Flag != 1)
499                Draw_Histo()
500        endif
501       
502        return(0)
503End
504
505//draws the histogram of the 2d data as specified by AvgCounts and Position
506//both wave are assumed to exist in the data folder. The SANS_Histogram() marquee
507//operation is responsible for creating them.
508//
509Function Draw_Histo()
510        Display /W=(197,329,567,461)/K=1 AvgCounts vs Position
511        DoWindow/C SANS_Histo
512        DoWindow/T SANS_Histo,"Histogram"
513        ModifyGraph mode=5,grid=1,mirror=2
514        ModifyGraph rgb=(21845,21845,21845)
515        ModifyGraph standoff=0
516        ModifyGraph hbFill=2
517        ModifyGraph useNegPat=1
518        ModifyGraph usePlusRGB=1
519        ModifyGraph useNegRGB=1
520        ModifyGraph hBarNegFill=2
521        ModifyGraph negRGB=(0,0,65535)
522        SetAxis/A/N=2 left
523        Label left "Counts"
524        Label bottom "Pixel (detector coordinates)"
525End
526
527//function will print marquee coordinates in axis terms, not detector terms
528//since IGOR is [0][127] and detector is (1,128)
529Function PrintMarqueeCoords() :  GraphMarquee
530        GetMarquee left,bottom
531        if(V_flag == 0)
532                Print "There is no Marquee"
533        else
534                Variable count,x1,x2,y1,y2
535                x1 = V_left
536                x2 = V_right
537                y1 = V_bottom
538                y2 = V_top
539                printf "marquee left in bottom axis terms: %g\r",round(V_left)
540                printf "marquee right in bottom axis terms: %g\r",round(V_right)
541                printf "marquee bottom in left axis terms: %g\r",round(V_bottom)
542                printf "marquee top in left axis terms: %g\r",round(V_top)
543                printf "**note that you must add 1 to each axis coordinate to get detector coordinates\r"
544               
545                KeepSelectionInBounds(x1,x2,y1,y2)
546                SVAR cur_folder=root:myGlobals:gDataDisplayType
547                count = SumCountsInBox(x1,x2,y1,y2,cur_folder)
548                Print "counts = "+ num2str(count)
549        endif
550End
Note: See TracBrowser for help on using the repository browser.