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

Last change on this file since 987 was 987, checked in by srkline, 7 years ago

lots of changes to how VCALC simulations can be written to Nexus files to effectively make test data with different callues before the virtual machine is ready. many changes for visualization to effectively handle zeros in log scaled images

File size: 29.2 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//
14// - copying workfiles to another folder
15//
16// - absolute scaling
17//
18// - (no) the WorkFile Math panel for simple image math (not done - maybe in the future?)
19// -
20// - (no) adding work.drk data without normalizing to monitor counts (the case not currently handled)
21//***************************
22
23//
24// Functions used for manipulation of the local Igor "WORK" folder
25// structure as raw data is displayed and processed.
26//
27//
28Strconstant ksDetectorListNoB = "FT;FB;FL;FR;MT;MB;ML;MR;"
29Strconstant ksDetectorListAll = "FT;FB;FL;FR;MT;MB;ML;MR;B;"
30
31
32//
33//Entry procedure from main panel
34//
35Proc CopyWorkFolder(oldType,newType)
36        String oldType,newType
37        Prompt oldType,"Source WORK data type",popup,"SAM;EMP;BGD;DIV;COR;CAL;RAW;ABS;STO;SUB;DRK;"
38        Prompt newType,"Destination WORK data type",popup,"SAM;EMP;BGD;DIV;COR;CAL;RAW;ABS;STO;SUB;DRK;"
39
40        // data folder "old" will be copied to "new" (either kills/copies or will overwrite)
41        CopyHDFToWorkFolder(oldtype,newtype)
42End
43
44//
45// copy what is needed for data processing (not the DAS_logs)
46// from the RawVSANS storage folder to the local WORK folder as needed
47//
48// TODO -- at what stage do I make copies of data in linear/log forms for data display?
49//                      -- when do I make the 2D error waves?
50//
51// TODO - decide what exactly I need to copy over. May be best to copy all, and delete
52//       what I know that I don't need
53//
54// TODO !!! DuplicateDataFolder will FAIL - in the base case of RAW data files, the
55//  data is actually in use - so it will fail every time. need an alternate solution. in SANS,
56// there are a limited number of waves to carry over, so Dupliate/O is used for rw, tw, data, etc.
57//
58// TODO : I also need a list of what is generated during processing that may be hanging around - that I need to
59//     be sure to get rid of - like the calibration waves, solidAngle, etc.
60//
61// hdfDF is the name only of the data in storage. May be full file name with extension (clean as needed)
62// type is the destination WORK folder for the copy
63//
64Function CopyHDFToWorkFolder(fromStr,toStr)
65        String fromStr,toStr
66       
67        String fromDF, toDF
68       
69        // make the DF paths - source and destination
70        fromDF = "root:Packages:NIST:VSANS:"+fromStr
71        toDF = "root:Packages:NIST:VSANS:"+toStr
72       
73//      // make a copy of the file name for my own use, since it's not in the file
74//      String/G $(toDF+":file_name") = root:
75       
76        // copy the folders
77        KillDataFolder/Z $toDF                  //DuplicateDataFolder will not overwrite, so Kill
78       
79        if(V_flag == 0)         // kill DF was OK
80                DuplicateDataFolder $("root:Packages:NIST:VSANS:"+fromStr),$("root:Packages:NIST:VSANS:"+toStr)
81               
82                // I can delete these if they came along with RAW
83                //   DAS_logs
84                //   top-level copies of data (duplicate links)
85                KillDataFolder/Z $(toDF+":entry:entry:DAS_logs")
86                KillDataFolder/Z $(toDF+":entry:entry:data")
87                KillDataFolder/Z $(toDF+":entry:entry:data_B")
88                KillDataFolder/Z $(toDF+":entry:entry:data_ML")
89                KillDataFolder/Z $(toDF+":entry:entry:data_MR")
90                KillDataFolder/Z $(toDF+":entry:entry:data_MT")
91                KillDataFolder/Z $(toDF+":entry:entry:data_MB")
92                KillDataFolder/Z $(toDF+":entry:entry:data_FL")
93                KillDataFolder/Z $(toDF+":entry:entry:data_FR")
94                KillDataFolder/Z $(toDF+":entry:entry:data_FT")
95                KillDataFolder/Z $(toDF+":entry:entry:data_FB")
96
97                return(0)
98        else
99                // need to do this the hard way, duplicate/O recursively
100                // see V_CopyToWorkFolder()
101               
102                // everything on the top level
103                V_DuplicateDataFolder($(fromDF+":entry:entry"),fromStr,toStr,0,"",0)    //no recursion here
104                // control
105                V_DuplicateDataFolder($(fromDF+":entry:entry:control"),fromStr,toStr,0,"",1)    //yes recursion here
106                // instrument
107                V_DuplicateDataFolder($(fromDF+":entry:entry:instrument"),fromStr,toStr,0,"",1) //yes recursion here
108                // reduction
109                V_DuplicateDataFolder($(fromDF+":entry:entry:reduction"),fromStr,toStr,0,"",1)  //yes recursion here
110                // sample
111                V_DuplicateDataFolder($(fromDF+":entry:entry:sample"),fromStr,toStr,0,"",1)     //yes recursion here
112
113        endif   
114       
115        return(0)
116end
117
118
119////////
120// see the help entry for IndexedDir for help on (possibly) how to do this faster
121// -- see the function Function ScanDirectories(pathName, printDirNames)
122//
123
124
125// from IgorExchange On July 17th, 2011 jtigor
126// started from "Recursively List Data Folder Contents"
127// Posted July 15th, 2011 by hrodstein
128//
129//
130//
131Proc V_CopyWorkFolderTest(dataFolderStr, fromStr, toStr, level, sNBName, recurse)
132        String dataFolderStr="root:Packages:NIST:VSANS:RAW"
133        String fromStr = "RAW"
134        String toStr="SAM"
135        Variable level=0
136        String sNBName="DataFolderTree"
137        Variable recurse = 1
138       
139        V_DuplicateDataFolder($dataFolderStr, fromStr, toStr, level, sNBName, recurse)
140
141
142end
143
144// ListDataFolder(dfr, level)
145// Recursively lists objects in data folder.
146// Pass data folder path for dfr and 0 for level.
147// pass level == 0 for the first call
148//  sNBName = "" prints nothing. any name will generate a notebook
149//
150// recurse == 0 will do only the specified folder, anything else will recurse all levels
151// toStr is the string name of the top-level folder only, not the full path
152//
153//
154Function V_DuplicateDataFolder(dfr, fromStr, toStr, level, sNBName,recurse)
155        DFREF dfr
156        String fromStr
157        String toStr
158        Variable level                  // Pass 0 to start
159        String sNBName
160        Variable recurse
161 
162        String name
163        String dfName
164        String sString
165       
166        String toDF = ""
167 
168        if (level == 0)         // this is the data folder, generate if needed in the destination
169                name = GetDataFolder(1, dfr)
170//              sPrintf sString, "%s (data folder)\r", name
171                toDF = ReplaceString(fromStr,name,toStr,1)              // case-sensitive replace
172                sprintf sString, "NewDataFolder/O %s\r",toDF
173                NewDataFolder/O $(RemoveEnding(toDF,":"))                       // remove trailing semicolon if it's there
174               
175                V_WriteBrowserInfo(sString, 1, sNBName)
176        endif
177 
178        dfName = GetDataFolder(1, dfr)
179        toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
180        Variable i
181 
182        String indentStr = "\t"
183        for(i=0; i<level; i+=1)
184                indentStr += "\t"
185        endfor
186 
187        Variable numWaves = CountObjectsDFR(dfr, 1)
188        for(i=0; i<numWaves; i+=1)
189                name = GetIndexedObjNameDFR(dfr, 1, i)
190                //
191                // wave type does not matter now. Duplicate does not care
192                //
193                sPrintf sString, "Duplicate/O  %s,  %s\r",dfName+name,toDF+name
194                Duplicate/O $(dfName+name),$(toDF+name)
195               
196                V_WriteBrowserInfo(sString, 2, sNBName)
197        endfor 
198 
199        Variable numNumericVariables = CountObjectsDFR(dfr, 2) 
200        for(i=0; i<numNumericVariables; i+=1)
201                name = GetIndexedObjNameDFR(dfr, 2, i)
202                sPrintf sString, "%s%s (numeric variable)\r", indentStr, name
203                V_WriteBrowserInfo(sString, 3, sNBName)
204        endfor 
205 
206        Variable numStringVariables = CountObjectsDFR(dfr, 3)   
207        for(i=0; i<numStringVariables; i+=1)
208                name = GetIndexedObjNameDFR(dfr, 3, i)
209                sPrintf sString, "%s%s (string variable)\r", indentStr, name
210                V_WriteBrowserInfo(sString, 4, sNBName)
211        endfor 
212
213        if(recurse)
214                Variable numDataFolders = CountObjectsDFR(dfr, 4)       
215                for(i=0; i<numDataFolders; i+=1)
216                        name = GetIndexedObjNameDFR(dfr, 4, i)
217//                      sPrintf sString, "%s%s (data folder)\r", indentStr, name
218                         dfName = GetDataFolder(1, dfr)
219                         
220                        toDF = ReplaceString(fromStr,dfName,toStr,1)            // case-sensitive replace
221                        sprintf sString, "NewDataFolder/O %s\r",toDF+name
222                        NewDataFolder/O $(toDF+name)
223                       
224                       
225                        V_WriteBrowserInfo(sString, 1, sNBName)
226                        DFREF childDFR = dfr:$(name)
227                        V_DuplicateDataFolder(childDFR, fromStr, toStr, level+1, sNBName, recurse)
228                endfor 
229        endif
230         
231//when finished walking tree, save as RTF with dialog   
232//      if(level == 0 && strlen(sNBName) != 0)
233//              SaveNotebook /I /S=4  $sNBName
234//      endif
235End
236 
237Function V_WriteBrowserInfo(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...
304Macro 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
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                        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        // -- currently, the "B" detector is skipped
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        // -- still not sure whether to duplicate these calculations as the RAW data is loaded. It would allow the RAW
420        //    data to be properly displayed, but without all of the (complete) set of 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 -other corrections modify the data
423        NVAR gDoNonLinearCor = root:Packages:NIST:VSANS:Globals:gDoNonLinearCor
424        // generate a distance matrix for each of the detectors
425        if (gDoNonLinearCor == 1)
426                Print "Doing Non-linear correction"// for "+ detStr
427                for(ii=0;ii<ItemsInList(ksDetectorListNoB);ii+=1)
428                        detStr = StringFromList(ii, ksDetectorListNoB, ";")
429                        Wave w = V_getDetectorDataW(fname,detStr)
430//                      Wave w_err = V_getDetectorDataErrW(fname,detStr)
431                        Wave w_calib = V_getDetTube_spatialCalib(fname,detStr)
432                        Variable tube_width = V_getDet_tubeWidth(fname,detStr)
433                        NonLinearCorrection(w,w_calib,tube_width,detStr,destPath)
434
435                        //(2.4) Convert the beam center values from pixels to mm
436                        // TODO -- there needs to be a permanent location for these values??
437                        //
438                        ConvertBeamCtr_to_mm(fname,detStr,destPath)
439                                                       
440                        // (2.5) Calculate the q-values
441                        // calculating q-values can't be done unless the non-linear corrections are calculated
442                        // so go ahead and put it in this loop.
443                        // TODO :
444                        // -- make sure that everything is present before the calculation
445                        // -- beam center must be properly defined in terms of real distance
446                        // -- distances/zero location/ etc. must be clearly documented for each detector
447                        //      ** this assumes that NonLinearCorrection() has been run to generate data_RealDistX and Y
448                        // ** this routine Makes the waves QTot, qx, qy, qz in each detector folder.
449                        //
450                        V_Detector_CalcQVals(fname,detStr,destPath)
451                       
452                endfor
453               
454                //"B" is separate
455                NonLinearCorrection_B(fname,"B",destPath)
456                ConvertBeamCtr_to_mmB(fname,"B",destPath)
457                V_Detector_CalcQVals(fname,"B",destPath)
458               
459        else
460                Print "Non-linear correction not done"
461        endif
462
463        // (3) solid angle correction
464        NVAR gDoSolidAngleCor = root:Packages:NIST:VSANS:Globals:gDoSolidAngleCor
465        if (gDoSolidAngleCor == 1)
466                Print "Doing Solid Angle correction"// for "+ detStr
467                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
468                        detStr = StringFromList(ii, ksDetectorListAll, ";")
469                        Wave w = V_getDetectorDataW(fname,detStr)
470                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
471                        Wave w_dt = V_getDetector_deadtime(fname,detStr)
472//                      SolidAngleCorrection(fill this in)
473                       
474                endfor
475        else
476                Print "Solid Angle correction not done"
477        endif   
478       
479        // (4) dead time correction
480        // TODO: -- remove the hard-wired test
481        // -- test for correct operation
482        // -- loop over all of the detectors
483        // -- B detector is a special case (do separately, then loop over NoB)
484        NVAR gDoDeadTimeCor = root:Packages:NIST:VSANS:Globals:gDoDeadTimeCor
485        ctTime = V_getCount_time(fname)
486        if (gDoDeadTimeCor == 1)
487                Print "Doing DeadTime correction"// for "+ detStr
488                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
489                        detStr = StringFromList(ii, ksDetectorListAll, ";")
490                        Wave w = V_getDetectorDataW(fname,detStr)
491                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
492                        Wave w_dt = V_getDetector_deadtime(fname,detStr)
493//                      DeadTimeCorrectionTubes(w,w_err,w_dt,ctTime)
494                                //deadtime corrections
495//      itim = integersread[2]
496//      cntrate = sum(data,-inf,inf)/itim               //use sum of detector counts rather than scaler value
497//      //TODO - do correct dead time correction for tubes
498//      deadtime = 1//DetectorDeadtime(textread[3],textread[9],dateAndTimeStr=textRead[1],dtime=realsRead[48])  //pick the correct deadtime
499//      dscale = 1/(1-deadTime*cntrate)
500//     
501       
502// dead time correction
503//      data *= dscale          //deadtime correction for everyone else, including NCNR
504//      data_err *= dscale
505
506                endfor
507        else
508                Print "Dead Time correction not done"
509        endif   
510       
511        // (5) angle-dependent tube shadowing
512        NVAR gDoTubeShadowCor = root:Packages:NIST:VSANS:Globals:gDoTubeShadowCor
513        if (gDoTubeShadowCor == 1)
514       
515        else
516                Print "Tube shadowing correction not done"
517        endif   
518               
519        // (6) angle dependent transmission correction
520        NVAR gDoTrans = root:Packages:NIST:VSANS:Globals:gDoTransmissionCor
521        if (gDoTrans == 1)
522                Print "Doing Large-angle transmission correction"// for "+ detStr
523                for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
524                        detStr = StringFromList(ii, ksDetectorListAll, ";")
525                        Wave w = V_getDetectorDataW(fname,detStr)
526                        Wave w_err = V_getDetectorDataErrW(fname,detStr)
527                        Wave w_dt = V_getDetector_deadtime(fname,detStr)
528//                      TransmissionCorrection(fill this in)
529                       
530                endfor
531        else
532                Print "Sample Transmission correction not done"
533        endif   
534       
535        // (7) normalize to default monitor counts
536        // TODO -- each detector is rescaled separately, but the rescaling factor is global (only one monitor!)
537        // TODO -- but there are TWO monitors - so how to switch?
538        // TODO -- what do I really need to save?
539        Print "Doing monitor normalization"// for "+ detStr
540
541        defmon=1e8                      //default monitor counts
542        for(ii=0;ii<ItemsInList(ksDetectorListAll);ii+=1)
543                detStr = StringFromList(ii, ksDetectorListAll, ";")
544                Wave w = V_getDetectorDataW(fname,detStr)
545                Wave w_err = V_getDetectorDataErrW(fname,detStr)
546                Variable monCt = V_getBeamMonNormData(fname)
547//                      MonitorNormalization(fill this in)
548        //scale the data to the default montor counts
549       
550        // TODO -- un-comment these three lines once monitor counts are reasonable - currently monCt = 9!!!
551//              scale = defmon/monCt
552//              w *= scale
553//              w_err *= scale          //assumes total monitor count is so large there is essentially no error
554
555// TODO
556// -- to write back to the local value, get the wave reference rather than the value, then I can
557//    re-assign the value directly, rather than this method (which is not terrible)     
558                // V_getBeamMonNormSaved_count()
559                // save the true monitor counts? save the scaling factor?
560                String path = "entry:instrument:beam_monitor_norm:saved_count"
561                Wave/Z savW = $("root:Packages:NIST:VSANS:"+fname+":entry:"+path)
562                savW[0] = scale
563        endfor
564       
565       
566        // (not done) angle dependent efficiency correction
567        NVAR doEfficiency = root:Packages:NIST:VSANS:Globals:gDoDetectorEffCor
568
569       
570// this function, in the past did the non-linear, solid angle, transmission, and efficiency corrections all at once
571//      DetCorr(data,data_err,realsread,doEfficiency,doTrans)           //the parameters are waves, and will be changed by the function
572
573
574       
575        //update totals to put in the work header (at the end of the function)
576//      total_mon += realsread[0]
577//
578//      total_det += dscale*realsread[2]
579//
580//      total_trn += realsread[39]
581//      total_rtime += integersread[2]
582//      total_numruns +=1
583//     
584
585        //all is done, except for the bookkeeping, updating the header information in the work folder
586
587//      integersread[3] = total_numruns                                         //numruns = 1
588//      realsread[1] = total_mon                        //save the true monitor count
589//      realsread[0] = defmon                                   //monitor ct = defmon
590//      realsread[2] = scale*total_det                  //scaled detector counts
591//     
592        //reset the current displaytype to "newtype"
593        String/G root:Packages:NIST:VSANS:Globals:gCurDispType=newType
594       
595        //return to root folder (redundant)
596        SetDataFolder root:
597       
598        Return(0)
599End
600
601
602//will "ADD" the current contents of the RAW folder to the newType work folder
603//and will ADD the RAW contents to the existing content of the newType folder
604// - used when adding multiple runs together
605//(the function Raw_to_work(type) makes a fresh workfile)
606//
607//the current display type is updated to newType (global)
608Function Add_raw_to_work(newType)
609        String newType
610       
611        // NEW OCT 2014
612        // this corrects for adding raw data files with different attenuation   
613        // does nothing if the attenuation of RAW and destination are the same
614        NVAR doAdjustRAW_Atten = root:Packages:NIST:gDoAdjustRAW_Atten
615        if(doAdjustRAW_Atten)
616                Adjust_RAW_Attenuation(newType)
617        endif
618       
619        String destPath=""
620       
621        // if the desired workfile doesn't exist, let the user know, and just make a new one
622        if(WaveExists($("root:Packages:NIST:"+newType + ":data")) == 0)
623                Print "There is no old work file to add to - a new one will be created"
624                //call Raw_to_work(), then return from this function
625                Raw_to_Work(newType)
626                Return(0)               //does not generate an error - a single file was converted to work.newtype
627        Endif
628       
629        NVAR pixelsX = root:myGlobals:gNPixelsX
630        NVAR pixelsY = root:myGlobals:gNPixelsY
631       
632        //now make references to data in newType folder
633        DestPath="root:Packages:NIST:"+newType 
634        WAVE data=$(destPath +":linear_data")                   // these wave references point to the EXISTING work data
635        WAVE data_copy=$(destPath +":data")                     // these wave references point to the EXISTING work data
636        WAVE dest_data_err=$(destPath +":linear_data_error")                    // these wave references point to the EXISTING work data
637        WAVE/T textread=$(destPath + ":textread")
638        WAVE integersread=$(destPath + ":integersread")
639        WAVE realsread=$(destPath + ":realsread")
640       
641        Variable deadTime,defmon,total_mon,total_det,total_trn,total_numruns,total_rtime
642        Variable ii,jj,itim,cntrate,dscale,scale,uscale,wrk_beamx,wrk_beamy,xshift,yshift
643
644
645        defmon=1e8                      //default monitor counts
646       
647        //Yes, add to previous run(s) in work, that does exist
648        //use the actual monitor count run.savmon rather than the normalized monitor count
649        //in run.moncnt and unscale the work data
650       
651        total_mon = realsread[1]        //saved monitor count
652        uscale = total_mon/defmon               //unscaling factor
653        total_det = uscale*realsread[2]         //unscaled detector count
654        total_trn = uscale*realsread[39]        //unscaled trans det count
655        total_numruns = integersread[3] //number of runs in workfile
656        total_rtime = integersread[2]           //total counting time in workfile
657        //retrieve workfile beamcenter
658        wrk_beamx = realsread[16]
659        wrk_beamy = realsread[17]
660        //unscale the workfile data in "newType"
661        //
662        //check for log-scaling and adjust if necessary
663        // should not be needed now - using display flag instead
664//      ConvertFolderToLinearScale(newType)
665        //
666        //then unscale the data array
667        data *= uscale
668        dest_data_err *= uscale
669       
670        //DetCorr() has not been applied to the data in RAW , do it now in a local reference to the raw data
671        WAVE raw_data = $"root:Packages:NIST:RAW:linear_data"
672        WAVE raw_data_err = $"root:Packages:NIST:RAW:linear_data_error"
673        WAVE raw_reals =  $"root:Packages:NIST:RAW:realsread"
674        WAVE/T raw_text = $"root:Packages:NIST:RAW:textread"
675        WAVE raw_ints = $"root:Packages:NIST:RAW:integersread"
676       
677        //check for log-scaling of the raw data - make sure it's linear
678        // should not be needed now - using display flag instead
679//      ConvertFolderToLinearScale("RAW")
680       
681        // switches to control what is done, don't do the transmission correction for the BGD measurement
682        NVAR doEfficiency = root:Packages:NIST:gDoDetectorEffCorr
683        NVAR gDoTrans = root:Packages:NIST:gDoTransmissionCorr
684        Variable doTrans = gDoTrans
685        if(cmpstr("BGD",newtype) == 0)
686                doTrans = 0             //skip the trans correction for the BGD file but don't change the value of the global
687        endif   
688       
689        DetCorr(raw_data,raw_data_err,raw_reals,doEfficiency,doTrans)   //applies correction to raw_data, and overwrites it
690       
691        //deadtime corrections to raw data
692        // TODO - do the tube correction for dead time now
693        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
694        itim = raw_ints[2]
695        cntrate = sum(raw_data,-inf,inf)/itim           //080802 use data sum, rather than scaler value
696        dscale = 1/(1-deadTime*cntrate)
697
698#if (exists("ILL_D22")==6)
699        Variable tubeSum
700        // for D22 detector might need to use cntrate/128 as it is the tube response
701        for(ii=0;ii<pixelsX;ii+=1)
702                //sum the counts in each tube
703                tubeSum = 0
704                for(jj=0;jj<pixelsY;jj+=1)
705                        tubeSum += data[jj][ii]
706                endfor
707                // countrate in tube ii
708                cntrate = tubeSum/itim
709                // deadtime scaling in tube ii
710                dscale = 1/(1-deadTime*cntrate)
711                // multiply data[ii][] by the dead time
712                raw_data[][ii] *= dscale
713                raw_data_err[][ii] *= dscale
714        endfor
715#else
716        // dead time correction on all other RAW data, including NCNR
717        raw_data *= dscale
718        raw_data_err *= dscale
719#endif
720
721        //update totals by adding RAW values to the local ones (write to work header at end of function)
722        total_mon += raw_reals[0]
723
724        total_det += dscale*raw_reals[2]
725
726        total_trn += raw_reals[39]
727        total_rtime += raw_ints[2]
728        total_numruns +=1
729       
730        //do the beamcenter shifting if there is a mismatch
731        //and then add the two data sets together, changing "data" since it is the workfile data
732        xshift = raw_reals[16] - wrk_beamx
733        yshift = raw_reals[17] - wrk_beamy
734       
735        If((xshift != 0) || (yshift != 0))
736                DoAlert 1,"Do you want to ignore the beam center mismatch?"
737                if(V_flag==1)
738                        xshift=0
739                        yshift=0
740                endif
741        endif
742       
743        If((xshift == 0) && (yshift == 0))              //no shift, just add them
744                data += raw_data                //deadtime correction has already been done to the raw data
745                dest_data_err = sqrt(dest_data_err^2 + raw_data_err^2)                  // error of the sum
746        Endif
747       
748        //scale the data to the default montor counts
749        scale = defmon/total_mon
750        data *= scale
751        dest_data_err *= scale
752       
753        // keep "data" and linear_data in sync in the destination folder
754        data_copy = data
755       
756        //all is done, except for the bookkeeping of updating the header info in the work folder
757        textread[1] = date() + " " + time()             //date + time stamp
758        integersread[3] = total_numruns                                         //numruns = more than one
759        realsread[1] = total_mon                        //save the true monitor count
760        realsread[0] = defmon                                   //monitor ct = defmon
761        integersread[2] = total_rtime                   // total counting time
762        realsread[2] = scale*total_det                  //scaled detector counts
763        realsread[39] = scale*total_trn                 //scaled transmission counts
764       
765        //Add the added raw filename to the list of files in the workfile
766        String newfile = ";" + raw_text[0]
767        SVAR oldList = $(destPath + ":fileList")
768        String/G $(destPath + ":fileList") = oldList + newfile
769       
770        //reset the current displaytype to "newtype"
771        String/G root:myGlobals:gDataDisplayType=newType
772       
773        //return to root folder (redundant)
774        SetDataFolder root:
775       
776        Return(0)
777End
778
779
780//used for adding DRK (beam shutter CLOSED) data to a workfile
781//force the monitor count to 1, since it's irrelevant
782// run data through normal "add" step, then unscale default monitor counts
783//to get the data back on a simple time basis
784//
785Function Raw_to_Work_NoNorm(type)
786        String type
787       
788        WAVE reals=$("root:Packages:NIST:RAW:realsread")
789        reals[1]=1              //true monitor counts, still in raw
790        Raw_to_work(type)
791        //data is now in "type" folder
792        WAVE data=$("root:Packages:NIST:"+type+":linear_data")
793        WAVE data_copy=$("root:Packages:NIST:"+type+":data")
794        WAVE data_err=$("root:Packages:NIST:"+type+":linear_data_error")
795        WAVE new_reals=$("root:Packages:NIST:"+type+":realsread")
796       
797        Variable norm_mon,tot_mon,scale
798       
799        norm_mon = new_reals[0]         //should be 1e8
800        tot_mon = new_reals[1]          //should be 1
801        scale= norm_mon/tot_mon
802       
803        data /= scale           //unscale the data
804        data_err /= scale
805       
806        // to keep "data" and linear_data in sync
807        data_copy = data
808       
809        return(0)
810End
811
812//used for adding DRK (beam shutter CLOSED) data to a workfile
813//force the monitor count to 1, since it's irrelevant
814// run data through normal "add" step, then unscale default monitor counts
815//to get the data back on a simple time basis
816//
817Function Add_Raw_to_Work_NoNorm(type)
818        String type
819       
820        WAVE reals=$("root:Packages:NIST:RAW:realsread")
821        reals[1]=1              //true monitor counts, still in raw
822        Add_Raw_to_work(type)
823        //data is now in "type" folder
824        WAVE data=$("root:Packages:NIST:"+type+":linear_data")
825        WAVE data_copy=$("root:Packages:NIST:"+type+":data")
826        WAVE data_err=$("root:Packages:NIST:"+type+":linear_data_error")
827        WAVE new_reals=$("root:Packages:NIST:"+type+":realsread")
828       
829        Variable norm_mon,tot_mon,scale
830       
831        norm_mon = new_reals[0]         //should be 1e8
832        tot_mon = new_reals[1]          //should be equal to the number of runs (1 count per run)
833        scale= norm_mon/tot_mon
834       
835        data /= scale           //unscale the data
836        data_err /= scale
837       
838        // to keep "data" and linear_data in sync
839        data_copy = data
840       
841        return(0)
842End
843
Note: See TracBrowser for help on using the repository browser.