source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_Sector_Average.ipf @ 1148

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

added a batch patch routine to correct the sample aperture shape and size since it may be incorrectly written by NICE

defined the sector angles to match the SANS definition and updated the associated logic for the averaging and sector display

updated sorting of the file catalog to allow sorting by SDD_F and countRate_F. Also added a second sort key to keep the run numbers in order whenever sorting any other column.

File size: 35.3 KB
Line 
1#pragma TextEncoding = "MacRoman"
2#pragma rtGlobals=3             // Use modern global access method and strict wave access.
3
4//
5// routines to do a sector average
6//
7// possible ideas:
8// - start with calculation of the phi matrix for each panel
9// - with the definition of the angle +/-, I can decide which points to keep during the average.
10// - I can also make this an actual mask (sum with the protocol mask) and use it this way
11//  so that I can still run it through the circular average routine.
12// - I can also then overlay the sector mask onto the data (once I figure out how to overlay masks on
13//  the regular data display
14
15
16
17Function/WAVE V_MakePhiMatrix(qTotal,folderStr,detStr,folderPath)
18        Wave qTotal
19        String folderStr,detStr,folderPath
20
21
22        Variable xctr = V_getDet_beam_center_x_pix(folderStr,detStr)
23        Variable yctr = V_getDet_beam_center_y_pix(folderStr,detStr)
24
25        Duplicate/O qTotal,$(folderPath+":phi")
26        Wave phi = $(folderPath+":phi")
27        Variable pixSizeX,pixSizeY
28        pixSizeX = V_getDet_x_pixel_size(folderStr,detStr)
29        pixSizeY = V_getDet_y_pixel_size(folderStr,detStr)
30        phi = V_FindPhi( pixSizeX*((p+1)-xctr) , pixSizeY*((q+1)-yctr))         //(dx,dy)
31       
32        return phi     
33End
34
35
36//
37// x- I want to mask out everything that is "out" of the sector
38//
39// 0 = keep the point
40// 1 = yes, mask the point
41//
42//
43// phiCtr is in the range (-90,90) degrees
44// delta is in the range (0,90) for a total width of 2*delta = 180 degrees
45//
46Function V_MarkSectorOverlayPixels(phi,overlay,phiCtr,delta,side)
47        Wave phi,overlay
48        Variable phiCtr,delta
49        String side
50       
51        Variable phiVal
52
53// convert the imput from degrees to radians    , since phi is in radians
54        phiCtr *= pi/180
55        delta *= pi/180         
56       
57        Variable xDim=DimSize(phi, 0)
58        Variable yDim=DimSize(phi, 1)
59
60        Variable ii,jj,exclude,mirror_phiCtr,crossZero,keepPix
61       
62// initialize the mask to == 1 == exclude everything
63        overlay = 1
64
65// now give every opportunity to keep pixel in
66// comparisons use a modified phiCtr to match the definition of the phi field (0= +x-axis)
67//
68        for(ii=0;ii<xDim;ii+=1)
69                for(jj=0;jj<yDim;jj+=1)
70                        //qTot = sqrt(qx[ii]^2 + qy[ii]^2+ qz[ii]^2)
71                        phiVal = phi[ii][jj]
72                        keepPix = 0             //start with not keeping
73
74                        // if within the right or left, flag to keep the pixel
75                        if(cmpstr(side,"right")==0)
76                                //right, when 0->pi/2
77                                if(V_CloseEnough(phiVal,phiCtr,delta))
78                                        keepPix = 1
79                                endif
80                                // condition here to get the 3pi/2 -> 2pi region
81                                if(V_CloseEnough(phiVal,phiCtr+2*pi,delta))
82                                        keepPix = 1
83                                endif
84                        endif
85                       
86                        if(cmpstr(side,"left")==0)
87                                if(V_CloseEnough(phiVal,phiCtr+pi,delta))
88                                        keepPix = 1
89                                endif
90                        endif
91                                               
92                //      both sides, duplicates the conditions above
93                        if(cmpstr(side,"both")==0)     
94                                //right, when 0->pi/2
95                                if(V_CloseEnough(phiVal,phiCtr,delta))
96                                        keepPix = 1
97                                endif
98                                // right, when 3pi/2 -> 2pi
99                                if(V_CloseEnough(phiVal,phiCtr+2*pi,delta))
100                                        keepPix = 1
101                                endif                           
102                               
103                                //left
104                                if(V_CloseEnough(phiVal,phiCtr+pi,delta))
105                                        keepPix = 1
106                                endif
107                               
108                        endif
109                               
110                        // set the mask value (entire overlay initialized to 1 to start)
111                        if(keepPix > 0)
112                                overlay[ii][jj] = 0
113                        endif
114                       
115                endfor
116        endfor
117
118
119        return(0)
120End
121
122       
123
124//
125// TODO -- binType == 4 (slit mode) should never end up here
126// -- new logic in calling routines to dispatch to proper routine
127// -- AND need to write the routine for binning_SlitMode
128//
129// side = one of "left;right;both;"
130// phi_rad = center of sector in radians
131// dphi_rad = half-width of sector, also in radians
132Function V_QBinAllPanels_Sector(folderStr,binType,collimationStr,side,phi_rad,dphi_rad)
133        String folderStr
134        Variable binType
135        String collimationStr,side
136        Variable phi_rad,dphi_rad
137
138        // do the back, middle, and front separately
139       
140//      figure out the binning type (where is it set?)
141        Variable ii,delQ
142        String detStr
143
144//      binType = V_GetBinningPopMode()
145
146        // set delta Q for binning (used later inside VC_fDoBinning_QxQy2D)
147        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
148                detStr = StringFromList(ii, ksDetectorListAll, ";")
149               
150                delQ = V_SetDeltaQ(folderStr,detStr)            // this sets (overwrites) the global value for each panel type
151        endfor
152       
153
154        switch(binType)
155                case 1:
156                        V_fDoSectorBin_QxQy2D(folderStr,"FL",collimationStr,side,phi_rad,dphi_rad)
157                        V_fDoSectorBin_QxQy2D(folderStr,"FR",collimationStr,side,phi_rad,dphi_rad)
158                        V_fDoSectorBin_QxQy2D(folderStr,"FT",collimationStr,side,phi_rad,dphi_rad)
159                        V_fDoSectorBin_QxQy2D(folderStr,"FB",collimationStr,side,phi_rad,dphi_rad)
160                        V_fDoSectorBin_QxQy2D(folderStr,"ML",collimationStr,side,phi_rad,dphi_rad)
161                        V_fDoSectorBin_QxQy2D(folderStr,"MR",collimationStr,side,phi_rad,dphi_rad)
162                        V_fDoSectorBin_QxQy2D(folderStr,"MT",collimationStr,side,phi_rad,dphi_rad)
163                        V_fDoSectorBin_QxQy2D(folderStr,"MB",collimationStr,side,phi_rad,dphi_rad)                     
164                        V_fDoSectorBin_QxQy2D(folderStr, "B",collimationStr,side,phi_rad,dphi_rad)             
165
166                        break
167                case 2:
168                        V_fDoSectorBin_QxQy2D(folderStr,"FLR",collimationStr,side,phi_rad,dphi_rad)
169                        V_fDoSectorBin_QxQy2D(folderStr,"FTB",collimationStr,side,phi_rad,dphi_rad)
170                        V_fDoSectorBin_QxQy2D(folderStr,"MLR",collimationStr,side,phi_rad,dphi_rad)
171                        V_fDoSectorBin_QxQy2D(folderStr,"MTB",collimationStr,side,phi_rad,dphi_rad)
172                        V_fDoSectorBin_QxQy2D(folderStr, "B",collimationStr,side,phi_rad,dphi_rad)             
173
174                        break
175                case 3:
176                        V_fDoSectorBin_QxQy2D(folderStr,"MLRTB",collimationStr,side,phi_rad,dphi_rad)
177                        V_fDoSectorBin_QxQy2D(folderStr,"FLRTB",collimationStr,side,phi_rad,dphi_rad)
178                        V_fDoSectorBin_QxQy2D(folderStr, "B",collimationStr,side,phi_rad,dphi_rad)             
179                       
180                        break
181                case 4:                         /// this is for a tall, narrow slit mode       
182                        VC_fBinDetector_byRows(folderStr,"FL")
183                        VC_fBinDetector_byRows(folderStr,"FR")
184                        VC_fBinDetector_byRows(folderStr,"ML")
185                        VC_fBinDetector_byRows(folderStr,"MR")
186                        VC_fBinDetector_byRows(folderStr,"B")
187
188                        break
189                case 5:
190                        V_fDoSectorBin_QxQy2D(folderStr,"FTB",collimationStr,side,phi_rad,dphi_rad)
191                        V_fDoSectorBin_QxQy2D(folderStr,"FLR",collimationStr,side,phi_rad,dphi_rad)
192                        V_fDoSectorBin_QxQy2D(folderStr,"MLRTB",collimationStr,side,phi_rad,dphi_rad)
193                        V_fDoSectorBin_QxQy2D(folderStr, "B",collimationStr,side,phi_rad,dphi_rad)             
194               
195                        break
196                case 6:
197                        V_fDoSectorBin_QxQy2D(folderStr,"FLRTB",collimationStr,side,phi_rad,dphi_rad)
198                        V_fDoSectorBin_QxQy2D(folderStr,"MLR",collimationStr,side,phi_rad,dphi_rad)
199                        V_fDoSectorBin_QxQy2D(folderStr, "B",collimationStr,side,phi_rad,dphi_rad)             
200               
201                        break
202                case 7:
203                        V_fDoSectorBin_QxQy2D(folderStr,"FTB",collimationStr,side,phi_rad,dphi_rad)
204                        V_fDoSectorBin_QxQy2D(folderStr,"FLR",collimationStr,side,phi_rad,dphi_rad)
205                        V_fDoSectorBin_QxQy2D(folderStr,"MLR",collimationStr,side,phi_rad,dphi_rad)
206                        V_fDoSectorBin_QxQy2D(folderStr, "B",collimationStr,side,phi_rad,dphi_rad)             
207               
208                        break
209                       
210                default:
211                        Abort "Binning mode not found in V_QBinAllPanels_Circular"// when no case matches       
212        endswitch
213       
214
215        return(0)
216End
217
218
219//////////
220//
221//              Function that bins a 2D detctor panel into I(q) based on the q-value of the pixel
222//              - each pixel QxQyQz has been calculated beforehand
223//              - if multiple panels are selected to be combined, it is done here during the binning
224//              - the setting of deltaQ step is still a little suspect (TODO)
225//
226//
227// see the equivalent function in PlotUtils2D_v40.ipf
228//
229//Function fDoBinning_QxQy2D(inten,qx,qy,qz)
230//
231// this has been modified to accept different detector panels and to take arrays
232// -- type = FL or FR or...other panel identifiers
233//
234// TODO "iErr" is not always defined correctly since it doesn't really apply here for data that is not 2D simulation
235//
236//
237// updated Feb2016 to take new folder structure
238// TODO
239// -- VERIFY
240// -- figure out what the best location is to put the averaged data? currently @ top level of WORK folder
241//    but this is a lousy choice.
242// x- binning is now Mask-aware. If mask is not present, all data is used. If data is from VCALC, all data is used
243// x- Where do I put the solid angle correction? In here as a weight for each point, or later on as
244//    a blanket correction (matrix multiply) for an entire panel? (Solid Angle correction is done in the
245//    step where data is added to a WORK file (see Raw_to_Work())
246//
247//
248// TODO:
249// -- some of the input parameters for the resolution calcuation are either assumed (apOff) or are currently
250//    hard-wired. these need to be corrected before even the pinhole resolution is correct
251// x- resolution calculation is in the correct place. The calculation is done per-panel (specified by TYPE),
252//    and then the unwanted points can be discarded (all 6 columns) as the data is trimmed and concatenated
253//    is separate functions that are resolution-aware.
254//
255//
256// folderStr = WORK folder, type = the binning type (may include multiple detectors)
257//
258// side = one of "left;right;both;"
259// phi_rad = center of sector in radians
260// dphi_rad = half-width of sector, also in radians
261//
262Function V_fDoSectorBin_QxQy2D(folderStr,type,collimationStr,side,phi_rad,dphi_rad)
263        String folderStr,type,collimationStr,side
264        Variable phi_rad,dphi_rad
265       
266        Variable nSets = 0
267        Variable xDim,yDim
268        Variable ii,jj
269        Variable qVal,nq,var,avesq,aveisq
270        Variable binIndex,val,isVCALC=0,maskMissing
271
272        String folderPath = "root:Packages:NIST:VSANS:"+folderStr
273        String instPath = ":entry:instrument:detector_"
274        String detStr
275               
276        if(cmpstr(folderStr,"VCALC") == 0)
277                isVCALC = 1
278        endif
279       
280        detStr = type
281
282// now switch on the type to determine which waves to declare and create
283// since there may be more than one panel to step through. There may be two, there may be four
284//
285
286// TODO:
287// -- Solid_Angle -- waves will be present for WORK data other than RAW, but not for RAW
288//
289// assume that the mask files are missing unless we can find them. If VCALC data,
290//  then the Mask is missing by definition
291        maskMissing = 1
292
293        strswitch(type) // string switch
294
295// only one panel, simply pick that panel and move on out of the switch
296                case "FL":
297                case "FR":
298                case "FT":
299                case "FB":
300                case "ML":
301                case "MR":
302                case "MT":
303                case "MB":                     
304                case "B":       
305                        if(isVCALC)
306                                WAVE inten = $(folderPath+instPath+detStr+":det_"+detStr)
307                                WAVE/Z iErr = $("iErr_"+detStr)                 // 2D errors -- may not exist, especially for simulation               
308                        else
309                                Wave inten = V_getDetectorDataW(folderStr,detStr)
310                                Wave iErr = V_getDetectorDataErrW(folderStr,detStr)
311                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+detStr+":data")
312                                if(WaveExists(mask) == 1)
313                                        maskMissing = 0
314                                endif
315                        endif   
316                        NVAR delQ = $(folderPath+instPath+detStr+":gDelQ_"+detStr)
317                        Wave qTotal = $(folderPath+instPath+detStr+":qTot_"+detStr)                     // 2D q-values 
318                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,detStr,folderPath+instPath+detStr)
319                        nSets = 1
320                        break   
321                       
322                case "FLR":
323                // detStr has multiple values now, so unfortuntely, I'm hard-wiring things...
324                // TODO
325                // -- see if I can un-hard-wire some of this below when more than one panel is combined
326                        if(isVCALC)
327                                WAVE inten = $(folderPath+instPath+"FL"+":det_"+"FL")
328                                WAVE/Z iErr = $("iErr_"+"FL")                   // 2D errors -- may not exist, especially for simulation               
329                                WAVE inten2 = $(folderPath+instPath+"FR"+":det_"+"FR")
330                                WAVE/Z iErr2 = $("iErr_"+"FR")                  // 2D errors -- may not exist, especially for simulation       
331                        else
332                                Wave inten = V_getDetectorDataW(folderStr,"FL")
333                                Wave iErr = V_getDetectorDataErrW(folderStr,"FL")
334                                Wave inten2 = V_getDetectorDataW(folderStr,"FR")
335                                Wave iErr2 = V_getDetectorDataErrW(folderStr,"FR")
336                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FL"+":data")
337                                Wave/Z mask2 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FR"+":data")
338                                if(WaveExists(mask) == 1 && WaveExists(mask2) == 1)
339                                        maskMissing = 0
340                                endif
341                        endif   
342                        NVAR delQ = $(folderPath+instPath+"FL"+":gDelQ_FL")
343                       
344                        Wave qTotal = $(folderPath+instPath+"FL"+":qTot_"+"FL")                 // 2D q-values 
345                        Wave qTotal2 = $(folderPath+instPath+"FR"+":qTot_"+"FR")                        // 2D q-values 
346                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,"FL",folderPath+instPath+"FL")
347                        Wave phi2 = V_MakePhiMatrix(qTotal2,folderStr,"FR",folderPath+instPath+"FR")
348
349                        nSets = 2
350                        break                   
351               
352                case "FTB":
353                        if(isVCALC)
354                                WAVE inten = $(folderPath+instPath+"FT"+":det_"+"FT")
355                                WAVE/Z iErr = $("iErr_"+"FT")                   // 2D errors -- may not exist, especially for simulation               
356                                WAVE inten2 = $(folderPath+instPath+"FB"+":det_"+"FB")
357                                WAVE/Z iErr2 = $("iErr_"+"FB")                  // 2D errors -- may not exist, especially for simulation       
358                        else
359                                Wave inten = V_getDetectorDataW(folderStr,"FT")
360                                Wave iErr = V_getDetectorDataErrW(folderStr,"FT")
361                                Wave inten2 = V_getDetectorDataW(folderStr,"FB")
362                                Wave iErr2 = V_getDetectorDataErrW(folderStr,"FB")
363                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FT"+":data")
364                                Wave/Z mask2 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FB"+":data")
365                                if(WaveExists(mask) == 1 && WaveExists(mask2) == 1)
366                                        maskMissing = 0
367                                endif
368                        endif   
369                        NVAR delQ = $(folderPath+instPath+"FT"+":gDelQ_FT")
370                       
371                        Wave qTotal = $(folderPath+instPath+"FT"+":qTot_"+"FT")                 // 2D q-values 
372                        Wave qTotal2 = $(folderPath+instPath+"FB"+":qTot_"+"FB")                        // 2D q-values 
373
374                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,"FT",folderPath+instPath+"FT")
375                        Wave phi2 = V_MakePhiMatrix(qTotal2,folderStr,"FB",folderPath+instPath+"FB")
376                               
377                        nSets = 2
378                        break           
379               
380                case "FLRTB":
381                        if(isVCALC)
382                                WAVE inten = $(folderPath+instPath+"FL"+":det_"+"FL")
383                                WAVE/Z iErr = $("iErr_"+"FL")                   // 2D errors -- may not exist, especially for simulation               
384                                WAVE inten2 = $(folderPath+instPath+"FR"+":det_"+"FR")
385                                WAVE/Z iErr2 = $("iErr_"+"FR")                  // 2D errors -- may not exist, especially for simulation       
386                                WAVE inten3 = $(folderPath+instPath+"FT"+":det_"+"FT")
387                                WAVE/Z iErr3 = $("iErr_"+"FT")                  // 2D errors -- may not exist, especially for simulation               
388                                WAVE inten4 = $(folderPath+instPath+"FB"+":det_"+"FB")
389                                WAVE/Z iErr4 = $("iErr_"+"FB")                  // 2D errors -- may not exist, especially for simulation       
390                        else
391                                Wave inten = V_getDetectorDataW(folderStr,"FL")
392                                Wave iErr = V_getDetectorDataErrW(folderStr,"FL")
393                                Wave inten2 = V_getDetectorDataW(folderStr,"FR")
394                                Wave iErr2 = V_getDetectorDataErrW(folderStr,"FR")
395                                Wave inten3 = V_getDetectorDataW(folderStr,"FT")
396                                Wave iErr3 = V_getDetectorDataErrW(folderStr,"FT")
397                                Wave inten4 = V_getDetectorDataW(folderStr,"FB")
398                                Wave iErr4 = V_getDetectorDataErrW(folderStr,"FB")
399                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FL"+":data")
400                                Wave/Z mask2 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FR"+":data")
401                                Wave/Z mask3 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FT"+":data")
402                                Wave/Z mask4 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"FB"+":data")
403                                if(WaveExists(mask) == 1 && WaveExists(mask2) == 1 && WaveExists(mask3) == 1 && WaveExists(mask4) == 1)
404                                        maskMissing = 0
405                                endif
406                        endif   
407                        NVAR delQ = $(folderPath+instPath+"FL"+":gDelQ_FL")
408                       
409                        Wave qTotal = $(folderPath+instPath+"FL"+":qTot_"+"FL")                 // 2D q-values 
410                        Wave qTotal2 = $(folderPath+instPath+"FR"+":qTot_"+"FR")                        // 2D q-values 
411                        Wave qTotal3 = $(folderPath+instPath+"FT"+":qTot_"+"FT")                        // 2D q-values 
412                        Wave qTotal4 = $(folderPath+instPath+"FB"+":qTot_"+"FB")                        // 2D q-values 
413
414                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,"FL",folderPath+instPath+"FL")
415                        Wave phi2 = V_MakePhiMatrix(qTotal2,folderStr,"FR",folderPath+instPath+"FR")
416                        Wave phi3 = V_MakePhiMatrix(qTotal3,folderStr,"FT",folderPath+instPath+"FT")
417                        Wave phi4 = V_MakePhiMatrix(qTotal4,folderStr,"FB",folderPath+instPath+"FB")   
418                               
419                        nSets = 4
420                        break           
421                       
422                case "MLR":
423                        if(isVCALC)
424                                WAVE inten = $(folderPath+instPath+"ML"+":det_"+"ML")
425                                WAVE/Z iErr = $("iErr_"+"ML")                   // 2D errors -- may not exist, especially for simulation               
426                                WAVE inten2 = $(folderPath+instPath+"MR"+":det_"+"MR")
427                                WAVE/Z iErr2 = $("iErr_"+"MR")                  // 2D errors -- may not exist, especially for simulation       
428                        else
429                                Wave inten = V_getDetectorDataW(folderStr,"ML")
430                                Wave iErr = V_getDetectorDataErrW(folderStr,"ML")
431                                Wave inten2 = V_getDetectorDataW(folderStr,"MR")
432                                Wave iErr2 = V_getDetectorDataErrW(folderStr,"MR")
433                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"ML"+":data")
434                                Wave/Z mask2 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"MR"+":data")
435                                if(WaveExists(mask) == 1 && WaveExists(mask2) == 1)
436                                        maskMissing = 0
437                                endif
438                        endif   
439                        NVAR delQ = $(folderPath+instPath+"ML"+":gDelQ_ML")
440                       
441                        Wave qTotal = $(folderPath+instPath+"ML"+":qTot_"+"ML")                 // 2D q-values 
442                        Wave qTotal2 = $(folderPath+instPath+"MR"+":qTot_"+"MR")                        // 2D q-values 
443
444                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,"ML",folderPath+instPath+"ML")
445                        Wave phi2 = V_MakePhiMatrix(qTotal2,folderStr,"MR",folderPath+instPath+"MR")
446                                       
447                        nSets = 2
448                        break                   
449               
450                case "MTB":
451                        if(isVCALC)
452                                WAVE inten = $(folderPath+instPath+"MT"+":det_"+"MT")
453                                WAVE/Z iErr = $("iErr_"+"MT")                   // 2D errors -- may not exist, especially for simulation               
454                                WAVE inten2 = $(folderPath+instPath+"MB"+":det_"+"MB")
455                                WAVE/Z iErr2 = $("iErr_"+"MB")                  // 2D errors -- may not exist, especially for simulation       
456                        else
457                                Wave inten = V_getDetectorDataW(folderStr,"MT")
458                                Wave iErr = V_getDetectorDataErrW(folderStr,"MT")
459                                Wave inten2 = V_getDetectorDataW(folderStr,"MB")
460                                Wave iErr2 = V_getDetectorDataErrW(folderStr,"MB")
461                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"MT"+":data")
462                                Wave/Z mask2 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"MB"+":data")
463                                if(WaveExists(mask) == 1 && WaveExists(mask2) == 1)
464                                        maskMissing = 0
465                                endif
466                        endif   
467                        NVAR delQ = $(folderPath+instPath+"MT"+":gDelQ_MT")
468                       
469                        Wave qTotal = $(folderPath+instPath+"MT"+":qTot_"+"MT")                 // 2D q-values 
470                        Wave qTotal2 = $(folderPath+instPath+"MB"+":qTot_"+"MB")                        // 2D q-values 
471
472                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,"MT",folderPath+instPath+"MT")
473                        Wave phi2 = V_MakePhiMatrix(qTotal2,folderStr,"MB",folderPath+instPath+"MB")
474                                       
475                        nSets = 2
476                        break                           
477               
478                case "MLRTB":
479                        if(isVCALC)
480                                WAVE inten = $(folderPath+instPath+"ML"+":det_"+"ML")
481                                WAVE/Z iErr = $("iErr_"+"ML")                   // 2D errors -- may not exist, especially for simulation               
482                                WAVE inten2 = $(folderPath+instPath+"MR"+":det_"+"MR")
483                                WAVE/Z iErr2 = $("iErr_"+"MR")                  // 2D errors -- may not exist, especially for simulation       
484                                WAVE inten3 = $(folderPath+instPath+"MT"+":det_"+"MT")
485                                WAVE/Z iErr3 = $("iErr_"+"MT")                  // 2D errors -- may not exist, especially for simulation               
486                                WAVE inten4 = $(folderPath+instPath+"MB"+":det_"+"MB")
487                                WAVE/Z iErr4 = $("iErr_"+"MB")                  // 2D errors -- may not exist, especially for simulation       
488                        else
489                                Wave inten = V_getDetectorDataW(folderStr,"ML")
490                                Wave iErr = V_getDetectorDataErrW(folderStr,"ML")
491                                Wave inten2 = V_getDetectorDataW(folderStr,"MR")
492                                Wave iErr2 = V_getDetectorDataErrW(folderStr,"MR")
493                                Wave inten3 = V_getDetectorDataW(folderStr,"MT")
494                                Wave iErr3 = V_getDetectorDataErrW(folderStr,"MT")
495                                Wave inten4 = V_getDetectorDataW(folderStr,"MB")
496                                Wave iErr4 = V_getDetectorDataErrW(folderStr,"MB")
497                                Wave/Z mask = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"ML"+":data")
498                                Wave/Z mask2 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"MR"+":data")
499                                Wave/Z mask3 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"MT"+":data")
500                                Wave/Z mask4 = $("root:Packages:NIST:VSANS:MSK:entry:instrument:detector_"+"MB"+":data")
501                                if(WaveExists(mask) == 1 && WaveExists(mask2) == 1 && WaveExists(mask3) == 1 && WaveExists(mask4) == 1)
502                                        maskMissing = 0
503                                endif
504                        endif   
505                        NVAR delQ = $(folderPath+instPath+"ML"+":gDelQ_ML")
506                       
507                        Wave qTotal = $(folderPath+instPath+"ML"+":qTot_"+"ML")                 // 2D q-values 
508                        Wave qTotal2 = $(folderPath+instPath+"MR"+":qTot_"+"MR")                        // 2D q-values 
509                        Wave qTotal3 = $(folderPath+instPath+"MT"+":qTot_"+"MT")                        // 2D q-values 
510                        Wave qTotal4 = $(folderPath+instPath+"MB"+":qTot_"+"MB")                        // 2D q-values 
511
512                        Wave phi = V_MakePhiMatrix(qTotal,folderStr,"ML",folderPath+instPath+"ML")
513                        Wave phi2 = V_MakePhiMatrix(qTotal2,folderStr,"MR",folderPath+instPath+"MR")
514                        Wave phi3 = V_MakePhiMatrix(qTotal3,folderStr,"MT",folderPath+instPath+"MT")
515                        Wave phi4 = V_MakePhiMatrix(qTotal4,folderStr,"MB",folderPath+instPath+"MB")
516                                       
517                        nSets = 4
518                        break                                                                   
519                                       
520                default:
521                        nSets = 0                                                       
522                        Print "ERROR   ---- type is not recognized "
523        endswitch
524
525//      Print "delQ = ",delQ," for ",type
526
527        if(nSets == 0)
528                SetDataFolder root:
529                return(0)
530        endif
531
532
533// RAW data is currently read in and the 2D error wave is correctly generated
534// 2D error is propagated through all reduction steps, but I have not
535// verified that it is an exact duplication of the 1D error
536//
537//
538//
539// IF ther is no 2D error wave present for some reason, make a fake one
540        if(WaveExists(iErr)==0  && WaveExists(inten) != 0)
541                Duplicate/O inten,iErr
542                Wave iErr=iErr
543//              iErr = 1+sqrt(inten+0.75)                       // can't use this -- it applies to counts, not intensity (already a count rate...)
544                iErr = sqrt(inten+0.75)                 // TODO -- here I'm just using some fictional value
545        endif
546        if(WaveExists(iErr2)==0 && WaveExists(inten2) != 0)
547                Duplicate/O inten2,iErr2
548                Wave iErr2=iErr2
549//              iErr2 = 1+sqrt(inten2+0.75)                     // can't use this -- it applies to counts, not intensity (already a count rate...)
550                iErr2 = sqrt(inten2+0.75)                       // TODO -- here I'm just using some fictional value
551        endif
552        if(WaveExists(iErr3)==0  && WaveExists(inten3) != 0)
553                Duplicate/O inten3,iErr3
554                Wave iErr3=iErr3
555//              iErr3 = 1+sqrt(inten3+0.75)                     // can't use this -- it applies to counts, not intensity (already a count rate...)
556                iErr3 = sqrt(inten3+0.75)                       // TODO -- here I'm just using some fictional value
557        endif
558        if(WaveExists(iErr4)==0  && WaveExists(inten4) != 0)
559                Duplicate/O inten4,iErr4
560                Wave iErr4=iErr4
561//              iErr4 = 1+sqrt(inten4+0.75)                     // can't use this -- it applies to counts, not intensity (already a count rate...)
562                iErr4 = sqrt(inten4+0.75)                       // TODO -- here I'm just using some fictional value
563        endif
564
565        // TODO -- nq will need to be larger, once the back detector is installed
566        //
567        if(cmpstr(type,"B") == 0)
568                nq = 8000
569        else
570                nq=600
571        endif
572
573//******TODO****** -- where to put the averaged data -- right now, folderStr is forced to ""   
574//      SetDataFolder $("root:"+folderStr)              //should already be here, but make sure...     
575        Make/O/D/N=(nq)  $(folderPath+":"+"iBin_qxqy"+"_"+type)
576        Make/O/D/N=(nq)  $(folderPath+":"+"qBin_qxqy"+"_"+type)
577        Make/O/D/N=(nq)  $(folderPath+":"+"nBin_qxqy"+"_"+type)
578        Make/O/D/N=(nq)  $(folderPath+":"+"iBin2_qxqy"+"_"+type)
579        Make/O/D/N=(nq)  $(folderPath+":"+"eBin_qxqy"+"_"+type)
580        Make/O/D/N=(nq)  $(folderPath+":"+"eBin2D_qxqy"+"_"+type)
581       
582        Wave iBin_qxqy = $(folderPath+":"+"iBin_qxqy_"+type)
583        Wave qBin_qxqy = $(folderPath+":"+"qBin_qxqy"+"_"+type)
584        Wave nBin_qxqy = $(folderPath+":"+"nBin_qxqy"+"_"+type)
585        Wave iBin2_qxqy = $(folderPath+":"+"iBin2_qxqy"+"_"+type)
586        Wave eBin_qxqy = $(folderPath+":"+"eBin_qxqy"+"_"+type)
587        Wave eBin2D_qxqy = $(folderPath+":"+"eBin2D_qxqy"+"_"+type)
588       
589       
590//      delQ = abs(sqrt(qx[2]^2+qy[2]^2+qz[2]^2) - sqrt(qx[1]^2+qy[1]^2+qz[1]^2))               //use bins of 1 pixel width
591// TODO: not sure if I want to set dQ in x or y direction...
592        // the short dimension is the 8mm tubes, use this direction as dQ?
593        // but don't use the corner of the detector, since dQ will be very different on T/B or L/R due to the location of [0,0]
594        // WRT the beam center. use qx or qy directly. Still not happy with this way...
595
596
597        qBin_qxqy[] =  p*delQ   
598        SetScale/P x,0,delQ,"",qBin_qxqy                //allows easy binning
599
600        iBin_qxqy = 0
601        iBin2_qxqy = 0
602        eBin_qxqy = 0
603        eBin2D_qxqy = 0
604        nBin_qxqy = 0   //number of intensities added to each bin
605
606// now there are situations of:
607// 1 panel
608// 2 panels
609// 4 panels
610//
611// this needs to be a double loop now...
612// TODO:
613// -- the iErr (=2D) wave and accumulation of error is NOT CALCULATED CORRECTLY YET
614// -- verify the 2D error propagation by reducing it to 1D error
615//
616//
617// The 1D error does not use iErr, and IS CALCULATED CORRECTLY
618//
619// x- the solid angle per pixel will be present for WORK data other than RAW, but not for RAW
620
621//
622// if any of the masks don't exist, display the error, and proceed with the averaging, using all data
623        if(maskMissing == 1)
624                Print "Mask file not found for at least one detector - so all data is used"
625        endif
626       
627        NVAR gIgnoreDetB = root:Packages:NIST:VSANS:Globals:gIgnoreDetB
628        if(gIgnoreDetB && cmpstr(type,"B") == 0)
629                maskMissing = 1
630                Print "Mask skipped for B due to possible mismatch (Panel B ignored in preferences)"
631        endif
632
633        Variable mask_val,phiVal,isIn
634// use set 1 (no number) only
635        if(nSets >= 1)
636                xDim=DimSize(inten,0)
637                yDim=DimSize(inten,1)
638       
639                for(ii=0;ii<xDim;ii+=1)
640                        for(jj=0;jj<yDim;jj+=1)
641                                //qTot = sqrt(qx[ii]^2 + qy[ii]^2+ qz[ii]^2)
642                                qVal = qTotal[ii][jj]
643                                binIndex = trunc(x2pnt(qBin_qxqy, qVal))
644                                val = inten[ii][jj]
645                               
646                                if(isVCALC || maskMissing)              // mask_val == 0 == keep, mask_val == 1 = YES, mask out the point
647                                        mask_val = 0
648                                else
649                                        mask_val = mask[ii][jj]
650                                endif
651                               
652                                phiVal = phi[ii][jj]
653                                isIn = 0                        // start with exclude, now see if it's a keeper
654                                               
655                                // if within the right or left, flag to keep the pixel
656                                if(cmpstr(side,"right")==0)
657                                        //right, when 0->pi/2
658                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
659                                                isIn = 1
660                                        endif
661                                        // condition here to get the 3pi/2 -> 2pi region
662                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
663                                                isIn = 1
664                                        endif
665                                endif
666                               
667                                if(cmpstr(side,"left")==0)
668                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
669                                                isIn = 1
670                                        endif
671                                endif
672                                                       
673                        //      both sides, duplicates the conditions above
674                                if(cmpstr(side,"both")==0)     
675                                        //right, when 0->pi/2
676                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
677                                                isIn = 1
678                                        endif
679                                        // right, when 3pi/2 -> 2pi
680                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
681                                                isIn = 1
682                                        endif                           
683                                       
684                                        //left
685                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
686                                                isIn = 1
687                                        endif
688                                       
689                                endif           //end the check of phiVal within sector and side
690                               
691               
692                                if (numType(val)==0 && mask_val == 0 && isIn > 0)               //count only the good points, in the sector, and ignore Nan or Inf
693                                        iBin_qxqy[binIndex] += val
694                                        iBin2_qxqy[binIndex] += val*val
695                                        eBin2D_qxqy[binIndex] += iErr[ii][jj]*iErr[ii][jj]
696                                        nBin_qxqy[binIndex] += 1
697                                endif
698                               
699                               
700                        endfor
701                endfor
702               
703        endif
704
705// add in set 2 (set 1 already done)
706        if(nSets >= 2)
707                xDim=DimSize(inten2,0)
708                yDim=DimSize(inten2,1)
709       
710                for(ii=0;ii<xDim;ii+=1)
711                        for(jj=0;jj<yDim;jj+=1)
712                                //qTot = sqrt(qx[ii]^2 + qy[ii]^2+ qz[ii]^2)
713                                qVal = qTotal2[ii][jj]
714                                binIndex = trunc(x2pnt(qBin_qxqy, qVal))
715                                val = inten2[ii][jj]
716                               
717                                if(isVCALC || maskMissing)
718                                        mask_val = 0
719                                else
720                                        mask_val = mask2[ii][jj]
721                                endif
722                               
723                                phiVal = phi2[ii][jj]                   
724                                isIn = 0                        // start with exclude, now see if it's a keeper
725                                               
726                                // if within the right or left, flag to keep the pixel
727                                if(cmpstr(side,"right")==0)
728                                        //right, when 0->pi/2
729                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
730                                                isIn = 1
731                                        endif
732                                        // condition here to get the 3pi/2 -> 2pi region
733                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
734                                                isIn = 1
735                                        endif
736                                endif
737                               
738                                if(cmpstr(side,"left")==0)
739                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
740                                                isIn = 1
741                                        endif
742                                endif
743                                                       
744                        //      both sides, duplicates the conditions above
745                                if(cmpstr(side,"both")==0)     
746                                        //right, when 0->pi/2
747                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
748                                                isIn = 1
749                                        endif
750                                        // right, when 3pi/2 -> 2pi
751                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
752                                                isIn = 1
753                                        endif                           
754                                       
755                                        //left
756                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
757                                                isIn = 1
758                                        endif
759                                       
760                                endif           //end the check of phiVal within sector and side
761                               
762                               
763                                if (numType(val)==0 && mask_val == 0 && isIn > 0)               //count only the good points, ignore Nan or Inf
764                                        iBin_qxqy[binIndex] += val
765                                        iBin2_qxqy[binIndex] += val*val
766                                        eBin2D_qxqy[binIndex] += iErr2[ii][jj]*iErr2[ii][jj]
767                                        nBin_qxqy[binIndex] += 1
768                                endif
769                        endfor
770                endfor
771               
772        endif
773
774// add in set 3 and 4 (set 1 and 2 already done)
775        if(nSets == 4)
776                xDim=DimSize(inten3,0)
777                yDim=DimSize(inten3,1)
778       
779                for(ii=0;ii<xDim;ii+=1)
780                        for(jj=0;jj<yDim;jj+=1)
781                                //qTot = sqrt(qx[ii]^2 + qy[ii]^2+ qz[ii]^2)
782                                qVal = qTotal3[ii][jj]
783                                binIndex = trunc(x2pnt(qBin_qxqy, qVal))
784                                val = inten3[ii][jj]
785                               
786                                if(isVCALC || maskMissing)
787                                        mask_val = 0
788                                else
789                                        mask_val = mask3[ii][jj]
790                                endif
791                               
792                                phiVal = phi3[ii][jj]                   
793                                isIn = 0                        // start with exclude, now see if it's a keeper
794                                               
795                                // if within the right or left, flag to keep the pixel
796                                if(cmpstr(side,"right")==0)
797                                        //right, when 0->pi/2
798                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
799                                                isIn = 1
800                                        endif
801                                        // condition here to get the 3pi/2 -> 2pi region
802                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
803                                                isIn = 1
804                                        endif
805                                endif
806                               
807                                if(cmpstr(side,"left")==0)
808                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
809                                                isIn = 1
810                                        endif
811                                endif
812                                                       
813                        //      both sides, duplicates the conditions above
814                                if(cmpstr(side,"both")==0)     
815                                        //right, when 0->pi/2
816                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
817                                                isIn = 1
818                                        endif
819                                        // right, when 3pi/2 -> 2pi
820                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
821                                                isIn = 1
822                                        endif                           
823                                       
824                                        //left
825                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
826                                                isIn = 1
827                                        endif
828                                       
829                                endif           //end the check of phiVal within sector and side
830                               
831                                if (numType(val)==0 && mask_val == 0 && isIn > 0)               //count only the good points, ignore Nan or Inf
832                                        iBin_qxqy[binIndex] += val
833                                        iBin2_qxqy[binIndex] += val*val
834                                        eBin2D_qxqy[binIndex] += iErr3[ii][jj]*iErr3[ii][jj]
835                                        nBin_qxqy[binIndex] += 1
836                                endif
837                        endfor
838                endfor
839               
840               
841                xDim=DimSize(inten4,0)
842                yDim=DimSize(inten4,1)
843       
844                for(ii=0;ii<xDim;ii+=1)
845                        for(jj=0;jj<yDim;jj+=1)
846                                //qTot = sqrt(qx[ii]^2 + qy[ii]^2+ qz[ii]^2)
847                                qVal = qTotal4[ii][jj]
848                                binIndex = trunc(x2pnt(qBin_qxqy, qVal))
849                                val = inten4[ii][jj]
850                               
851                                if(isVCALC || maskMissing)
852                                        mask_val = 0
853                                else
854                                        mask_val = mask4[ii][jj]
855                                endif
856                               
857                                phiVal = phi4[ii][jj]                   
858                                isIn = 0                        // start with exclude, now see if it's a keeper
859                                               
860                                // if within the right or left, flag to keep the pixel
861                                if(cmpstr(side,"right")==0)
862                                        //right, when 0->pi/2
863                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
864                                                isIn = 1
865                                        endif
866                                        // condition here to get the 3pi/2 -> 2pi region
867                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
868                                                isIn = 1
869                                        endif
870                                endif
871                               
872                                if(cmpstr(side,"left")==0)
873                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
874                                                isIn = 1
875                                        endif
876                                endif
877                                                       
878                        //      both sides, duplicates the conditions above
879                                if(cmpstr(side,"both")==0)     
880                                        //right, when 0->pi/2
881                                        if(V_CloseEnough(phiVal,phi_rad,dphi_rad))
882                                                isIn = 1
883                                        endif
884                                        // right, when 3pi/2 -> 2pi
885                                        if(V_CloseEnough(phiVal,phi_rad+2*pi,dphi_rad))
886                                                isIn = 1
887                                        endif                           
888                                       
889                                        //left
890                                        if(V_CloseEnough(phiVal,phi_rad+pi,dphi_rad))
891                                                isIn = 1
892                                        endif
893                                       
894                                endif           //end the check of phiVal within sector and side
895                               
896                                if (numType(val)==0 && mask_val == 0 && isIn > 0)               //count only the good points, ignore Nan or Inf
897                                        iBin_qxqy[binIndex] += val
898                                        iBin2_qxqy[binIndex] += val*val
899                                        eBin2D_qxqy[binIndex] += iErr4[ii][jj]*iErr4[ii][jj]
900                                        nBin_qxqy[binIndex] += 1
901                                endif
902                        endfor
903                endfor
904               
905        endif
906
907
908// after looping through all of the data on the panels, calculate errors on I(q),
909// just like in CircSectAve.ipf
910// TODO:
911// -- 2D Errors were (maybe) properly acculumated through reduction, so this loop of calculations is NOT VERIFIED (yet)
912// x- the error on the 1D intensity, is correctly calculated as the standard error of the mean.
913        for(ii=0;ii<nq;ii+=1)
914                if(nBin_qxqy[ii] == 0)
915                        //no pixels in annuli, data unknown
916                        iBin_qxqy[ii] = 0
917                        eBin_qxqy[ii] = 1
918                        eBin2D_qxqy[ii] = NaN
919                else
920                        if(nBin_qxqy[ii] <= 1)
921                                //need more than one pixel to determine error
922                                iBin_qxqy[ii] /= nBin_qxqy[ii]
923                                eBin_qxqy[ii] = 1
924                                eBin2D_qxqy[ii] /= (nBin_qxqy[ii])^2
925                        else
926                                //assume that the intensity in each pixel in annuli is normally distributed about mean...
927                                //  -- this is correctly calculating the error as the standard error of the mean, as
928                                //    was always done for SANS as well.
929                                iBin_qxqy[ii] /= nBin_qxqy[ii]
930                                avesq = iBin_qxqy[ii]^2
931                                aveisq = iBin2_qxqy[ii]/nBin_qxqy[ii]
932                                var = aveisq-avesq
933                                if(var<=0)
934                                        eBin_qxqy[ii] = 1e-6
935                                else
936                                        eBin_qxqy[ii] = sqrt(var/(nBin_qxqy[ii] - 1))
937                                endif
938                                // and calculate as it is propagated pixel-by-pixel
939                                eBin2D_qxqy[ii] /= (nBin_qxqy[ii])^2
940                        endif
941                endif
942        endfor
943       
944        eBin2D_qxqy = sqrt(eBin2D_qxqy)         // as equation (3) of John's memo
945       
946        // find the last non-zero point, working backwards
947        val=nq
948        do
949                val -= 1
950        while((nBin_qxqy[val] == 0) && val > 0)
951       
952//      print val, nBin_qxqy[val]
953        DeletePoints val, nq-val, iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy,eBin2D_qxqy
954
955        if(val == 0)
956                // all the points were deleted, make dummy waves for resolution
957                Make/O/D/N=0  $(folderPath+":"+"sigmaQ_qxqy"+"_"+type)
958                Make/O/D/N=0  $(folderPath+":"+"qBar_qxqy"+"_"+type)
959                Make/O/D/N=0  $(folderPath+":"+"fSubS_qxqy"+"_"+type)
960                return(0)
961        endif
962       
963       
964        // since the beam center is not always on the detector, many of the low Q bins will have zero pixels
965        // find the first non-zero point, working forwards
966        val = -1
967        do
968                val += 1
969        while(nBin_qxqy[val] == 0)     
970        DeletePoints 0, val, iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy,eBin2D_qxqy
971
972        // ?? there still may be a point in the q-range that gets zero pixel contribution - so search this out and get rid of it
973        val = numpnts(nBin_qxqy)-1
974        do
975                if(nBin_qxqy[val] == 0)
976                        DeletePoints val, 1, iBin_qxqy,qBin_qxqy,nBin_qxqy,iBin2_qxqy,eBin_qxqy,eBin2D_qxqy
977                endif
978                val -= 1
979        while(val>0)
980
981// utility function to remove NaN values from the waves
982        V_RemoveNaNsQIS(qBin_qxqy, iBin_qxqy, eBin_qxqy)
983
984       
985        // TODO:
986        // -- This is where I calculate the resolution in SANS (see CircSectAve)
987        // -- use the isVCALC flag to exclude VCALC from the resolution calculation if necessary
988        // -- from the top of the function, folderStr = work folder, type = "FLRTB" or other type of averaging
989        //
990        nq = numpnts(qBin_qxqy)
991        Make/O/D/N=(nq)  $(folderPath+":"+"sigmaQ_qxqy"+"_"+type)
992        Make/O/D/N=(nq)  $(folderPath+":"+"qBar_qxqy"+"_"+type)
993        Make/O/D/N=(nq)  $(folderPath+":"+"fSubS_qxqy"+"_"+type)
994        Wave sigmaq = $(folderPath+":"+"sigmaQ_qxqy_"+type)
995        Wave qbar = $(folderPath+":"+"qBar_qxqy_"+type)
996        Wave fsubs = $(folderPath+":"+"fSubS_qxqy_"+type)
997        Variable ret1,ret2,ret3                 
998
999// all of the different collimation conditions are handled within the V_getResolution function
1000// which is responsible for switching based on the different collimation types (white beam, slit, Xtal, etc)
1001// to calculate the correct resolution, or fill the waves with the correct "flags"
1002//
1003
1004// For white beam data, the wavelength distribution can't be represented as a gaussian, but all of the other
1005//  geometric corrections still apply. Passing zero for the lambdaWidth will return the geometry contribution,
1006//  as long as the wavelength can be handled separately. It appears to be correct to do as a double integral,
1007//  with the inner(lambda) calculated first, then the outer(geometry).
1008//
1009
1010// possible values are:
1011//
1012// pinhole
1013// pinhole_whiteBeam
1014// convergingPinholes
1015//
1016// *slit data should be reduced using the slit routine, not here, proceed but warn
1017// narrowSlit
1018// narrowSlit_whiteBeam
1019//
1020        ii=0
1021        do
1022                V_getResolution(qBin_qxqy[ii],folderStr,type,collimationStr,ret1,ret2,ret3)
1023                sigmaq[ii] = ret1       
1024                qbar[ii] = ret2
1025                fsubs[ii] = ret3       
1026                ii+=1
1027        while(ii<nq)
1028
1029
1030        SetDataFolder root:
1031       
1032        return(0)
1033End
1034
Note: See TracBrowser for help on using the repository browser.