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

Last change on this file since 1051 was 1051, checked in by srkline, 5 years ago

a lot of little changes:

changed the name of the Raw Data display procedure file (removed test)

lots of bug fixes, moving items from the macros menu to proper locations, getting the file status to display properly, some error checking, and cleaning up a few TODO items.

File size: 33.9 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 V_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        V_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 V_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, these should not be present in a proper NICE file)
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       
103                V_KillWavesFullTree($toDF,toStr,0,"",1)                 // this will traverse the whole tree, trying to kill what it can
104
105                if(DataFolderExists("root:Packages:NIST:VSANS:"+toStr) == 0)            //if the data folder (RAW, SAM, etc.) was just killed?
106                        NewDataFolder/O $("root:Packages:NIST:VSANS:"+toStr)
107                endif   
108                       
109                // need to do this the hard way, duplicate/O recursively
110                // see V_CopyToWorkFolder()
111               
112                // everything on the top level
113                V_DuplicateDataFolder($(fromDF+":entry"),fromStr,toStr,0,"",0)  //no recursion here
114                // control
115                V_DuplicateDataFolder($(fromDF+":entry:control"),fromStr,toStr,0,"",1)  //yes recursion here
116                // instrument
117                V_DuplicateDataFolder($(fromDF+":entry:instrument"),fromStr,toStr,0,"",1)       //yes recursion here
118                // reduction
119                V_DuplicateDataFolder($(fromDF+":entry:reduction"),fromStr,toStr,0,"",1)        //yes recursion here
120                // sample
121                V_DuplicateDataFolder($(fromDF+":entry:sample"),fromStr,toStr,0,"",1)   //yes recursion here
122                // user
123                V_DuplicateDataFolder($(fromDF+":entry:user"),fromStr,toStr,0,"",1)     //yes recursion here
124               
125        endif   
126       
127        return(0)
128end
129
130
131
132Function V_KillWavesInFolder(folderStr)
133        String folderStr
134       
135        if(DataFolderExists(folderStr) && strlen(folderStr) != 0)
136                SetDataFolder $folderStr
137                KillWaves/A/Z
138        endif
139       
140        SetDataFolder root:
141        return(0)
142end
143
144////////
145// see the help entry for IndexedDir for help on (possibly) how to do this faster
146// -- see  Function ScanDirectories(pathName, printDirNames)
147//
148
149
150// from IgorExchange On July 17th, 2011 jtigor
151// started from "Recursively List Data Folder Contents"
152// Posted July 15th, 2011 by hrodstein
153//
154//
155//
156Proc V_CopyWorkFolderProc(dataFolderStr, fromStr, toStr, level, sNBName, recurse)
157        String dataFolderStr="root:Packages:NIST:VSANS:RAW"
158        String fromStr = "RAW"
159        String toStr="SAM"
160        Variable level=0
161        String sNBName="DataFolderTree"
162        Variable recurse = 1
163       
164        V_DuplicateDataFolder($dataFolderStr, fromStr, toStr, level, sNBName, recurse)
165
166
167end
168
169// ListDataFolder(dfr, level)
170// Recursively lists objects in data folder.
171// Pass data folder path for dfr and 0 for level.
172// pass level == 0 for the first call
173//  sNBName = "" prints nothing. any name will generate a notebook
174//
175// recurse == 0 will do only the specified folder, anything else will recurse all levels
176// toStr is the string name of the top-level folder only, not the full path
177//
178//
179Function V_DuplicateDataFolder(dfr, fromStr, toStr, level, sNBName,recurse)
180        DFREF dfr
181        String fromStr
182        String toStr
183        Variable level                  // Pass 0 to start
184        String sNBName
185        Variable recurse
186 
187        String name
188        String dfName
189        String sString
190       
191        String toDF = ""
192 
193        if (level == 0)         // this is the data folder, generate if needed in the destination
194                name = GetDataFolder(1, dfr)
195//              sPrintf sString, "%s (data folder)\r", name
196                toDF = ReplaceString(fromStr,name,toStr,1)              // case-sensitive replace
197                sprintf sString, "NewDataFolder/O %s\r",toDF
198                NewDataFolder/O $(RemoveEnding(toDF,":"))                       // remove trailing semicolon if it's there
199               
200                V_WriteBrowserInfo_test(sString, 1, sNBName)
201        endif
202 
203        dfName = GetDataFolder(1, dfr)
204        toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
205        Variable i
206 
207        String indentStr = "\t"
208        for(i=0; i<level; i+=1)
209                indentStr += "\t"
210        endfor
211 
212        Variable numWaves = CountObjectsDFR(dfr, 1)
213        for(i=0; i<numWaves; i+=1)
214                name = GetIndexedObjNameDFR(dfr, 1, i)
215                //
216                // wave type does not matter now. Duplicate does not care
217                //
218                sPrintf sString, "Duplicate/O  %s,  %s\r",dfName+name,toDF+name
219                Duplicate/O $(dfName+name),$(toDF+name)
220               
221                V_WriteBrowserInfo_test(sString, 2, sNBName)
222        endfor 
223 
224        Variable numNumericVariables = CountObjectsDFR(dfr, 2) 
225        for(i=0; i<numNumericVariables; i+=1)
226                name = GetIndexedObjNameDFR(dfr, 2, i)
227                sPrintf sString, "%s%s (numeric variable)\r", indentStr, name
228                V_WriteBrowserInfo_test(sString, 3, sNBName)
229        endfor 
230 
231        Variable numStringVariables = CountObjectsDFR(dfr, 3)   
232        for(i=0; i<numStringVariables; i+=1)
233                name = GetIndexedObjNameDFR(dfr, 3, i)
234                sPrintf sString, "%s%s (string variable)\r", indentStr, name
235                V_WriteBrowserInfo_test(sString, 4, sNBName)
236        endfor 
237
238        if(recurse)
239                Variable numDataFolders = CountObjectsDFR(dfr, 4)       
240                for(i=0; i<numDataFolders; i+=1)
241                        name = GetIndexedObjNameDFR(dfr, 4, i)
242//                      sPrintf sString, "%s%s (data folder)\r", indentStr, name
243                         dfName = GetDataFolder(1, dfr)
244                         
245                        toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
246                        sprintf sString, "NewDataFolder/O %s\r",toDF+name
247                        NewDataFolder/O $(toDF+name)
248                       
249                       
250                        V_WriteBrowserInfo_test(sString, 1, sNBName)
251                        DFREF childDFR = dfr:$(name)
252                        V_DuplicateDataFolder(childDFR, fromStr, toStr, level+1, sNBName, recurse)
253                endfor 
254        endif
255         
256
257End
258
259
260// ListDataFolder(dfr, level)
261// Recursively lists objects in data folder.
262// Pass data folder path for dfr and 0 for level.
263// pass level == 0 for the first call
264//  sNBName = "" prints nothing. any name will generate a notebook
265//
266// recurse == 0 will do only the specified folder, anything else will recurse all levels
267// toStr is the string name of the top-level folder only, not the full path
268//
269//
270Function V_KillWavesFullTree(dfr, fromStr, level, sNBName,recurse)
271        DFREF dfr
272        String fromStr
273//      String toStr
274        Variable level                  // Pass 0 to start
275        String sNBName
276        Variable recurse
277 
278        String name
279        String dfName
280        String sString
281       
282        String toDF = ""
283 
284        if (level == 0)         // this is the data folder, generate if needed in the destination
285                name = GetDataFolder(1, dfr)
286                sPrintf sString, "%s (data folder)\r", name
287//              toDF = ReplaceString(fromStr,name,toStr,1)              // case-sensitive replace
288//              sprintf sString, "NewDataFolder/O %s\r",toDF
289//              NewDataFolder/O $(RemoveEnding(toDF,":"))                       // remove trailing semicolon if it's there
290               
291                V_WriteBrowserInfo_test(sString, 1, sNBName)
292        endif
293 
294        dfName = GetDataFolder(1, dfr)
295//      toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
296        Variable i
297 
298        String indentStr = "\t"
299        for(i=0; i<level; i+=1)
300                indentStr += "\t"
301        endfor
302 
303        Variable numWaves = CountObjectsDFR(dfr, 1)
304        for(i=0; i<numWaves; i+=1)
305                name = GetIndexedObjNameDFR(dfr, 1, i)
306                //
307                // wave type does not matter now. Kill does not care
308                //
309                sPrintf sString, "Killing  %s\r",dfName+name
310                KillWaves/Z $(dfName+name)
311               
312                V_WriteBrowserInfo_test(sString, 2, sNBName)
313        endfor 
314 
315 // now kill the data folder if possible
316        KillDataFolder/Z $dfName
317       
318       
319        Variable numNumericVariables = CountObjectsDFR(dfr, 2) 
320        for(i=0; i<numNumericVariables; i+=1)
321                name = GetIndexedObjNameDFR(dfr, 2, i)
322                sPrintf sString, "%s%s (numeric variable)\r", indentStr, name
323                V_WriteBrowserInfo_test(sString, 3, sNBName)
324        endfor 
325 
326        Variable numStringVariables = CountObjectsDFR(dfr, 3)   
327        for(i=0; i<numStringVariables; i+=1)
328                name = GetIndexedObjNameDFR(dfr, 3, i)
329                sPrintf sString, "%s%s (string variable)\r", indentStr, name
330                V_WriteBrowserInfo_test(sString, 4, sNBName)
331        endfor 
332
333        if(recurse)
334                Variable numDataFolders = CountObjectsDFR(dfr, 4)       
335                for(i=0; i<numDataFolders; i+=1)
336                        name = GetIndexedObjNameDFR(dfr, 4, i)
337                        sPrintf sString, "%s%s (data folder)\r", indentStr, name
338                         dfName = GetDataFolder(1, dfr)
339                         
340//                      toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
341//                      sprintf sString, "NewDataFolder/O %s\r",toDF+name
342//                      NewDataFolder/O $(toDF+name)
343                       
344                       
345                        V_WriteBrowserInfo_test(sString, 1, sNBName)
346                        DFREF childDFR = dfr:$(name)
347                        V_KillWavesFullTree(childDFR, fromStr, level+1, sNBName, recurse)
348                endfor 
349        endif
350         
351
352End
353 
354
355
356Function V_WriteBrowserInfo_test(sString, vType, sNBName)
357        String sString
358        Variable vType
359        String sNBName
360 
361        if(strlen(sNBName) == 0)
362//              print sString
363                return 0
364        endif
365        DoWindow $sNBName
366        if(V_flag != 1)
367                NewNoteBook/F=0 /N=$sNBName /V=1 as sNBName
368        else
369                DoWindow/F $sNBName
370        endif
371        Notebook $sNBName selection={endOfFile, endOfFile}
372        if(vType == 1)          // a data folder
373//              Notebook $sNBName fstyle=1
374                Notebook $sNBName text=sString
375//              Notebook $sNBName fstyle=-1
376        else
377                Notebook $sNBName text=sString 
378        endif
379 
380End
381
382///////////////////////////////
383
384
385//
386// given the folder, duplicate the data -> linear_data and generate the error
387// TODO
388// -- do I want to use different names here? If it turns out that I don't need to drag a copy of
389//    the data around as "linear_data", then I can eliminate that, and rename the error wave
390// -- be sure the data is either properly written as 2D in the file, or converted to 2D before
391//    duplicating here
392// -- ? do I recast to DP here. Probably necessary since I'm doing a DP calculation, but Redimension
393//    is done in the Raw_to_Work step too. very confusing.
394Function V_MakeDataError(folderStr)
395        String folderStr
396       
397        SetDataFolder $folderStr
398        Wave data=data
399        Duplicate/O data linear_data                    // at this point, the data is still the raw data, and is linear_data
400       
401        // proper error for counting statistics, good for low count values too
402        // rather than just sqrt(n)
403        // see N. Gehrels, Astrophys. J., 303 (1986) 336-346, equation (7)
404        // for S = 1 in eq (7), this corresponds to one sigma error bars
405        Duplicate/O linear_data linear_data_error
406        linear_data_error = 1 + sqrt(linear_data + 0.75)                               
407        //
408       
409        SetDataFolder root:
410        return(0)
411End
412
413
414/////////////////////
415
416
417
418//testing procedure
419// TODO -- can't duplicate this with another proceudre, but if I change the name of the variable
420//   "newType" to "type", then when Raw_to_work() gets to CopyHDFToWorkFolder(), the KillDataFolder/Z
421//   line fails (but reports no error), then DuplicateDataFolder fails, and reports an error. Trying
422//   to simplify this condition, I can't reproduce the bug for WM...
423Proc V_Convert_to_Workfile(newtype, doadd)
424        String newtype,doadd
425        Prompt newtype,"WORK data type",popup,"SAM;EMP;BGD;ADJ;"
426        Prompt doadd,"Add to current WORK contents?",popup,"No;Yes;"
427       
428        //macro will take whatever is in RAW folder and "ADD" it to the folder specified
429        //in the popup menu
430       
431        //"add" = yes/no, don't add to previous runs
432        //switch here - two separate functions to avoid (my) confusion
433        Variable err// = Raw_to_work(newtype)
434        if(cmpstr(doadd,"No")==0)
435                //don't add to prev work contents, copy RAW contents to work and convert
436                err = V_Raw_to_work(newtype)
437        else
438                //yes, add RAW to the current work folder contents
439                Abort "Adding RAW data files is currently unsupported"
440                err = V_Add_raw_to_work(newtype)
441        endif
442       
443        String newTitle = "WORK_"+newtype
444        DoWindow/F VSANS_Data
445        DoWindow/T VSANS_Data, newTitle
446        KillStrings/Z newTitle
447       
448        //need to update the display with "data" from the correct dataFolder
449        V_UpdateDisplayInformation(newtype)
450       
451End
452
453
454//
455// THIS IS THE MAJOR ROUTINE TO APPLY DATA CORRECTIONS
456//
457//will copy the current contents of the RAW folder to the newType work folder
458//and do the geometric corrections and normalization to monitor counts
459//(the function Add_Raw_to_work(type) adds multiple runs together - and is LOW priority)
460//
461//the current display type is updated to newType (global)
462//
463Function V_Raw_to_work(newType)
464        String newType
465       
466        Variable deadTime,defmon,total_mon,total_det,total_trn,total_numruns,total_rtime
467        Variable ii,jj,itim,cntrate,dscale,scale,uscale
468        String destPath
469       
470        String fname = newType
471        String detStr
472        Variable ctTime
473
474        //initialize values before normalization
475        total_mon=0
476        total_det=0
477        total_trn=0
478        total_numruns=0
479        total_rtime=0
480       
481        //Not adding multiple runs, so wipe out the old contents of the work folder and
482        // replace with the contents of raw
483
484        destPath = "root:Packages:NIST:VSANS:" + newType
485       
486        //copy from current dir (RAW) to work, defined by newType
487        V_CopyHDFToWorkFolder("RAW",newType)
488       
489        // now work with the waves from the destination folder.
490       
491        // apply corrections ---
492        // switches to control what is done, don't do the transmission correction for the BGD measurement
493        // start with the DIV correction, before conversion to mm
494        // then do all of the other corrections, order doesn't matter.
495        // rescaling to default monitor counts however, must be LAST.
496
497// each correction must loop over each detector. tedious.
498
499        // (0) Redimension the data waves in the destination folder
500        //     so that they are DP, not integer
501        // TODO
502        // -- currently only redimensioning the data and linear_data_error - What else???
503        // -- ?? some of this is done at load time for RAW data. shouldn't be an issue to re-do the redimension
504        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
505                detStr = StringFromList(ii, ksDetectorListAll, ";")
506                Wave w = V_getDetectorDataW(fname,detStr)
507                Wave w_err = V_getDetectorDataErrW(fname,detStr)
508                Redimension/D w,w_err
509        endfor
510       
511       
512        // (1) DIV correction
513        // do this in terms of pixels.
514        // TODO : This must also exist at the time the first work folder is generated.
515        //   So it must be in the user folder at the start of the experiment, and defined.
516        NVAR gDoDIVCor = root:Packages:NIST:VSANS:Globals:gDoDIVCor
517        if (gDoDIVCor == 1)
518                // need extra check here for file existence
519                // if not in DIV folder, load.
520                // if unable to load, skip correction and report error (Alert?) (Ask to Load?)
521                Print "Doing DIV correction"// for "+ detStr
522                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
523                        detStr = StringFromList(ii, ksDetectorListAll, ";")
524                        Wave w = V_getDetectorDataW(fname,detStr)
525                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
526                       
527                        V_DIVCorrection(w,w_err,detStr,newType)
528                endfor
529        else
530                Print "DIV correction NOT DONE"         // not an error since correction was unchecked
531        endif
532       
533        // (2) non-linear correction   
534        // TODO:
535        // x-  the "B" detector is calculated in its own routines
536        // -- document what is generated here:
537        //    **in each detector folder: data_realDistX and data_realDistY (2D waves of the mm? position of the pixel)
538        // x- these spatial calculations ARE DONE as the RAW data is loaded. It allows the RAW
539        //    data to be properly displayed, but without all of the (complete) set of detector corrections
540        // * the corrected distances are calculated into arrays, but nothing is done with them yet
541        // * there is enough information now to calculate the q-arrays, so it is done now
542        // - other corrections may modify the data, this calculation does NOT modify the data
543        NVAR gDoNonLinearCor = root:Packages:NIST:VSANS:Globals:gDoNonLinearCor
544        // generate a distance matrix for each of the detectors
545        if (gDoNonLinearCor == 1)
546                Print "Doing Non-linear correction"// for "+ detStr
547                for(ii=0;ii<ItemsInList(ksDetectorListNoB);ii+=1)
548                        detStr = StringFromList(ii, ksDetectorListNoB, ";")
549                        Wave w = V_getDetectorDataW(fname,detStr)
550//                      Wave w_err = V_getDetectorDataErrW(fname,detStr)
551                        Wave w_calib = V_getDetTube_spatialCalib(fname,detStr)
552                        Variable tube_width = V_getDet_tubeWidth(fname,detStr)
553                        V_NonLinearCorrection(w,w_calib,tube_width,detStr,destPath)
554
555                        //(2.4) Convert the beam center values from pixels to mm
556                        // TODO -- there needs to be a permanent location for these values??
557                        //
558                        V_ConvertBeamCtr_to_mm(fname,detStr,destPath)
559                                                       
560                        // (2.5) Calculate the q-values
561                        // calculating q-values can't be done unless the non-linear corrections are calculated
562                        // so go ahead and put it in this loop.
563                        // TODO :
564                        // -- make sure that everything is present before the calculation
565                        // -- beam center must be properly defined in terms of real distance
566                        // -- distances/zero location/ etc. must be clearly documented for each detector
567                        //      ** this assumes that NonLinearCorrection() has been run to generate data_RealDistX and Y
568                        // ** this routine Makes the waves QTot, qx, qy, qz in each detector folder.
569                        //
570                        V_Detector_CalcQVals(fname,detStr,destPath)
571                       
572                endfor
573               
574                //"B" is separate
575                V_NonLinearCorrection_B(fname,"B",destPath)
576                V_ConvertBeamCtr_to_mmB(fname,"B",destPath)
577                V_Detector_CalcQVals(fname,"B",destPath)
578               
579        else
580                Print "Non-linear correction NOT DONE"
581        endif
582
583
584        // (3) dead time correction
585        // TODO:
586        // x- remove the hard-wired test - done
587        // -- test for correct operation
588        // x- loop over all of the detectors
589        // x- B detector is a special case (do separately, then loop over NoB)
590        // -- this DOES alter the data
591        // -- verify the error propagation (not done yet)
592        //
593        Variable countRate
594        NVAR gDoDeadTimeCor = root:Packages:NIST:VSANS:Globals:gDoDeadTimeCor
595        if (gDoDeadTimeCor == 1)
596                Print "Doing DeadTime correction"// for "+ detStr
597                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
598                        detStr = StringFromList(ii, ksDetectorListAll, ";")
599                        Wave w = V_getDetectorDataW(fname,detStr)
600                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
601                        ctTime = V_getCount_time(fname)
602
603                        if(cmpstr(detStr,"B") == 0)
604                                Variable b_dt = V_getDetector_deadtime_B(fname,detStr)
605                                // do the correction for the back panel
606                                countRate = sum(w,-inf,inf)/ctTime              //use sum of detector counts
607
608                                w = w/(1-countRate*b_dt)
609                                w_err = w_err/(1-countRate*b_dt)
610                                                               
611                        else
612                                // do the corrections for 8 tube panels
613                                Wave w_dt = V_getDetector_deadtime(fname,detStr)
614                                V_DeadTimeCorrectionTubes(w,w_err,w_dt,ctTime)
615                               
616                        endif
617                endfor
618               
619        else
620                Print "Dead Time correction NOT DONE"
621        endif   
622       
623
624        // (4) solid angle correction
625        // TODO -- this currently calculates the correction factor AND applys it to the data
626        //  -- as a result, the data values are very large since they are divided by a very small
627        //     solid angle per pixel. But all of the count values are now on the basis of
628        //    counts/(solid angle) --- meaning that they can all be binned together for I(q)
629        //    -and- TODO - this will need to be taken into account for absolute scaling (this part is already done)
630        //
631        // the solid angle correction is calculated for ALL detector panels.
632        NVAR gDoSolidAngleCor = root:Packages:NIST:VSANS:Globals:gDoSolidAngleCor
633        if (gDoSolidAngleCor == 1)
634                Print "Doing Solid Angle correction"// for "+ detStr
635                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
636                        detStr = StringFromList(ii, ksDetectorListAll, ";")
637                        Wave w = V_getDetectorDataW(fname,detStr)
638                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
639                        // any other dimensions to pass in?
640                        V_SolidAngleCorrection(w,w_err,fname,detStr,destPath)
641                       
642                endfor
643        else
644                Print "Solid Angle correction NOT DONE"
645        endif   
646       
647               
648        // (5) angle-dependent tube shadowing
649        // TODO:
650        // -- not sure about this correction yet...
651        //
652        NVAR gDoTubeShadowCor = root:Packages:NIST:VSANS:Globals:gDoTubeShadowCor
653        if (gDoTubeShadowCor == 1)
654                Print "(stub)Tube Shadow correction"
655        else
656                Print "Tube shadowing correction NOT DONE"
657        endif   
658               
659        // (6) angle dependent transmission correction
660        // TODO:
661        // x- this still needs to be filled in
662        // -- still some debate of when/where in the corrections that this is best applied
663        //    - do it here, and it's done whether the output is 1D or 2D
664        //    - do it later (where SAMPLE information is used) since this section is ONLY instrument-specific
665        // -- verify that the calculation is correct
666        // -- verify that the error propagation (in 2D) is correct
667        //
668        NVAR gDoTrans = root:Packages:NIST:VSANS:Globals:gDoTransmissionCor
669        if (gDoTrans == 1)
670                Print "Doing Large-angle transmission correction"// for "+ detStr
671                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
672                        detStr = StringFromList(ii, ksDetectorListAll, ";")
673                        Wave w = V_getDetectorDataW(fname,detStr)
674                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
675                       
676                        V_LargeAngleTransmissionCorr(w,w_err,fname,detStr,destPath)
677                       
678                endfor
679        else
680                Print "Sample Transmission correction NOT DONE"
681        endif   
682       
683        // (7) normalize to default monitor counts
684        // TODO(DONE) x- each detector is rescaled separately, but the rescaling factor is global (only one monitor!)
685        // TODO -- but there are TWO monitors - so how to switch?
686        //  --- AND, there is also /entry/control/monitor_counts !!! Which one is the correct value? Which will NICE write
687        // TODO -- what do I really need to save?
688       
689        NVAR gDoMonitorNormalization = root:Packages:NIST:VSANS:Globals:gDoMonitorNormalization
690        if (gDoMonitorNormalization == 1)
691               
692                Variable monCount,savedMonCount
693                defmon=1e8                      //default monitor counts
694                monCount = V_getMonitorCount(fname)                     // TODO -- this is read in since VCALC fakes this on output
695//              monCount = V_getBeamMonNormData(fname)          // TODO -- I think this is the *real* one to read
696                savedMonCount   = monCount
697                scale = defMon/monCount         // scale factor to MULTIPLY data by to rescale to defmon
698
699                // write to newType=fname will put new values in the destination WORK folder
700                V_writeBeamMonNormSaved_count(fname,savedMonCount)                      // save the true count
701                V_writeBeamMonNormData(fname,defMon)            // mon ct is now 10^8
702                                       
703//                      // TODO
704//                      // the low efficiency monitor, expect to use this for white beam mode
705//                              -- need code switch here to determine which monitor to use.
706//
707//                      V_getBeamMonLowData(fname)
708//                      V_getBeamMonLowSaved_count(fname)       
709                       
710                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
711                        detStr = StringFromList(ii, ksDetectorListAll, ";")
712                        Wave w = V_getDetectorDataW(fname,detStr)
713                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
714
715                        // do the calculation right here. It's a simple scaling and not worth sending to another function.     
716                        //scale the data and error to the default montor counts
717               
718//
719                        w *= scale
720                        w_err *= scale          //assumes total monitor count is so large there is essentially no error
721
722                        // TODO
723                        // -- do I want to update and save the integrated detector count?
724                        Variable integratedCount = V_getDet_IntegratedCount(fname,detStr)
725                        V_writeDet_IntegratedCount(fname,detStr,integratedCount*scale)
726
727                endfor
728        else
729                Print "Monitor Normalization correction NOT DONE"
730        endif
731       
732       
733        // (not done) angle dependent efficiency correction
734        NVAR doEfficiency = root:Packages:NIST:VSANS:Globals:gDoDetectorEffCor
735
736//
737///// these lines are if files are added together, not done (yet) for VSANS
738//
739//update totals to put in the work header (at the end of the function)
740//      total_mon += realsread[0]
741//
742//      total_det += dscale*realsread[2]
743//
744//      total_trn += realsread[39]
745//      total_rtime += integersread[2]
746//      total_numruns +=1
747//     
748
749//     
750        //reset the current displaytype to "newtype"
751        String/G root:Packages:NIST:VSANS:Globals:gCurDispType=newType
752       
753        //return to root folder (redundant)
754        SetDataFolder root:
755       
756        Return(0)
757End
758
759
760//will "ADD" the current contents of the RAW folder to the newType work folder
761//and will ADD the RAW contents to the existing content of the newType folder
762// - used when adding multiple runs together
763//(the function Raw_to_work(type) makes a fresh workfile)
764//
765//the current display type is updated to newType (global)
766Function V_Add_raw_to_work(newType)
767        String newType
768       
769        // NEW OCT 2014
770        // this corrects for adding raw data files with different attenuation   
771        // does nothing if the attenuation of RAW and destination are the same
772        NVAR doAdjustRAW_Atten = root:Packages:NIST:gDoAdjustRAW_Atten
773        if(doAdjustRAW_Atten)
774                V_Adjust_RAW_Attenuation(newType)
775        endif
776       
777        String destPath=""
778       
779        // if the desired workfile doesn't exist, let the user know, and just make a new one
780        if(WaveExists($("root:Packages:NIST:"+newType + ":data")) == 0)
781                Print "There is no old work file to add to - a new one will be created"
782                //call Raw_to_work(), then return from this function
783                V_Raw_to_Work(newType)
784                Return(0)               //does not generate an error - a single file was converted to work.newtype
785        Endif
786       
787        NVAR pixelsX = root:myGlobals:gNPixelsX
788        NVAR pixelsY = root:myGlobals:gNPixelsY
789       
790        //now make references to data in newType folder
791        DestPath="root:Packages:NIST:"+newType 
792        WAVE data=$(destPath +":linear_data")                   // these wave references point to the EXISTING work data
793        WAVE data_copy=$(destPath +":data")                     // these wave references point to the EXISTING work data
794        WAVE dest_data_err=$(destPath +":linear_data_error")                    // these wave references point to the EXISTING work data
795        WAVE/T textread=$(destPath + ":textread")
796        WAVE integersread=$(destPath + ":integersread")
797        WAVE realsread=$(destPath + ":realsread")
798       
799        Variable deadTime,defmon,total_mon,total_det,total_trn,total_numruns,total_rtime
800        Variable ii,jj,itim,cntrate,dscale,scale,uscale,wrk_beamx,wrk_beamy,xshift,yshift
801
802
803        defmon=1e8                      //default monitor counts
804       
805        //Yes, add to previous run(s) in work, that does exist
806        //use the actual monitor count run.savmon rather than the normalized monitor count
807        //in run.moncnt and unscale the work data
808       
809        total_mon = realsread[1]        //saved monitor count
810        uscale = total_mon/defmon               //unscaling factor
811        total_det = uscale*realsread[2]         //unscaled detector count
812        total_trn = uscale*realsread[39]        //unscaled trans det count
813        total_numruns = integersread[3] //number of runs in workfile
814        total_rtime = integersread[2]           //total counting time in workfile
815        //retrieve workfile beamcenter
816        wrk_beamx = realsread[16]
817        wrk_beamy = realsread[17]
818        //unscale the workfile data in "newType"
819        //
820        //check for log-scaling and adjust if necessary
821        // should not be needed now - using display flag instead
822//      ConvertFolderToLinearScale(newType)
823        //
824        //then unscale the data array
825        data *= uscale
826        dest_data_err *= uscale
827       
828        //DetCorr() has not been applied to the data in RAW , do it now in a local reference to the raw data
829        WAVE raw_data = $"root:Packages:NIST:RAW:linear_data"
830        WAVE raw_data_err = $"root:Packages:NIST:RAW:linear_data_error"
831        WAVE raw_reals =  $"root:Packages:NIST:RAW:realsread"
832        WAVE/T raw_text = $"root:Packages:NIST:RAW:textread"
833        WAVE raw_ints = $"root:Packages:NIST:RAW:integersread"
834       
835        //check for log-scaling of the raw data - make sure it's linear
836        // should not be needed now - using display flag instead
837//      ConvertFolderToLinearScale("RAW")
838       
839        // switches to control what is done, don't do the transmission correction for the BGD measurement
840        NVAR doEfficiency = root:Packages:NIST:gDoDetectorEffCor
841        NVAR gDoTrans = root:Packages:NIST:gDoTransmissionCor
842        Variable doTrans = gDoTrans
843        if(cmpstr("BGD",newtype) == 0)
844                doTrans = 0             //skip the trans correction for the BGD file but don't change the value of the global
845        endif   
846       
847        V_DetCorr(raw_data,raw_data_err,raw_reals,doEfficiency,doTrans) //applies correction to raw_data, and overwrites it
848       
849        //deadtime corrections to raw data
850        // TODO - do the tube correction for dead time now
851        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
852        itim = raw_ints[2]
853        cntrate = sum(raw_data,-inf,inf)/itim           //080802 use data sum, rather than scaler value
854        dscale = 1/(1-deadTime*cntrate)
855
856#if (exists("ILL_D22")==6)
857        Variable tubeSum
858        // for D22 detector might need to use cntrate/128 as it is the tube response
859        for(ii=0;ii<pixelsX;ii+=1)
860                //sum the counts in each tube
861                tubeSum = 0
862                for(jj=0;jj<pixelsY;jj+=1)
863                        tubeSum += data[jj][ii]
864                endfor
865                // countrate in tube ii
866                cntrate = tubeSum/itim
867                // deadtime scaling in tube ii
868                dscale = 1/(1-deadTime*cntrate)
869                // multiply data[ii][] by the dead time
870                raw_data[][ii] *= dscale
871                raw_data_err[][ii] *= dscale
872        endfor
873#else
874        // dead time correction on all other RAW data, including NCNR
875        raw_data *= dscale
876        raw_data_err *= dscale
877#endif
878
879        //update totals by adding RAW values to the local ones (write to work header at end of function)
880        total_mon += raw_reals[0]
881
882        total_det += dscale*raw_reals[2]
883
884        total_trn += raw_reals[39]
885        total_rtime += raw_ints[2]
886        total_numruns +=1
887       
888        //do the beamcenter shifting if there is a mismatch
889        //and then add the two data sets together, changing "data" since it is the workfile data
890        xshift = raw_reals[16] - wrk_beamx
891        yshift = raw_reals[17] - wrk_beamy
892       
893        If((xshift != 0) || (yshift != 0))
894                DoAlert 1,"Do you want to ignore the beam center mismatch?"
895                if(V_flag==1)
896                        xshift=0
897                        yshift=0
898                endif
899        endif
900       
901        If((xshift == 0) && (yshift == 0))              //no shift, just add them
902                data += raw_data                //deadtime correction has already been done to the raw data
903                dest_data_err = sqrt(dest_data_err^2 + raw_data_err^2)                  // error of the sum
904        Endif
905       
906        //scale the data to the default montor counts
907        scale = defmon/total_mon
908        data *= scale
909        dest_data_err *= scale
910       
911        // keep "data" and linear_data in sync in the destination folder
912        data_copy = data
913       
914        //all is done, except for the bookkeeping of updating the header info in the work folder
915        textread[1] = date() + " " + time()             //date + time stamp
916        integersread[3] = total_numruns                                         //numruns = more than one
917        realsread[1] = total_mon                        //save the true monitor count
918        realsread[0] = defmon                                   //monitor ct = defmon
919        integersread[2] = total_rtime                   // total counting time
920        realsread[2] = scale*total_det                  //scaled detector counts
921        realsread[39] = scale*total_trn                 //scaled transmission counts
922       
923        //Add the added raw filename to the list of files in the workfile
924        String newfile = ";" + raw_text[0]
925        SVAR oldList = $(destPath + ":fileList")
926        String/G $(destPath + ":fileList") = oldList + newfile
927       
928        //reset the current display type to "newtype"
929        SVAR gCurDispType = root:Packages:NIST:VSANS:Globals:gCurDispType
930        gCurDispType = newType
931       
932        //return to root folder (redundant)
933        SetDataFolder root:
934       
935        Return(0)
936End
937
938
939//used for adding DRK (beam shutter CLOSED) data to a workfile
940//force the monitor count to 1, since it's irrelevant
941// run data through normal "add" step, then unscale default monitor counts
942//to get the data back on a simple time basis
943//
944Function V_Raw_to_Work_NoNorm(type)
945        String type
946       
947        WAVE reals=$("root:Packages:NIST:RAW:realsread")
948        reals[1]=1              //true monitor counts, still in raw
949        V_Raw_to_work(type)
950        //data is now in "type" folder
951        WAVE data=$("root:Packages:NIST:"+type+":linear_data")
952        WAVE data_copy=$("root:Packages:NIST:"+type+":data")
953        WAVE data_err=$("root:Packages:NIST:"+type+":linear_data_error")
954        WAVE new_reals=$("root:Packages:NIST:"+type+":realsread")
955       
956        Variable norm_mon,tot_mon,scale
957       
958        norm_mon = new_reals[0]         //should be 1e8
959        tot_mon = new_reals[1]          //should be 1
960        scale= norm_mon/tot_mon
961       
962        data /= scale           //unscale the data
963        data_err /= scale
964       
965        // to keep "data" and linear_data in sync
966        data_copy = data
967       
968        return(0)
969End
970
971//used for adding DRK (beam shutter CLOSED) data to a workfile
972//force the monitor count to 1, since it's irrelevant
973// run data through normal "add" step, then unscale default monitor counts
974//to get the data back on a simple time basis
975//
976Function V_Add_Raw_to_Work_NoNorm(type)
977        String type
978       
979        WAVE reals=$("root:Packages:NIST:RAW:realsread")
980        reals[1]=1              //true monitor counts, still in raw
981        V_Add_Raw_to_work(type)
982        //data is now in "type" folder
983        WAVE data=$("root:Packages:NIST:"+type+":linear_data")
984        WAVE data_copy=$("root:Packages:NIST:"+type+":data")
985        WAVE data_err=$("root:Packages:NIST:"+type+":linear_data_error")
986        WAVE new_reals=$("root:Packages:NIST:"+type+":realsread")
987       
988        Variable norm_mon,tot_mon,scale
989       
990        norm_mon = new_reals[0]         //should be 1e8
991        tot_mon = new_reals[1]          //should be equal to the number of runs (1 count per run)
992        scale= norm_mon/tot_mon
993       
994        data /= scale           //unscale the data
995        data_err /= scale
996       
997        // to keep "data" and linear_data in sync
998        data_copy = data
999       
1000        return(0)
1001End
1002
Note: See TracBrowser for help on using the repository browser.