source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_WorkFolderUtils.ipf @ 1022

Last change on this file since 1022 was 1022, checked in by srkline, 6 years ago

ADDED:

included common NCNR procedures for the PlotManager?, to allow plotting of 1D data sets using a familiar interface

greatly expanded Patch functionality to include input in the multiple sections of the Nexus file, including separate panels to handle patching of waves to the file - needed for non-linear coefficients, dead time, and XY beam centers. All patch operations are expandable as more fields become necessary to patch.

removed bug of group_id being defined in /reduction and in /sample (removed R/W that referenced /reduction)

added panel to "isolate" a single detector panel, allowing the corrections to be applied/removed/recalculated as needed to directly see their effects.

linked new procedures to their appropriate action buttons

Added more data fields (label, intent, etc.) to the VCALC to Nexus data writer to get more realistic values into the fake data files for testing

Added VCALC simulation functions with EMP and BGD in anticipation of testing the CORRECT step

more little bug and documentation fixes which I can't remember, but they are all important...

File size: 30.1 KB
Line 
1#pragma rtGlobals=3             // Use modern global access method and strict wave access.
2#pragma version=1.0
3#pragma IgorVersion=6.1
4
5//*******************
6// Vers 1.0 JAN2016
7//
8//*******************
9//  VSANS Utility procedures for handling of workfiles (each is housed in a separate datafolder)
10//
11// - adding RAW data to a workfile
12// -- **this conversion applies the detector corrections**
13// -- Raw_to_work(newType) IS THE MAJOR ROUTINE TO APPLY DETECTOR CORRECTIONS
14//
15//
16// - copying workfiles to another folder
17//
18// - absolute scaling
19//
20// - (no) the WorkFile Math panel for simple image math (not done - maybe in the future?)
21// -
22// - (no) adding work.drk data without normalizing to monitor counts (the case not currently handled)
23//***************************
24
25//
26// Functions used for manipulation of the local Igor "WORK" folder
27// structure as raw data is displayed and processed.
28//
29//
30Strconstant ksDetectorListNoB = "FL;FR;FT;FB;ML;MR;MT;MB;"
31Strconstant ksDetectorListAll = "FL;FR;FT;FB;ML;MR;MT;MB;B;"
32
33
34//
35//Entry procedure from main panel
36//
37Proc CopyWorkFolder(oldType,newType)
38        String oldType,newType
39        Prompt oldType,"Source WORK data type",popup,"SAM;EMP;BGD;DIV;COR;CAL;RAW;ABS;STO;SUB;DRK;"
40        Prompt newType,"Destination WORK data type",popup,"SAM;EMP;BGD;DIV;COR;CAL;RAW;ABS;STO;SUB;DRK;"
41
42        // data folder "old" will be copied to "new" (either kills/copies or will overwrite)
43        CopyHDFToWorkFolder(oldtype,newtype)
44End
45
46//
47// copy what is needed for data processing (not the DAS_logs)
48// from the RawVSANS storage folder to the local WORK folder as needed
49//
50// TODO -- at what stage do I make copies of data in linear/log forms for data display?
51//      -- do I need to make copies, if I'm displaying with the lookup wave (no copy needed if this works)
52//                      -- when do I make the 2D error waves?
53//
54// TODO - decide what exactly I need to copy over. May be best to copy all, and delete
55//       what I know that I don't need
56//
57// TODO !!! DuplicateDataFolder will FAIL - in the base case of RAW data files, the
58//  data is actually in use - so it will fail every time. need an alternate solution. in SANS,
59// there are a limited number of waves to carry over, so Dupliate/O is used for rw, tw, data, etc.
60//
61// TODO : I also need a list of what is generated during processing that may be hanging around - that I need to
62//     be sure to get rid of - like the calibration waves, solidAngle, etc.
63//
64// hdfDF is the name only of the data in storage. May be full file name with extension (clean as needed)
65// type is the destination WORK folder for the copy
66//
67Function CopyHDFToWorkFolder(fromStr,toStr)
68        String fromStr,toStr
69       
70        String fromDF, toDF
71       
72        // make the DF paths - source and destination
73        fromDF = "root:Packages:NIST:VSANS:"+fromStr
74        toDF = "root:Packages:NIST:VSANS:"+toStr
75       
76//      // make a copy of the file name for my own use, since it's not in the file
77//      String/G $(toDF+":file_name") = root:
78       
79        // copy the folders
80        KillDataFolder/Z $toDF                  //DuplicateDataFolder will not overwrite, so Kill
81       
82        if(V_flag == 0)         // kill DF was OK
83                DuplicateDataFolder $("root:Packages:NIST:VSANS:"+fromStr),$("root:Packages:NIST:VSANS:"+toStr)
84               
85                // I can delete these if they came along with RAW
86                //   DAS_logs
87                //   top-level copies of data (duplicate links)
88                KillDataFolder/Z $(toDF+":entry:DAS_logs")
89                KillDataFolder/Z $(toDF+":entry:data")
90                KillDataFolder/Z $(toDF+":entry:data_B")
91                KillDataFolder/Z $(toDF+":entry:data_ML")
92                KillDataFolder/Z $(toDF+":entry:data_MR")
93                KillDataFolder/Z $(toDF+":entry:data_MT")
94                KillDataFolder/Z $(toDF+":entry:data_MB")
95                KillDataFolder/Z $(toDF+":entry:data_FL")
96                KillDataFolder/Z $(toDF+":entry:data_FR")
97                KillDataFolder/Z $(toDF+":entry:data_FT")
98                KillDataFolder/Z $(toDF+":entry:data_FB")
99
100                return(0)
101        else
102                // need to do this the hard way, duplicate/O recursively
103                // see V_CopyToWorkFolder()
104               
105                // everything on the top level
106                V_DuplicateDataFolder($(fromDF+":entry"),fromStr,toStr,0,"",0)  //no recursion here
107                // control
108                V_DuplicateDataFolder($(fromDF+":entry:control"),fromStr,toStr,0,"",1)  //yes recursion here
109                // instrument
110                V_DuplicateDataFolder($(fromDF+":entry:instrument"),fromStr,toStr,0,"",1)       //yes recursion here
111                // reduction
112                V_DuplicateDataFolder($(fromDF+":entry:reduction"),fromStr,toStr,0,"",1)        //yes recursion here
113                // sample
114                V_DuplicateDataFolder($(fromDF+":entry:sample"),fromStr,toStr,0,"",1)   //yes recursion here
115
116        endif   
117       
118        return(0)
119end
120
121
122////////
123// see the help entry for IndexedDir for help on (possibly) how to do this faster
124// -- see  Function ScanDirectories(pathName, printDirNames)
125//
126
127
128// from IgorExchange On July 17th, 2011 jtigor
129// started from "Recursively List Data Folder Contents"
130// Posted July 15th, 2011 by hrodstein
131//
132//
133//
134Proc V_CopyWorkFolderProc(dataFolderStr, fromStr, toStr, level, sNBName, recurse)
135        String dataFolderStr="root:Packages:NIST:VSANS:RAW"
136        String fromStr = "RAW"
137        String toStr="SAM"
138        Variable level=0
139        String sNBName="DataFolderTree"
140        Variable recurse = 1
141       
142        V_DuplicateDataFolder($dataFolderStr, fromStr, toStr, level, sNBName, recurse)
143
144
145end
146
147// ListDataFolder(dfr, level)
148// Recursively lists objects in data folder.
149// Pass data folder path for dfr and 0 for level.
150// pass level == 0 for the first call
151//  sNBName = "" prints nothing. any name will generate a notebook
152//
153// recurse == 0 will do only the specified folder, anything else will recurse all levels
154// toStr is the string name of the top-level folder only, not the full path
155//
156//
157Function V_DuplicateDataFolder(dfr, fromStr, toStr, level, sNBName,recurse)
158        DFREF dfr
159        String fromStr
160        String toStr
161        Variable level                  // Pass 0 to start
162        String sNBName
163        Variable recurse
164 
165        String name
166        String dfName
167        String sString
168       
169        String toDF = ""
170 
171        if (level == 0)         // this is the data folder, generate if needed in the destination
172                name = GetDataFolder(1, dfr)
173//              sPrintf sString, "%s (data folder)\r", name
174                toDF = ReplaceString(fromStr,name,toStr,1)              // case-sensitive replace
175                sprintf sString, "NewDataFolder/O %s\r",toDF
176                NewDataFolder/O $(RemoveEnding(toDF,":"))                       // remove trailing semicolon if it's there
177               
178                V_WriteBrowserInfo_test(sString, 1, sNBName)
179        endif
180 
181        dfName = GetDataFolder(1, dfr)
182        toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
183        Variable i
184 
185        String indentStr = "\t"
186        for(i=0; i<level; i+=1)
187                indentStr += "\t"
188        endfor
189 
190        Variable numWaves = CountObjectsDFR(dfr, 1)
191        for(i=0; i<numWaves; i+=1)
192                name = GetIndexedObjNameDFR(dfr, 1, i)
193                //
194                // wave type does not matter now. Duplicate does not care
195                //
196                sPrintf sString, "Duplicate/O  %s,  %s\r",dfName+name,toDF+name
197                Duplicate/O $(dfName+name),$(toDF+name)
198               
199                V_WriteBrowserInfo_test(sString, 2, sNBName)
200        endfor 
201 
202        Variable numNumericVariables = CountObjectsDFR(dfr, 2) 
203        for(i=0; i<numNumericVariables; i+=1)
204                name = GetIndexedObjNameDFR(dfr, 2, i)
205                sPrintf sString, "%s%s (numeric variable)\r", indentStr, name
206                V_WriteBrowserInfo_test(sString, 3, sNBName)
207        endfor 
208 
209        Variable numStringVariables = CountObjectsDFR(dfr, 3)   
210        for(i=0; i<numStringVariables; i+=1)
211                name = GetIndexedObjNameDFR(dfr, 3, i)
212                sPrintf sString, "%s%s (string variable)\r", indentStr, name
213                V_WriteBrowserInfo_test(sString, 4, sNBName)
214        endfor 
215
216        if(recurse)
217                Variable numDataFolders = CountObjectsDFR(dfr, 4)       
218                for(i=0; i<numDataFolders; i+=1)
219                        name = GetIndexedObjNameDFR(dfr, 4, i)
220//                      sPrintf sString, "%s%s (data folder)\r", indentStr, name
221                         dfName = GetDataFolder(1, dfr)
222                         
223                        toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
224                        sprintf sString, "NewDataFolder/O %s\r",toDF+name
225                        NewDataFolder/O $(toDF+name)
226                       
227                       
228                        V_WriteBrowserInfo_test(sString, 1, sNBName)
229                        DFREF childDFR = dfr:$(name)
230                        V_DuplicateDataFolder(childDFR, fromStr, toStr, level+1, sNBName, recurse)
231                endfor 
232        endif
233         
234
235End
236 
237Function V_WriteBrowserInfo_test(sString, vType, sNBName)
238        String sString
239        Variable vType
240        String sNBName
241 
242        if(strlen(sNBName) == 0)
243//              print sString
244                return 0
245        endif
246        DoWindow $sNBName
247        if(V_flag != 1)
248                NewNoteBook/F=0 /N=$sNBName /V=1 as sNBName
249        else
250                DoWindow/F $sNBName
251        endif
252        Notebook $sNBName selection={endOfFile, endOfFile}
253        if(vType == 1)          // a data folder
254//              Notebook $sNBName fstyle=1
255                Notebook $sNBName text=sString
256//              Notebook $sNBName fstyle=-1
257        else
258                Notebook $sNBName text=sString 
259        endif
260 
261End
262
263///////////////////////////////
264
265
266//
267// given the folder, duplicate the data -> linear_data and generate the error
268// TODO
269// -- do I want to use different names here? If it turns out that I don't need to drag a copy of
270//    the data around as "linear_data", then I can eliminate that, and rename the error wave
271// -- be sure the data is either properly written as 2D in the file, or converted to 2D before
272//    duplicating here
273// -- ? do I recast to DP here. Probably necessary since I'm doing a DP calculation, but Redimension
274//    is done in the Raw_to_Work step too. very confusing.
275Function V_MakeDataError(folderStr)
276        String folderStr
277       
278        SetDataFolder $folderStr
279        Wave data=data
280        Duplicate/O data linear_data                    // at this point, the data is still the raw data, and is linear_data
281       
282        // proper error for counting statistics, good for low count values too
283        // rather than just sqrt(n)
284        // see N. Gehrels, Astrophys. J., 303 (1986) 336-346, equation (7)
285        // for S = 1 in eq (7), this corresponds to one sigma error bars
286        Duplicate/O linear_data linear_data_error
287        linear_data_error = 1 + sqrt(linear_data + 0.75)                               
288        //
289       
290        SetDataFolder root:
291        return(0)
292End
293
294
295/////////////////////
296
297
298
299//testing procedure
300// TODO -- can't duplicate this with another proceudre, but if I change the name of the variable
301//   "newType" to "type", then when Raw_to_work() gets to CopyHDFToWorkFolder(), the KillDataFolder/Z
302//   line fails (but reports no error), then DuplicateDataFolder fails, and reports an error. Trying
303//   to simplify this condition, I can't reproduce the bug for WM...
304Proc Convert_to_Workfile(newtype, doadd)
305        String newtype,doadd
306        Prompt newtype,"WORK data type",popup,"SAM;EMP;BGD;ADJ;"
307        Prompt doadd,"Add to current WORK contents?",popup,"No;Yes;"
308       
309        //macro will take whatever is in RAW folder and "ADD" it to the folder specified
310        //in the popup menu
311       
312        //"add" = yes/no, don't add to previous runs
313        //switch here - two separate functions to avoid (my) confusion
314        Variable err// = Raw_to_work(newtype)
315        if(cmpstr(doadd,"No")==0)
316                //don't add to prev work contents, copy RAW contents to work and convert
317                err = Raw_to_work(newtype)
318        else
319                //yes, add RAW to the current work folder contents
320                Abort "Adding RAW data files is currently unsupported"
321                err = Add_raw_to_work(newtype)
322        endif
323       
324        String newTitle = "WORK_"+newtype
325        DoWindow/F VSANS_Data
326        DoWindow/T VSANS_Data, newTitle
327        KillStrings/Z newTitle
328       
329        //need to update the display with "data" from the correct dataFolder
330        UpdateDisplayInformation(newtype)
331       
332End
333
334
335//
336// THIS IS THE MAJOR ROUTINE TO APPLY DATA CORRECTIONS
337//
338//will copy the current contents of the RAW folder to the newType work folder
339//and do the geometric corrections and normalization to monitor counts
340//(the function Add_Raw_to_work(type) adds multiple runs together - and is LOW priority)
341//
342//the current display type is updated to newType (global)
343//
344Function Raw_to_work(newType)
345        String newType
346       
347        Variable deadTime,defmon,total_mon,total_det,total_trn,total_numruns,total_rtime
348        Variable ii,jj,itim,cntrate,dscale,scale,uscale
349        String destPath
350       
351        String fname = newType
352        String detStr
353        Variable ctTime
354
355        //initialize values before normalization
356        total_mon=0
357        total_det=0
358        total_trn=0
359        total_numruns=0
360        total_rtime=0
361       
362        //Not adding multiple runs, so wipe out the old contents of the work folder and
363        // replace with the contents of raw
364
365        destPath = "root:Packages:NIST:VSANS:" + newType
366       
367        //copy from current dir (RAW) to work, defined by newType
368        CopyHDFToWorkFolder("RAW",newType)
369       
370        // now work with the waves from the destination folder.
371       
372        // apply corrections ---
373        // switches to control what is done, don't do the transmission correction for the BGD measurement
374        // start with the DIV correction, before conversion to mm
375        // then do all of the other corrections, order doesn't matter.
376        // rescaling to default monitor counts however, must be LAST.
377
378// each correction must loop over each detector. tedious.
379
380        // (0) Redimension the data waves in the destination folder
381        //     so that they are DP, not integer
382        // TODO
383        // -- currently only redimensioning the data and linear_data_error - What else???
384        // -- ?? some of this is done at load time for RAW data. shouldn't be an issue to re-do the redimension
385        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
386                detStr = StringFromList(ii, ksDetectorListAll, ";")
387                Wave w = V_getDetectorDataW(fname,detStr)
388                Wave w_err = V_getDetectorDataErrW(fname,detStr)
389                Redimension/D w,w_err
390        endfor
391       
392       
393        // (1) DIV correction
394        // do this in terms of pixels.
395        // TODO : This must also exist at the time the first work folder is generated.
396        //   So it must be in the user folder at the start of the experiment, and defined.
397        NVAR gDoDIVCor = root:Packages:NIST:VSANS:Globals:gDoDIVCor
398        if (gDoDIVCor == 1)
399                // need extra check here for file existence
400                // if not in DIV folder, load.
401                // if unable to load, skip correction and report error (Alert?) (Ask to Load?)
402                Print "Doing DIV correction"// for "+ detStr
403                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
404                        detStr = StringFromList(ii, ksDetectorListAll, ";")
405                        Wave w = V_getDetectorDataW(fname,detStr)
406                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
407                       
408                        V_DIVCorrection(w,w_err,detStr,newType)
409                endfor
410        else
411                Print "DIV correction NOT DONE"         // not an error since correction was unchecked
412        endif
413       
414        // (2) non-linear correction   
415        // TODO:
416        // x-  the "B" detector is calculated in its own routines
417        // -- document what is generated here:
418        //    **in each detector folder: data_realDistX and data_realDistY (2D waves of the mm? position of the pixel)
419        // x- these spatial calculations ARE DONE as the RAW data is loaded. It allows the RAW
420        //    data to be properly displayed, but without all of the (complete) set of detector corrections
421        // * the corrected distances are calculated into arrays, but nothing is done with them yet
422        // * there is enough information now to calculate the q-arrays, so it is done now
423        // - other corrections may modify the data, this calculation does NOT modify the data
424        NVAR gDoNonLinearCor = root:Packages:NIST:VSANS:Globals:gDoNonLinearCor
425        // generate a distance matrix for each of the detectors
426        if (gDoNonLinearCor == 1)
427                Print "Doing Non-linear correction"// for "+ detStr
428                for(ii=0;ii<ItemsInList(ksDetectorListNoB);ii+=1)
429                        detStr = StringFromList(ii, ksDetectorListNoB, ";")
430                        Wave w = V_getDetectorDataW(fname,detStr)
431//                      Wave w_err = V_getDetectorDataErrW(fname,detStr)
432                        Wave w_calib = V_getDetTube_spatialCalib(fname,detStr)
433                        Variable tube_width = V_getDet_tubeWidth(fname,detStr)
434                        NonLinearCorrection(w,w_calib,tube_width,detStr,destPath)
435
436                        //(2.4) Convert the beam center values from pixels to mm
437                        // TODO -- there needs to be a permanent location for these values??
438                        //
439                        ConvertBeamCtr_to_mm(fname,detStr,destPath)
440                                                       
441                        // (2.5) Calculate the q-values
442                        // calculating q-values can't be done unless the non-linear corrections are calculated
443                        // so go ahead and put it in this loop.
444                        // TODO :
445                        // -- make sure that everything is present before the calculation
446                        // -- beam center must be properly defined in terms of real distance
447                        // -- distances/zero location/ etc. must be clearly documented for each detector
448                        //      ** this assumes that NonLinearCorrection() has been run to generate data_RealDistX and Y
449                        // ** this routine Makes the waves QTot, qx, qy, qz in each detector folder.
450                        //
451                        V_Detector_CalcQVals(fname,detStr,destPath)
452                       
453                endfor
454               
455                //"B" is separate
456                NonLinearCorrection_B(fname,"B",destPath)
457                ConvertBeamCtr_to_mmB(fname,"B",destPath)
458                V_Detector_CalcQVals(fname,"B",destPath)
459               
460        else
461                Print "Non-linear correction NOT DONE"
462        endif
463
464        // (3) solid angle correction
465        // TODO -- this currently calculates the correction factor AND applys it to the data
466        //  -- as a result, the data values are very large since they are divided by a very small
467        //     solid angle per pixel. But all of the count values are now on the basis of
468        //    counts/(solid angle) --- meaning that they can all be binned together for I(q)
469        //    -and- TODO - this will need to be taken into account for absolute scaling (this part is already done)
470        //
471        // the solid angle correction is calculated for ALL detector panels.
472        NVAR gDoSolidAngleCor = root:Packages:NIST:VSANS:Globals:gDoSolidAngleCor
473        if (gDoSolidAngleCor == 1)
474                Print "Doing Solid Angle correction"// for "+ detStr
475                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
476                        detStr = StringFromList(ii, ksDetectorListAll, ";")
477                        Wave w = V_getDetectorDataW(fname,detStr)
478                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
479                        // any other dimensions to pass in?
480                        SolidAngleCorrection(w,w_err,fname,detStr,destPath)
481                       
482                endfor
483        else
484                Print "Solid Angle correction NOT DONE"
485        endif   
486       
487        // (4) dead time correction
488        // TODO:
489        // x- remove the hard-wired test - done
490        // -- test for correct operation
491        // x- loop over all of the detectors
492        // x- B detector is a special case (do separately, then loop over NoB)
493        // -- this DOES alter the data
494        // -- verify the error propagation (not done yet)
495        //
496        Variable countRate
497        NVAR gDoDeadTimeCor = root:Packages:NIST:VSANS:Globals:gDoDeadTimeCor
498        if (gDoDeadTimeCor == 1)
499                Print "Doing DeadTime correction"// for "+ detStr
500                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
501                        detStr = StringFromList(ii, ksDetectorListAll, ";")
502                        Wave w = V_getDetectorDataW(fname,detStr)
503                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
504                        ctTime = V_getCount_time(fname)
505
506                        if(cmpstr(detStr,"B") == 0)
507                                Variable b_dt = V_getDetector_deadtime_B(fname,detStr)
508                                // do the correction for the back panel
509                                countRate = sum(w,-inf,inf)/ctTime              //use sum of detector counts
510
511                                w = w/(1-countRate*b_dt)
512                                w_err = w_err/(1-countRate*b_dt)
513                                                               
514                        else
515                                // do the corrections for 8 tube panels
516                                Wave w_dt = V_getDetector_deadtime(fname,detStr)
517                                DeadTimeCorrectionTubes(w,w_err,w_dt,ctTime)
518                               
519                        endif
520                endfor
521               
522        else
523                Print "Dead Time correction NOT DONE"
524        endif   
525       
526       
527        // (5) angle-dependent tube shadowing
528        // TODO:
529        // -- not sure about this correction yet...
530        //
531        NVAR gDoTubeShadowCor = root:Packages:NIST:VSANS:Globals:gDoTubeShadowCor
532        if (gDoTubeShadowCor == 1)
533                Print "(stub)Tube Shadow correction"
534        else
535                Print "Tube shadowing correction NOT DONE"
536        endif   
537               
538        // (6) angle dependent transmission correction
539        // TODO:
540        // -- this still needs to be filled in
541        // -- still some debate of when/where in the corrections that this is best applied
542        //    - do it here, and it's done whether the output is 1D or 2D
543        //    - do it later (where SAMPLE information is used) since this section is ONLY instrument-specific
544        NVAR gDoTrans = root:Packages:NIST:VSANS:Globals:gDoTransmissionCor
545        if (gDoTrans == 1)
546                Print "(stub)Doing Large-angle transmission correction -- INCOMPLETE FUNCTION "// for "+ detStr
547                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
548                        detStr = StringFromList(ii, ksDetectorListAll, ";")
549                        Wave w = V_getDetectorDataW(fname,detStr)
550                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
551                       
552//                      TransmissionCorrection(fill this in)
553                       
554                endfor
555        else
556                Print "Sample Transmission correction NOT DONE"
557        endif   
558       
559        // (7) normalize to default monitor counts
560        // TODO -- each detector is rescaled separately, but the rescaling factor is global (only one monitor!)
561        // TODO -- but there are TWO monitors - so how to switch?
562        // TODO -- what do I really need to save?
563       
564        NVAR gDoMonitorNormalization = root:Packages:NIST:VSANS:Globals:gDoMonitorNormalization
565        if (gDoMonitorNormalization == 1)
566                Print "(stub)Doing monitor normalization -- INCOMPLETE --"// for "+ detStr
567               
568                defmon=1e8                      //default monitor counts
569                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
570                        detStr = StringFromList(ii, ksDetectorListAll, ";")
571                        Wave w = V_getDetectorDataW(fname,detStr)
572                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
573                        Variable monCt = V_getBeamMonNormData(fname)
574        //                      MonitorNormalization(fill this in)
575                //scale the data to the default montor counts
576               
577                // TODO -- un-comment these three lines once monitor counts are reasonable - currently monCt = 9!!!
578        //              scale = defmon/monCt
579        //              w *= scale
580        //              w_err *= scale          //assumes total monitor count is so large there is essentially no error
581       
582        // TODO
583        // -- to write back to the local value, get the wave reference rather than the value, then I can
584        //    re-assign the value directly, rather than this method (which is not terrible)     
585                        // V_getBeamMonNormSaved_count()
586                        // save the true monitor counts? save the scaling factor?
587                        String path = "entry:instrument:beam_monitor_norm:saved_count"
588                        Wave/Z savW = $("root:Packages:NIST:VSANS:"+fname+":"+path)
589                        savW[0] = scale
590                endfor
591        else
592                Print "Monitor Normalization correction NOT DONE"
593        endif
594       
595        // (not done) angle dependent efficiency correction
596        NVAR doEfficiency = root:Packages:NIST:VSANS:Globals:gDoDetectorEffCor
597
598//update totals to put in the work header (at the end of the function)
599//      total_mon += realsread[0]
600//
601//      total_det += dscale*realsread[2]
602//
603//      total_trn += realsread[39]
604//      total_rtime += integersread[2]
605//      total_numruns +=1
606//     
607
608        //all is done, except for the bookkeeping, updating the header information in the work folder
609
610//      integersread[3] = total_numruns                                         //numruns = 1
611//      realsread[1] = total_mon                        //save the true monitor count
612//      realsread[0] = defmon                                   //monitor ct = defmon
613//      realsread[2] = scale*total_det                  //scaled detector counts
614//     
615        //reset the current displaytype to "newtype"
616        String/G root:Packages:NIST:VSANS:Globals:gCurDispType=newType
617       
618        //return to root folder (redundant)
619        SetDataFolder root:
620       
621        Return(0)
622End
623
624
625//will "ADD" the current contents of the RAW folder to the newType work folder
626//and will ADD the RAW contents to the existing content of the newType folder
627// - used when adding multiple runs together
628//(the function Raw_to_work(type) makes a fresh workfile)
629//
630//the current display type is updated to newType (global)
631Function Add_raw_to_work(newType)
632        String newType
633       
634        // NEW OCT 2014
635        // this corrects for adding raw data files with different attenuation   
636        // does nothing if the attenuation of RAW and destination are the same
637        NVAR doAdjustRAW_Atten = root:Packages:NIST:gDoAdjustRAW_Atten
638        if(doAdjustRAW_Atten)
639                Adjust_RAW_Attenuation(newType)
640        endif
641       
642        String destPath=""
643       
644        // if the desired workfile doesn't exist, let the user know, and just make a new one
645        if(WaveExists($("root:Packages:NIST:"+newType + ":data")) == 0)
646                Print "There is no old work file to add to - a new one will be created"
647                //call Raw_to_work(), then return from this function
648                Raw_to_Work(newType)
649                Return(0)               //does not generate an error - a single file was converted to work.newtype
650        Endif
651       
652        NVAR pixelsX = root:myGlobals:gNPixelsX
653        NVAR pixelsY = root:myGlobals:gNPixelsY
654       
655        //now make references to data in newType folder
656        DestPath="root:Packages:NIST:"+newType 
657        WAVE data=$(destPath +":linear_data")                   // these wave references point to the EXISTING work data
658        WAVE data_copy=$(destPath +":data")                     // these wave references point to the EXISTING work data
659        WAVE dest_data_err=$(destPath +":linear_data_error")                    // these wave references point to the EXISTING work data
660        WAVE/T textread=$(destPath + ":textread")
661        WAVE integersread=$(destPath + ":integersread")
662        WAVE realsread=$(destPath + ":realsread")
663       
664        Variable deadTime,defmon,total_mon,total_det,total_trn,total_numruns,total_rtime
665        Variable ii,jj,itim,cntrate,dscale,scale,uscale,wrk_beamx,wrk_beamy,xshift,yshift
666
667
668        defmon=1e8                      //default monitor counts
669       
670        //Yes, add to previous run(s) in work, that does exist
671        //use the actual monitor count run.savmon rather than the normalized monitor count
672        //in run.moncnt and unscale the work data
673       
674        total_mon = realsread[1]        //saved monitor count
675        uscale = total_mon/defmon               //unscaling factor
676        total_det = uscale*realsread[2]         //unscaled detector count
677        total_trn = uscale*realsread[39]        //unscaled trans det count
678        total_numruns = integersread[3] //number of runs in workfile
679        total_rtime = integersread[2]           //total counting time in workfile
680        //retrieve workfile beamcenter
681        wrk_beamx = realsread[16]
682        wrk_beamy = realsread[17]
683        //unscale the workfile data in "newType"
684        //
685        //check for log-scaling and adjust if necessary
686        // should not be needed now - using display flag instead
687//      ConvertFolderToLinearScale(newType)
688        //
689        //then unscale the data array
690        data *= uscale
691        dest_data_err *= uscale
692       
693        //DetCorr() has not been applied to the data in RAW , do it now in a local reference to the raw data
694        WAVE raw_data = $"root:Packages:NIST:RAW:linear_data"
695        WAVE raw_data_err = $"root:Packages:NIST:RAW:linear_data_error"
696        WAVE raw_reals =  $"root:Packages:NIST:RAW:realsread"
697        WAVE/T raw_text = $"root:Packages:NIST:RAW:textread"
698        WAVE raw_ints = $"root:Packages:NIST:RAW:integersread"
699       
700        //check for log-scaling of the raw data - make sure it's linear
701        // should not be needed now - using display flag instead
702//      ConvertFolderToLinearScale("RAW")
703       
704        // switches to control what is done, don't do the transmission correction for the BGD measurement
705        NVAR doEfficiency = root:Packages:NIST:gDoDetectorEffCorr
706        NVAR gDoTrans = root:Packages:NIST:gDoTransmissionCorr
707        Variable doTrans = gDoTrans
708        if(cmpstr("BGD",newtype) == 0)
709                doTrans = 0             //skip the trans correction for the BGD file but don't change the value of the global
710        endif   
711       
712        DetCorr(raw_data,raw_data_err,raw_reals,doEfficiency,doTrans)   //applies correction to raw_data, and overwrites it
713       
714        //deadtime corrections to raw data
715        // TODO - do the tube correction for dead time now
716        deadTime = 1//DetectorDeadtime(raw_text[3],raw_text[9],dateAndTimeStr=raw_text[1],dtime=raw_reals[48])          //pick the correct detector deadtime, switch on date too
717        itim = raw_ints[2]
718        cntrate = sum(raw_data,-inf,inf)/itim           //080802 use data sum, rather than scaler value
719        dscale = 1/(1-deadTime*cntrate)
720
721#if (exists("ILL_D22")==6)
722        Variable tubeSum
723        // for D22 detector might need to use cntrate/128 as it is the tube response
724        for(ii=0;ii<pixelsX;ii+=1)
725                //sum the counts in each tube
726                tubeSum = 0
727                for(jj=0;jj<pixelsY;jj+=1)
728                        tubeSum += data[jj][ii]
729                endfor
730                // countrate in tube ii
731                cntrate = tubeSum/itim
732                // deadtime scaling in tube ii
733                dscale = 1/(1-deadTime*cntrate)
734                // multiply data[ii][] by the dead time
735                raw_data[][ii] *= dscale
736                raw_data_err[][ii] *= dscale
737        endfor
738#else
739        // dead time correction on all other RAW data, including NCNR
740        raw_data *= dscale
741        raw_data_err *= dscale
742#endif
743
744        //update totals by adding RAW values to the local ones (write to work header at end of function)
745        total_mon += raw_reals[0]
746
747        total_det += dscale*raw_reals[2]
748
749        total_trn += raw_reals[39]
750        total_rtime += raw_ints[2]
751        total_numruns +=1
752       
753        //do the beamcenter shifting if there is a mismatch
754        //and then add the two data sets together, changing "data" since it is the workfile data
755        xshift = raw_reals[16] - wrk_beamx
756        yshift = raw_reals[17] - wrk_beamy
757       
758        If((xshift != 0) || (yshift != 0))
759                DoAlert 1,"Do you want to ignore the beam center mismatch?"
760                if(V_flag==1)
761                        xshift=0
762                        yshift=0
763                endif
764        endif
765       
766        If((xshift == 0) && (yshift == 0))              //no shift, just add them
767                data += raw_data                //deadtime correction has already been done to the raw data
768                dest_data_err = sqrt(dest_data_err^2 + raw_data_err^2)                  // error of the sum
769        Endif
770       
771        //scale the data to the default montor counts
772        scale = defmon/total_mon
773        data *= scale
774        dest_data_err *= scale
775       
776        // keep "data" and linear_data in sync in the destination folder
777        data_copy = data
778       
779        //all is done, except for the bookkeeping of updating the header info in the work folder
780        textread[1] = date() + " " + time()             //date + time stamp
781        integersread[3] = total_numruns                                         //numruns = more than one
782        realsread[1] = total_mon                        //save the true monitor count
783        realsread[0] = defmon                                   //monitor ct = defmon
784        integersread[2] = total_rtime                   // total counting time
785        realsread[2] = scale*total_det                  //scaled detector counts
786        realsread[39] = scale*total_trn                 //scaled transmission counts
787       
788        //Add the added raw filename to the list of files in the workfile
789        String newfile = ";" + raw_text[0]
790        SVAR oldList = $(destPath + ":fileList")
791        String/G $(destPath + ":fileList") = oldList + newfile
792       
793        //reset the current display type to "newtype"
794        SVAR gCurDispType = root:Packages:NIST:VSANS:Globals:gCurDispType
795        gCurDispType = newType
796       
797        //return to root folder (redundant)
798        SetDataFolder root:
799       
800        Return(0)
801End
802
803
804//used for adding DRK (beam shutter CLOSED) data to a workfile
805//force the monitor count to 1, since it's irrelevant
806// run data through normal "add" step, then unscale default monitor counts
807//to get the data back on a simple time basis
808//
809Function Raw_to_Work_NoNorm(type)
810        String type
811       
812        WAVE reals=$("root:Packages:NIST:RAW:realsread")
813        reals[1]=1              //true monitor counts, still in raw
814        Raw_to_work(type)
815        //data is now in "type" folder
816        WAVE data=$("root:Packages:NIST:"+type+":linear_data")
817        WAVE data_copy=$("root:Packages:NIST:"+type+":data")
818        WAVE data_err=$("root:Packages:NIST:"+type+":linear_data_error")
819        WAVE new_reals=$("root:Packages:NIST:"+type+":realsread")
820       
821        Variable norm_mon,tot_mon,scale
822       
823        norm_mon = new_reals[0]         //should be 1e8
824        tot_mon = new_reals[1]          //should be 1
825        scale= norm_mon/tot_mon
826       
827        data /= scale           //unscale the data
828        data_err /= scale
829       
830        // to keep "data" and linear_data in sync
831        data_copy = data
832       
833        return(0)
834End
835
836//used for adding DRK (beam shutter CLOSED) data to a workfile
837//force the monitor count to 1, since it's irrelevant
838// run data through normal "add" step, then unscale default monitor counts
839//to get the data back on a simple time basis
840//
841Function Add_Raw_to_Work_NoNorm(type)
842        String type
843       
844        WAVE reals=$("root:Packages:NIST:RAW:realsread")
845        reals[1]=1              //true monitor counts, still in raw
846        Add_Raw_to_work(type)
847        //data is now in "type" folder
848        WAVE data=$("root:Packages:NIST:"+type+":linear_data")
849        WAVE data_copy=$("root:Packages:NIST:"+type+":data")
850        WAVE data_err=$("root:Packages:NIST:"+type+":linear_data_error")
851        WAVE new_reals=$("root:Packages:NIST:"+type+":realsread")
852       
853        Variable norm_mon,tot_mon,scale
854       
855        norm_mon = new_reals[0]         //should be 1e8
856        tot_mon = new_reals[1]          //should be equal to the number of runs (1 count per run)
857        scale= norm_mon/tot_mon
858       
859        data /= scale           //unscale the data
860        data_err /= scale
861       
862        // to keep "data" and linear_data in sync
863        data_copy = data
864       
865        return(0)
866End
867
Note: See TracBrowser for help on using the repository browser.