source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/ANSTO_DataReadWrite.ipf

Last change on this file was 940, checked in by srkline, 8 years ago

Updated the reduction Read/Write? and corresponding factility stubs to accomodate detector dead time written to the VAX header. It is currently not written to the header, but may be with NICE (hopefully).

With the move of NG3 SANS to CGB(upper), the NG3 designation in the account name [NGxSANSxx] has been replaced with CGB. Folders on charlotte and radio button on SASCALC are labeled "NGB30" to be more obvious which instrument it is.

FFT routines have minor typos cleaned up.

File size: 55.2 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma version=5.0
3#pragma IgorVersion = 6.1 //required to read and write files created with HDF 1.8 library
4
5//**************************
6//
7// Vers. 1.2 092101
8// Vers. 5.0 29MAR07 - branched from main reduction to split out facility
9//                     specific calls
10//
11// functions for reading raw data files from the VAX
12// - RAW data files are read into the RAW folder - integer data from the detector
13//   is decompressed and given the proper orientation
14// - header information is placed into real,integer, or text waves in the order they appear
15//   in the file header
16//
17// Work data (DIV File) is read into the DIV folder
18//
19//*****************************
20
21//simple, main entry procedure that will load a RAW sans data file (not a work file)
22//into the RAW dataFolder. It is up to the calling procedure to display the file
23//
24// called by MainPanel.ipf and ProtocolAsPanel.ipf
25//
26Function LoadRawSANSData(msgStr)
27        String msgStr
28
29        String filename=""
30
31        //each routine is responsible for checking the current (displayed) data folder
32        //selecting it, and returning to root when done
33        PathInfo/S catPathName          //should set the next dialog to the proper path...
34        //get the filename, then read it in
35        filename = PromptForPath(msgStr)                //in SANS_Utils.ipf
36        //check for cancel from dialog
37        if(strlen(filename)==0)
38                //user cancelled, abort
39                SetDataFolder root:
40                DoAlert 0, "No file selected, action aborted"
41                return(1)
42        Endif
43
44        ReadHeaderAndData(filename)     //this is the full Path+file
45       
46        Return(0)
47End
48
49
50//function that does the guts of reading the binary data file
51//fname is the full path:name;vers required to open the file
52//
53// The final root:Packages:NIST:RAW:data wave is the real
54//neutron counts and can be directly operated on
55//
56// header information is put into three waves: integersRead, realsRead, and textRead
57// logicals in the header are currently skipped, since they are no use in the
58// data reduction process.
59//
60// The output is the three R/T/I waves that are filled at least with minimal values
61// and the detector data loaded into an array named "data"
62//
63// see documentation for the expected information in each element of the R/T/I waves
64// and the minimum set of information. These waves can be increased in length so that
65// more information can be accessed as needed (propagating changes...)
66//
67// called by multiple .ipfs (when the file name is known)
68//
69//
70// THIS FUNCTION DOES NOT NEED TO BE MODIFIED. ONLY THE DATA ACCESSORS NEED TO BE CONSTRUCTED
71//
72Function ReadHeaderAndData(fname)
73        String fname
74        //this function is for reading in RAW data only, so it will always put data in RAW folder
75        String curPath = "root:Packages:NIST:RAW"
76        SetDataFolder curPath           //use the full path, so it will always work
77        variable/g root:Packages:NIST:RAW:gIsLogScale = 0               //initial state is linear, keep this in RAW folder
78       
79        Variable refNum,integer,realval
80        String sansfname,textstr
81       
82        Make/D/O/N=23 $"root:Packages:NIST:RAW:IntegersRead"
83        Make/D/O/N=52 $"root:Packages:NIST:RAW:RealsRead"
84        Make/O/T/N=11 $"root:Packages:NIST:RAW:TextRead"
85        Make/O/N=7 $"root:Packages:NIST:RAW:LogicalsRead"
86
87        Wave intw=$"root:Packages:NIST:RAW:IntegersRead"
88        Wave realw=$"root:Packages:NIST:RAW:RealsRead"
89        Wave/T textw=$"root:Packages:NIST:RAW:TextRead"
90        Wave logw=$"root:Packages:NIST:RAW:LogicalsRead"
91       
92        // FILL IN 3 ARRAYS WITH HEADER INFORMATION FOR LATER USE
93        // THESE ARE JUST THE MINIMALLY NECESSARY VALUES
94       
95        // filename as stored in the file header
96        textw[0]= ParseFilePath(0, fname, ":", 1, 0)   
97       
98        // date and time of collection
99        textw[1]= getFileCreationDate(fname)
100       
101        // run type identifier (this is a reader for RAW data)
102        textw[2]= "RAW"
103       
104        // user account identifier (currently used only for NCNR-specific operations)
105        textw[3]= ""
106
107        // sample label
108        textw[6]= getSampleLabel(fname)
109       
110        // identifier of detector type, useful for setting detector constants
111        //(currently used only for NCNR-specific operations)
112        textw[9]= ""
113
114        //total counting time in seconds
115        intw[2] = getCountTime(fname)
116       
117//      realw[1] = 0
118//      realw[6] = 0
119//      realw[7] = 0
120//      realw[8] = 0
121//      realw[9] = 0
122//      realw[19] = 0
123//      realw[22] = 0
124
125        // total monitor count
126        realw[0] = getMonitorCount(fname)
127
128       
129        // attenuator number (NCNR-specific, your stub returns 0)
130        // may also be used to hold attenuator transmission (< 1)
131        realw[3] = getAttenNumber(fname)
132       
133        // sample transmission
134        realw[4] = getSampleTrans(fname)
135       
136        //sample thickness (cm)
137        realw[5] = getSampleThickness(fname)
138       
139        // following 6 values are for non-linear spatial corrections to a detector (RC timing)
140        // these values are set for a detctor that has a linear correspondence between
141        // pixel number and real-space distance
142        // 10 and 13 are the X and Y pixel dimensions, respectively (in mm!)
143        //(11,12 and 13,14 are set to values for a linear response, as from a new Ordela detector)
144        realw[10] = getDetectorPixelXSize(fname)
145        realw[11] = 10000
146        realw[12] = 0
147        realw[13] = getDetectorPixelYSize(fname)
148        realw[14] = 10000
149        realw[15] = 0
150       
151        // beam center X,Y on the detector (in units of pixel coordinates (1,N))
152        realw[16] = getBeamXPos(fname)
153        realw[17] = getBeamYPos(fname)
154       
155        // sample to detector distance (meters)
156        realw[18] = getSDD(fname)
157
158        // detector physical width (right now assumes square...) (in cm)
159        //realw[20] = 65
160        realw[20] = getPhysicalWidth(fname)
161       
162        // beam stop diameter (assumes circular) (in mm)
163        realw[21] = getBSDiameter(fname)
164       
165        // source aperture diameter (mm)
166        realw[23] = getSourceApertureDiam(fname)
167       
168        // sample aperture diameter (mm)
169        realw[24] = getSampleApertureDiam(fname)
170       
171        // source aperture to sample aperture distance
172        realw[25] = getSourceToSampleDist(fname)
173       
174        // wavelength (A)
175        realw[26] = getWavelength(fname)
176       
177        // wavelength spread (FWHM)
178        realw[27] = getWavelengthSpread(fname)
179       
180        // beam stop X-position (motor reading, approximate cm from zero position)
181        // currently NCNR-specific use to identify transmission measurements
182        // you return 0
183        realw[37] = 0
184
185// the actual data array, dimensions are set as globals in
186// InitFacilityGlobals()
187        NVAR XPix = root:myGlobals:gNPixelsX
188        NVAR YPix = root:myGlobals:gNPixelsX
189       
190        Make/D/O/N=(XPix,YPix) $"root:Packages:NIST:RAW:data"
191        WAVE data=$"root:Packages:NIST:RAW:data"
192
193        // fill the data array with the detector values
194        getDetectorData(fname,data)
195       
196        // total detector count
197        //nha 21/5/10 moved here because it requires the detector data to already be written
198        //Result of issue with 0 counts being written for a while in metadata.
199        realw[2] = getDetCount(fname)
200       
201        //keep a string with the filename in the RAW folder
202        String/G root:Packages:NIST:RAW:fileList = textw[0]
203       
204        Duplicate/O data linear_data                    //data read in is on linear scale, copy it now
205
206        // proper error for counting statistics, good for low count values too
207        // rather than just sqrt(n)
208        // see N. Gehrels, Astrophys. J., 303 (1986) 336-346, equation (7)
209        // for S = 1 in eq (7), this corresponds to one sigma error bars
210        Duplicate/O linear_data linear_data_error
211        linear_data_error = 1 + sqrt(linear_data + 0.75)                               
212        //
213       
214        Return 0
215
216End
217
218//****************
219//main entry procedure for reading a "WORK.DIV" file
220//displays a quick image of the  file, to check that it's correct
221//data is deposited in root:Packages:NIST:DIV data folder
222//
223// local, just for testing
224//
225Proc ReadWork_DIV()
226       
227        String fname = PromptForPath("Select detector sensitivity file")
228        ReadHeaderAndWork("DIV",fname)          //puts what is read in work.div
229       
230        String waveStr = "root:Packages:NIST:DIV:data"
231        NewImage/F/K=1/S=2 $waveStr
232        ModifyImage '' ctab= {*,*,YellowHot,0}
233       
234        String/G root:Packages:NIST:DIV:fileList = "WORK.DIV"
235       
236        SetDataFolder root:             //(redundant)
237End
238
239
240
241// Detector sensitivity files have the same header structure as RAW SANS data
242// as NCNR, just with a different data array (real, rather than integer data)
243//
244// So for your facility, make this reader specific to the format of whatever
245// file you use for a pixel-by-pixel division of the raw detector data
246// to correct for non-uniform sensitivities of the detector. This is typically a
247// measurement of water, plexiglas, or another uniform scattering sample.
248//
249// The data must be normalized to a mean value of 1
250//
251// called from ProtocolAsPanel.ipf
252//
253// type is "DIV" on input
254Function ReadHeaderAndWork(type,fname)
255        String type,fname
256       
257        //type is the desired folder to read the workfile to
258        //this data will NOT be automatically displayed
259        // gDataDisplayType is unchanged
260
261        String cur_folder = type
262        String curPath = "root:Packages:NIST:"+cur_folder
263        SetDataFolder curPath           //use the full path, so it will always work
264       
265        Variable refNum,integer,realval
266        String sansfname,textstr
267        Variable/G $(curPath + ":gIsLogScale") = 0              //initial state is linear, keep this in DIV folder
268       
269        Make/O/D/N=23 $(curPath + ":IntegersRead")
270        Make/O/D/N=52 $(curPath + ":RealsRead")
271        Make/O/T/N=11 $(curPath + ":TextRead")
272       
273        WAVE intw=$(curPath + ":IntegersRead")
274        WAVE realw=$(curPath + ":RealsRead")
275        WAVE/T textw=$(curPath + ":TextRead")
276       
277        // the actual data array, dimensions are set as globals in
278        // InitFacilityGlobals()
279        NVAR XPix = root:myGlobals:gNPixelsX
280        NVAR YPix = root:myGlobals:gNPixelsX
281
282        Make/O/D/N=(XPix,YPix) $(curPath + ":data")
283        WAVE data = $(curPath + ":data")
284       
285        // (1)
286        // fill in your reader for the header here so that intw, realw, and textW are filled in
287        // ? possibly a duplication of the raw data reader
288               
289        //(2)
290        // then fill in a reader for the data array that will DIVIDE your data
291        // to get the corrected values.
292        String dfName=""
293        variable err
294        err = hdfRead(fname, dfName)
295        Wave tempData = $dfName+":data:div"
296        data = tempData
297       
298        //funky edge correction bodgy ???
299        //copy second column to first column
300        data[][0] = data[p][1]
301        //copy second last column to last column
302        data[][191] = data[p][190]
303        //copy second row to first row
304        data[0][] = data[1][q]
305        //copy second last row to last row
306        data[191][] = data[190][q]
307        //keep a string with the filename in the DIV folder
308        String/G $(curPath + ":fileList") = textw[0]
309       
310        Duplicate/O data linear_data                    //data read in is on linear scale, copy it now
311
312        //return the data folder to root
313        SetDataFolder root:
314       
315        Return(0)
316End
317
318Function WriteHeaderAndWork(type)
319        String type
320       
321        // your writer here
322        NVAR XPix = root:myGlobals:gNPixelsX
323        NVAR YPix = root:myGlobals:gNPixelsX   
324       
325        Wave wData=$("root:Packages:NIST:"+type+":data")
326       
327//      Variable refnum,ii=0,hdrBytes=516,a,b,offset
328        String fname=""
329//      Duplicate/O wData,tempData
330       
331        //changed for Quokka
332//      Redimension/S/N=(XPix*YPix) tempData
333//      tempData *= 4
334       
335        PathInfo/S catPathName
336        fname = DoSaveFileDialog("Save data as")          //won't actually open the file
337        If(cmpstr(fname,"")==0)
338                //user cancel, don't write out a file
339          Close/A
340          Abort "no data file was written"
341        Endif
342       
343        variable fileID
344        HDF5CreateFile /O fileID as fname
345       
346        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
347        //Make /N=(1,1) wTransmission
348        String groupName = "/reduce"
349        String varName = "div"
350        // your code returning value
351        variable err
352        err = hdfWrite(fname, groupName, varName, wData)
353
354        // KillWaves wData, tempData
355        return(0)
356End
357
358
359
360/////   ASC FORMAT READER  //////
361/////   FOR WORKFILE MATH PANEL //////
362//
363//function to read in the ASC output of SANS reduction
364// currently the file has 20 header lines, followed by a single column
365// of 16384 values, Data is written by row, starting with Y=1 and X=(1->128)
366//
367//returns 0 if read was ok
368//returns 1 if there was an error
369//
370// called by WorkFileUtils.ipf
371//
372//
373// If the ASC data was generated by the NCNR data writer, then
374// NOTHING NEEDS TO BE CHANGED HERE
375//
376Function ReadASCData(fname,destPath)
377        String fname, destPath
378        //this function is for reading in ASCII data so put data in user-specified folder
379        SetDataFolder "root:Packages:NIST:"+destPath
380
381        NVAR pixelsX = root:myGlobals:gNPixelsX
382        NVAR pixelsY = root:myGlobals:gNPixelsY
383        Variable refNum=0,ii,p1,p2,tot,num=pixelsX,numHdrLines=20
384        String str=""
385        //data is initially linear scale
386        Variable/G :gIsLogScale=0
387        Make/O/T/N=(numHdrLines) hdrLines
388        Make/O/D/N=(pixelsX*pixelsY) data                       //,linear_data
389       
390        //full filename and path is now passed in...
391        //actually open the file
392//      SetDataFolder destPath
393        Open/R/Z refNum as fname                // /Z flag means I must handle open errors
394        if(refnum==0)           //FNF error, get out
395                DoAlert 0,"Could not find file: "+fname
396                Close/A
397                SetDataFolder root:
398                return(1)
399        endif
400        if(V_flag!=0)
401                DoAlert 0,"File open error: V_flag="+num2Str(V_Flag)
402                Close/A
403                SetDataFolder root:
404                return(1)
405        Endif
406        //
407        for(ii=0;ii<numHdrLines;ii+=1)          //read (or skip) 18 header lines
408                FReadLine refnum,str
409                hdrLines[ii]=str
410        endfor
411        //     
412        Close refnum
413       
414//      SetDataFolder destPath
415        LoadWave/Q/G/D/N=temp fName
416        Wave/Z temp0=temp0
417        data=temp0
418        Redimension/N=(pixelsX,pixelsY) data            //,linear_data
419
420        Duplicate/O data linear_data_error
421        linear_data_error = 1 + sqrt(data + 0.75)
422       
423        //just in case there are odd inputs to this, like negative intensities
424        WaveStats/Q linear_data_error
425        linear_data_error = numtype(linear_data_error[p]) == 0 ? linear_data_error[p] : V_avg
426        linear_data_error = linear_data_error[p] != 0 ? linear_data_error[p] : V_avg
427       
428        //linear_data = data
429       
430        KillWaves/Z temp0
431       
432        //return the data folder to root
433        SetDataFolder root:
434       
435        Return(0)
436End
437
438// fills the "default" fake header so that the SANS Reduction machinery does not have to be altered
439// pay attention to what is/not to be trusted due to "fake" information.
440// uses what it can from the header lines from the ASC file (hdrLines wave)
441//
442// destFolder is of the form "myGlobals:WorkMath:AAA"
443//
444//
445// called by WorkFileUtils.ipf
446//
447// If the ASC data was generated by the NCNR data writer, then
448// NOTHING NEEDS TO BE CHANGED HERE
449//
450Function FillFakeHeader_ASC(destFolder)
451        String destFolder
452        Make/O/D/N=23 $("root:Packages:NIST:"+destFolder+":IntegersRead")
453        Make/O/D/N=52 $("root:Packages:NIST:"+destFolder+":RealsRead")
454        Make/O/T/N=11 $("root:Packages:NIST:"+destFolder+":TextRead")
455       
456        Wave intw=$("root:Packages:NIST:"+destFolder+":IntegersRead")
457        Wave realw=$("root:Packages:NIST:"+destFolder+":RealsRead")
458        Wave/T textw=$("root:Packages:NIST:"+destFolder+":TextRead")
459       
460        //Put in appropriate "fake" values
461        //parse values as needed from headerLines
462        Wave/T hdr=$("root:Packages:NIST:"+destFolder+":hdrLines")
463        Variable monCt,lam,offset,sdd,trans,thick
464        Variable xCtr,yCtr,a1,a2,a1a2Dist,dlam,bsDiam
465        String detTyp=""
466        String tempStr="",formatStr="",junkStr=""
467        formatStr = "%g %g %g %g %g %g"
468        tempStr=hdr[3]
469        sscanf tempStr, formatStr, monCt,lam,offset,sdd,trans,thick
470//      Print monCt,lam,offset,sdd,trans,thick,avStr,step
471        formatStr = "%g %g %g %g %g %g %g %s"
472        tempStr=hdr[5]
473        sscanf tempStr,formatStr,xCtr,yCtr,a1,a2,a1a2Dist,dlam,bsDiam,detTyp
474//      Print xCtr,yCtr,a1,a2,a1a2Dist,dlam,bsDiam,detTyp
475       
476        realw[16]=xCtr          //xCtr(pixels)
477        realw[17]=yCtr  //yCtr (pixels)
478        realw[18]=sdd           //SDD (m)
479        realw[26]=lam           //wavelength (A)
480        //
481        // necessary values
482        realw[10]=5                     //detector calibration constants, needed for averaging
483        realw[11]=10000
484        realw[12]=0
485        realw[13]=5
486        realw[14]=10000
487        realw[15]=0
488        //
489        // used in the resolution calculation, ONLY here to keep the routine from crashing
490        realw[20]=65            //det size
491        realw[27]=dlam  //delta lambda
492        realw[21]=bsDiam        //BS size
493        realw[23]=a1            //A1
494        realw[24]=a2    //A2
495        realw[25]=a1a2Dist      //A1A2 distance
496        realw[4]=trans          //trans
497        realw[3]=0              //atten
498        realw[5]=thick          //thick
499        //
500        //
501        realw[0]=monCt          //def mon cts
502
503        // fake values to get valid deadtime and detector constants
504        //
505        textw[9]=detTyp+"  "            //6 characters 4+2 spaces
506        textw[3]="[NGxSANS00]"  //11 chars, NGx will return default values for atten trans, deadtime...
507       
508        //set the string values
509        formatStr="FILE: %s CREATED: %s"
510        sscanf hdr[0],formatStr,tempStr,junkStr
511//      Print tempStr
512//      Print junkStr
513        String/G $("root:Packages:NIST:"+destFolder+":fileList") = tempStr
514        textw[0] = tempStr              //filename
515        textw[1] = junkStr              //run date-time
516       
517        //file label = hdr[1]
518        tempStr = hdr[1]
519        tempStr = tempStr[0,strlen(tempStr)-2]          //clean off the last LF
520//      Print tempStr
521        textW[6] = tempStr      //sample label
522       
523        return(0)
524End
525
526
527///////// ACCESSORS FOR WRITING DATA TO HEADER  /////////
528/////
529
530// Write* routines all must:
531// (1) open file "fname". fname is a full file path and name to the file on disk
532// (2) write the specified value to the header at the correct location in the file
533// (3) close the file
534
535// new, April 2011 for error propagation. fill these in with the facility-
536// specific versions, if desired.
537Function WriteTransmissionErrorToHeader(fname,transErr)
538        String fname
539        Variable transErr
540       
541
542        return(0)
543End
544
545Function WriteBoxCountsErrorToHeader(fname,rel_err)
546        String fname
547        Variable rel_err
548       
549        return(0)
550End
551
552Function getSampleTransError(fname)
553        String fname
554       
555        return(0)
556end
557
558Function getBoxCountsError(fname)
559        String fname
560       
561        return(0)
562end
563
564
565// end April 2011 additions
566
567// read the detector deadtime (in seconds)
568Function getDetectorDeadtime(fname)
569        String fname
570       
571        return(0)
572end
573
574
575//whole transmission is NCNR-specific right now
576// leave this stub empty
577Function WriteWholeTransToHeader(fname,wholeTrans)
578        String fname
579        Variable wholeTrans
580       
581        String groupName = "/reduce"
582        variable err
583       
584        Wave wCounts
585        Make /N=(1,1) wWholeTrans
586               
587        wWholeTrans[0] = wholeTrans
588       
589        String varName = "wholeTransmission"   
590        err = hdfWrite(fname, groupName, varName,wWholeTrans)
591
592        KillWaves wWholeTrans
593       
594        //err not handled here
595        return(0)       
596End
597
598//box sum counts is a real value
599// used for transmission calculation module
600Function WriteBoxCountsToHeader(fname,counts)
601        String fname
602        Variable counts
603       
604        // do nothing if not using NCNR Transmission module
605       
606        String groupName = "/reduce"
607        variable err
608       
609        Wave wCounts
610        Make /N=(1,1) wCounts
611               
612        wCounts[0] = counts
613       
614        String varName = "boxCounts"   
615        err = hdfWrite(fname, groupName, varName,wCounts)
616
617        KillWaves wCounts
618       
619        //err not handled here
620        return(0)       
621End
622
623//beam stop X-position
624// used for transmission module to manually tag transmission files
625Function WriteBSXPosToHeader(fname,xpos)
626        String fname
627        Variable xpos
628       
629        // do nothing if not using NCNR Transmission module
630       
631        return(0)
632End
633
634
635
636
637
638//attenuator number (not its transmission)
639// if your beam attenuation is indexed in some way, use that number here
640// if not, write a 1 to the file here as a default
641//
642Function WriteAttenNumberToHeader(fname,attenNumber)
643        String fname
644        Variable attenNumber
645       
646        // your writer here
647        Wave wAttenNumber
648        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
649        Make /N=(1,1) wAttenNumber
650        String groupName = "/instrument/collimator"
651        String varName = "att"
652        //convert number to a rotation angle
653        attenNumber = attenNumber * 30
654        wAttenNumber[0] = attenNumber //
655        // your code returning value
656        variable err
657        err = hdfWrite(fname, groupName, varName, wAttenNumber)
658        KillWaves wAttenNumber
659        //err not handled here
660               
661        return(0)
662
663End
664
665
666// total monitor count during data collection
667Function WriteMonitorCountToHeader(fname,num)
668        String fname
669        Variable num
670       
671        // your code here
672       
673        return(0)
674End
675
676//total detector count
677Function WriteDetectorCountToHeader(fname,num)
678        String fname
679        Variable num
680       
681        // your code here
682       
683        return(0)
684End
685
686//transmission detector count
687// (currently unused in data reduction)
688Function WriteTransDetCountToHeader(fname,num)
689        String fname
690        Variable num
691       
692        // do nothing for now
693       
694        return(0)
695End
696
697
698
699//temperature of the sample (C)
700Function WriteTemperatureToHeader(fname,num)
701        String fname
702        Variable num
703       
704        //  your code here
705       
706        return(0)
707End
708
709//magnetic field (Oe)
710Function WriteMagnFieldToHeader(fname,num)
711        String fname
712        Variable num
713       
714        // your code here
715       
716        return(0)
717End
718
719//lateral detector offset (centimeters)
720Function WriteDetectorOffsetToHeader(fname,DetectorOffset)
721        String fname
722        Variable DetectorOffset
723       
724        // your writer here
725        Wave wDetectorOffset
726        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
727        Make /N=(1,1) wDetectorOffset
728        String groupName = "/instrument/detector"
729        String varName = "detector_x"
730        //convert from cm (NIST standard) to mm (NeXus standard)
731        DetectorOffset = DetectorOffset * 10
732        wDetectorOffset[0] = DetectorOffset //
733        // your code returning value
734        variable err
735        err = hdfWrite(fname, groupName, varName, wDetectorOffset)
736        KillWaves wDetectorOffset
737        //err not handled here
738       
739        return(0)
740End
741
742
743
744
745// physical dimension of detector x-pixel (mm)
746Function WriteDetPixelXToHeader(fname,num)
747        String fname
748        Variable num
749       
750        //your code here
751       
752        return(0)
753end
754
755// physical dimension of detector y-pixel (mm)
756Function WriteDetPixelYToHeader(fname,num)
757        String fname
758        Variable num
759       
760        //your code here
761       
762        return(0)
763end
764
765// sample label string
766Function WriteSamLabelToHeader(fname,str)
767        String fname,str
768       
769        // your code here
770
771        return(0)
772End
773
774// total counting time (seconds)
775Function WriteCountTimeToHeader(fname,num)
776        String fname
777        Variable num
778       
779        // your code here
780       
781        return(0)
782End
783
784// Write the detector deadtime to the file header (in seconds)
785Function WriteDeadtimeToHeader(fname,num)
786        String fname
787        Variable num
788       
789        return(0)
790End
791
792
793//////// ACCESSORS FOR READING DATA FROM THE HEADER  //////////////
794//
795// read specific bits of information from the header
796// each of these operations MUST take care of open/close on their own
797//
798// fname is the full path and filname to the file on disk
799// return values are either strings or real values as appropriate
800//
801//////
802
803
804// function that reads in the 2D detector data and fills the array
805// data[nx][ny] with the data values
806// fname is the full name and path to the data file for opening and closing
807//
808//
809Function getDetectorData(fname,data)
810        String fname
811        Wave data
812        NVAR XPix = root:myGlobals:gNPixelsX
813       
814        // your reader here
815        variable err
816        string dfName = ""
817        err = hdfRead(fname, dfName)
818        if(err)
819                return 0
820        endif
821
822        Wave hmm_xy = $(dfName+":data:hmm_xy")
823       
824        //redimension /I /N = (dimsize(hmm_xy, 2), dimsize(hmm_xy, 1)), data
825        //nha. Count arrays need to be floating point, since the data will be divided, normalised etc.
826        redimension /N = (dimsize(hmm_xy, 2), dimsize(hmm_xy, 1)), data
827        data[][] = hmm_xy[0][q][p]
828       
829        //nha workaround. for wrongly dimensioned Quokka data 191x192
830        variable x_dim = dimsize(data,0)
831        if (x_dim!=XPix)
832                //redimension to add an extra row(s) to the data
833                redimension /I /N = (XPix,-1) data
834                //copy row 190 to row 191
835                data[191][] = data[190][q]
836        endif   
837        // end workaround
838       
839        //funky edge correction bodgy ???
840        //copy second column to first column
841        data[][0] = data[p][1]
842        //copy second last column to last column
843        data[][191] = data[p][190]
844        //copy second row to first row
845        data[0][] = data[1][q]
846        //copy second last row to last row
847        data[191][] = data[190][q]
848                       
849        KillWaves hmm_xy
850       
851       
852        return(0)
853End
854
855// file suffix (NCNR data file name specific)
856// return filename as suffix
857Function/S getSuffix(fname)
858        String fname
859       
860        return(ParseFilePath(0, fname, ":", 1, 0))
861End
862
863// associated file suffix (for transmission)
864// NCNR Transmission calculation specific
865// return null string
866Function/S getAssociatedFileSuffix(fname)
867        String fname
868       
869        return(getFileAssoc(fname))
870End
871
872// sample label
873Function/S getSampleLabel(fname)
874        String fname
875        String str
876       
877        // your code, returning str
878        variable err
879        string dfName = ""
880        err = hdfRead(fname, dfName)
881        //err not handled here
882       
883        Wave/T wSampleName = $dfname+":sample:name"
884        str = wSampleName[0]
885        KillWaves wSampleName
886       
887        return(str)
888End
889
890// file creation date
891Function/S getFileCreationDate(fname)
892        String fname
893        String str
894       
895        // your code, returning str
896        variable err
897        string dfName = ""
898        err = hdfRead(fname, dfName)
899        //err not handled here
900       
901        Wave/T wStartTime = $dfName+":start_time"
902        str = wStartTime[0]
903        KillWaves wStartTime
904       
905        return(str)
906End
907
908
909//monitor count
910Function getMonitorCount(fname)
911//not patched
912        String fname
913        Variable value
914       
915        // your code returning value
916        variable err
917        string dfName = ""
918        err = hdfRead(fname, dfName)
919        //err not handled here
920       
921        Wave wCounts = $dfName+":monitor:bm1_counts"
922        value = wCounts[0]
923        KillWaves wCounts
924       
925        return(value)
926end
927
928//saved monitor count
929// NCNR specific for pre-processed data, never used
930// return 0
931Function getSavMon(fname)
932        String fname
933       
934        Variable value
935       
936        // your code returning value
937        //??? to do. Is this required if 'never used'? nha
938       
939        return(0)
940end
941
942//total detector count
943Function getDetCount(fname)
944//not patched, but could be
945        String fname
946        Variable value
947       
948        // your code returning value
949        variable err
950        string dfName = ""
951        err = hdfRead(fname, dfName)
952        //err not handled here
953
954        //      Wave wCounts = $(dfName+":data:total_counts")
955        // changed 22/12/09 nha
956                if(WaveExists($(dfName+":data:total_counts")))
957                        Wave wCounts = $(dfName+":data:total_counts")
958                elseif(WaveExists($(dfName+":instrument:detector:total_counts")))
959                Wave wCounts = $(dfName+":instrument:detector:total_counts")
960               else
961                print "Can't find detector total_counts in " + fname
962               endif
963       
964        value = wCounts[0]
965       
966        //nha 21/5/10 temporary glitch wrote detector count to file as 0       
967                if (value<1)
968                        NVAR XPix = root:myGlobals:gNPixelsX
969                        NVAR YPix = root:myGlobals:gNPixelsX
970                        Make/D/O/N=(XPix,YPix) $"root:Packages:NIST:RAW:data"
971                        WAVE data=$"root:Packages:NIST:RAW:data"
972                        getDetectorData(fname,data)
973                        value = sum(data)
974                endif
975       
976        KillWaves wCounts
977       
978        return(value)
979end
980
981//Attenuator number, return 1 if your attenuators are not
982// indexed by number
983Function getAttenNumber(fname)
984        String fname
985        Variable value
986        Variable att, err
987        string dfName = ""
988               
989        err = hdfRead(fname, dfName)
990        //err not handled here
991
992        if(WaveExists($(dfName+":instrument:collimator:att")))
993                        Wave wAttrotdeg = $(dfName+":instrument:collimator:att")
994        elseif(WaveExists($(dfName+":instrument:parameters:derived_parameters:AttRotDeg")))
995                        Wave wAttrotdeg = $(dfName+":instrument:parameters:derived_parameters:AttRotDeg")
996        else
997                        print "Can't find attenuator in " + fname
998        endif   
999       
1000        value = wAttrotdeg[0]
1001        att = round(value)/30
1002        KillWaves wAttrotdeg
1003        return(att)
1004end
1005
1006
1007//transmission
1008Function getSampleTrans(fname)
1009        String fname
1010       
1011        Variable value
1012       
1013        // your code returning value
1014        variable err
1015        string dfName = ""
1016        err = hdfRead(fname, dfName)
1017        //err not handled here
1018
1019        if(WaveExists($(dfName+":reduce:Transmission"))) //canonical location
1020                        Wave wTransmission = $(dfName+":reduce:Transmission")
1021        elseif(WaveExists($(dfName+":instrument:parameters:Transmission")))
1022                        Wave wTransmission = $(dfName+":instrument:parameters:Transmission")
1023        else
1024                        print "Can't find Transmission in " + fname
1025        endif
1026        value = wTransmission[0]
1027        KillWaves wTransmission
1028               
1029        return(value)
1030end
1031
1032//sample transmission (0<T<=1)
1033Function WriteTransmissionToHeader(fname,trans)
1034        String fname
1035        Variable trans
1036       
1037        // your writer here
1038//      Wave wTransmission
1039        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1040//      Make /N=(1,1) wTransmission
1041        Make/O/D/N=1 wTransmission
1042        String groupName = "/reduce"
1043        String varName = "Transmission"
1044        wTransmission[0] = trans //
1045        // your code returning value
1046        variable err
1047        err = hdfWrite(fname, groupName, varName, wTransmission)
1048        KillWaves wTransmission
1049        //err not handled here
1050               
1051        return(0)
1052End
1053
1054//box counts from stored transmission calculation
1055// return 0 if not using NCNR transmission module
1056Function getBoxCounts(fname)
1057        String fname
1058        Variable value
1059       
1060        // do nothing if not using NCNR Transmission module
1061        variable err
1062        string dfName = ""
1063        err = hdfRead(fname, dfName)
1064        //err not handled here
1065
1066        Wave wBoxCounts = $(dfName+":reduce:boxCounts")
1067        if (waveexists(wBoxCounts) == 0)
1068                //boxcounts not yet set in  reduce group
1069                //return 0
1070                value = 0
1071        else
1072                value = wBoxCounts[0]
1073        endif
1074
1075        KillWaves/Z wBoxCounts
1076       
1077        return(value)
1078end
1079
1080//whole detector transmission
1081// return 0 if not using NCNR transmission module
1082Function getSampleTransWholeDetector(fname)
1083        String fname
1084        Variable value
1085       
1086        // your code returning value
1087        // ??? don't know what to put here. nha
1088        value=0
1089        return(value)
1090end
1091
1092//SampleThickness in centimeters
1093Function getSampleThickness(fname)
1094        String fname
1095        Variable value
1096
1097        // your code returning value
1098       
1099        variable err
1100        string dfName = ""
1101        err = hdfRead(fname, dfName)
1102        //err not handled here
1103
1104        if(WaveExists($(dfName+":sample:SampleThickness"))) //canonical location - a bit ugly and verbose, but that's just my opinion
1105                Wave wThickness = $(dfName+":sample:SampleThickness")
1106        elseif(WaveExists($(dfName+":sample:thickness")))
1107                Wave wThickness = $(dfName+":sample:thickness")
1108        else
1109                print "Can't find Sample Thickness in " + fname
1110        endif
1111                       
1112        value = wThickness[0]/10
1113        //value = 1 //??? temporary fix. nha
1114        KillWaves wThickness
1115       
1116        return(value)
1117end
1118
1119//sample thickness in cm
1120Function WriteThicknessToHeader(fname,thickness)
1121        String fname
1122        Variable thickness
1123       
1124        // your writer here
1125        Wave wThickness
1126        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1127        Make /N=(1,1) wThickness
1128        String groupName = "/sample"
1129        String varName = "SampleThickness"
1130        wThickness[0] = thickness*10 //
1131        // your code returning value
1132        variable err
1133        err = hdfWrite(fname, groupName, varName, wThickness) //does not exist ???
1134        KillWaves wThickness
1135        //err not handled here
1136               
1137        return(0)
1138End
1139
1140//Sample Rotation Angle (degrees)
1141Function getSampleRotationAngle(fname)
1142        String fname
1143        Variable value
1144       
1145        // your code returning value
1146        variable err
1147        string dfName = ""
1148        err = hdfRead(fname, dfName)
1149        //err not handled here
1150
1151        Wave wSample_rotation_angle = $(dfName+":sample:sample_theta") //is this correct
1152        value = wSample_rotation_angle[0]
1153        KillWaves wSample_rotation_angle
1154               
1155        return(value)
1156end
1157
1158//temperature (C)
1159Function getTemperature(fname)
1160        String fname
1161       
1162        Variable value
1163       
1164        // your code returning value
1165       
1166        return(value)
1167end
1168
1169//field strength (Oe)
1170Function getFieldStrength(fname)
1171        String fname
1172       
1173        Variable value
1174       
1175        // your code returning value
1176       
1177        return(value)
1178end
1179
1180//center of beam xPos in pixel coordinates
1181Function getBeamXPos(fname)
1182        String fname
1183        Variable value
1184        // your code returning value
1185        variable err
1186        string dfName = ""
1187        err = hdfRead(fname, dfName)
1188        //err not handled here
1189
1190        if(WaveExists($(dfName+":instrument:reduce:BeamCenterX"))) //canonical location
1191                        Wave wBeamXPos = $(dfName+":instrument:reduce:BeamCenterX")
1192        elseif(WaveExists($(dfName+":instrument:parameters:BeamCenterX")))
1193                        Wave wBeamXPos = $(dfName+":instrument:parameters:BeamCenterX")
1194        else
1195                        print "Can't find BeamCenterX in" $fname
1196        endif
1197        value = wBeamXPos[0]   
1198        KillWaves wBeamXPos
1199       
1200        return(value)   
1201
1202end
1203
1204//beam center X pixel location (detector coordinates)
1205Function WriteBeamCenterXToHeader(fname,beamCenterX)
1206        String fname
1207        Variable beamCenterX
1208       
1209        // your writer here
1210        Wave wBeamCenterX
1211        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1212        Make /N=(1,1) wBeamCenterX
1213        String groupName = "/instrument/reduce"
1214        String varName = "BeamCenterX"
1215        wBeamCenterX[0] = beamCenterX //
1216        // your code returning value
1217        variable err
1218        err = hdfWrite(fname, groupName, varName, wBeamCenterX)
1219        KillWaves wBeamCenterX
1220        //err not handled here
1221       
1222        return(0)
1223End
1224
1225//center of beam Y pos in pixel coordinates
1226Function getBeamYPos(fname)
1227        String fname
1228        Variable value
1229        // your code returning value
1230        variable err
1231        string dfName = ""
1232        err = hdfRead(fname, dfName)
1233        //err not handled here
1234
1235        if(WaveExists($(dfName+":instrument:reduce:BeamCenterZ"))) //canonical location
1236                        Wave wBeamYPos = $(dfName+":instrument:reduce:BeamCenterZ")
1237        elseif(WaveExists($(dfName+":instrument:parameters:BeamCenterZ")))
1238                        Wave wBeamYPos = $(dfName+":instrument:parameters:BeamCenterZ")
1239        else
1240                        print "Can't find BeamCenterZ in" $fname
1241        endif
1242        value = wBeamYPos[0]   
1243        KillWaves wBeamYPos
1244               
1245        return(value)
1246end
1247
1248//beam center Y pixel location (detector coordinates)
1249Function WriteBeamCenterYToHeader(fname,beamCenterY)
1250        String fname
1251        Variable beamCenterY
1252       
1253        // your writer here
1254        Wave wBeamCenterY
1255        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1256        Make /N=(1,1) wBeamCenterY
1257        String groupName = "/instrument/reduce"
1258        String varName = "BeamCenterZ"
1259        wBeamCenterY[0] = beamCenterY //
1260        // your code returning value
1261        variable err
1262        err = hdfWrite(fname, groupName, varName, wBeamCenterY)
1263        KillWaves wBeamCenterY
1264        //err not handled here
1265
1266        return(0)
1267End
1268
1269
1270//sample to detector distance (meters)
1271Function getSDD(fname)
1272        String fname
1273        Variable value
1274        // your code returning value
1275        variable err
1276        string dfName = ""
1277        err = hdfRead(fname, dfName)
1278        //err not handled here
1279
1280        //workaround for bad HDF5 dataset
1281        if(WaveExists($(dfName+":instrument:parameters:L2"))) //canonical location
1282                Wave wSourceToDetectorDist = $(dfName+":instrument:parameters:L2")
1283        elseif(WaveExists($(dfName+":instrument:parameters:L2mm")))
1284                Wave wSourceToDetectorDist = $(dfName+":instrument:parameters:L2mm")
1285        elseif(WaveExists($(dfName+":instrument:parameters:derived_parameters:L2mm"))) 
1286                Wave wSourceToDetectorDist = $(dfName+":instrument:parameters:derived_parameters:L2mm")
1287        else
1288                print "Can't find L2 in " + fname
1289        endif
1290       
1291        value = wSourceToDetectorDist[0]/1000   
1292        KillWaves wSourceToDetectorDist
1293               
1294        return(value)
1295end
1296
1297//sample to detector distance (meters)
1298Function WriteSDDToHeader(fname,sdd)
1299        String fname
1300        Variable sdd
1301       
1302// your writer here
1303        Wave wSDD
1304        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1305        Make /N=(1,1) wSDD
1306        String groupName = "/instrument/parameters"
1307        String varName = "L2"
1308        wSDD[0] = sdd * 1000 //
1309        // your code returning value
1310        variable err
1311        err = hdfWrite(fname, groupName, varName, wSDD)
1312        KillWaves wSDD
1313        //err not handled here
1314       
1315        return(0)
1316End
1317
1318//lateral detector offset (centimeters)
1319Function getDetectorOffset(fname)
1320        String fname
1321        Variable value
1322        // your code returning value
1323        variable err
1324        string dfName = ""
1325        err = hdfRead(fname, dfName)
1326        //err not handled here
1327
1328        Wave wDetectorOffset = $(dfName+":instrument:detector:detector_x") //is this correct
1329        value = wDetectorOffset[0]/10
1330        KillWaves wDetectorOffset
1331       
1332        return(value)
1333end
1334
1335//Beamstop diameter (millimeters)
1336Function getBSDiameter(fname)
1337        String fname
1338        Variable value
1339        // your code returning value
1340        variable err
1341        string dfName = ""
1342        err = hdfRead(fname, dfName)
1343        //err not handled here
1344
1345        if(WaveExists($(dfName+":instrument:parameters:BSdiam"))) //canonical location
1346                Wave wBSdiameter = $(dfName+":instrument:parameters:BSdiam")
1347        elseif(WaveExists($(dfName+":instrument:parameters:BSXmm")))
1348                Wave wBSdiameter = $(dfName+":instrument:parameters:BSXmm")
1349        else
1350                print "Can't find Beamstop Diameter in " + fname
1351        endif
1352        value = wBSdiameter[0]
1353        KillWaves wBSdiameter
1354       
1355        return(value)   
1356end
1357
1358//beam stop diameter (millimeters)
1359Function WriteBeamStopDiamToHeader(fname,bs)
1360        String fname
1361        Variable bs     
1362       
1363        // your writer here
1364        Wave wBS
1365        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1366        Make /N=(1,1) wBS
1367        String groupName = "/instrument/parameters"
1368        String varName = "BSdiam"
1369        wBS[0] = bs //
1370        // your code returning value
1371        variable err
1372        err = hdfWrite(fname, groupName, varName, wBS)
1373        KillWaves wBS
1374        //err not handled here
1375        return(0)
1376End
1377
1378//source aperture diameter (millimeters)
1379Function getSourceApertureDiam(fname)
1380        String fname
1381        Variable value
1382        // your code returning value
1383        variable err
1384        string dfName = ""
1385        err = hdfRead(fname, dfName)
1386        //err not handled here
1387
1388        if(WaveExists($(dfName+":instrument:parameters:EApX")))
1389                Wave wSourceApertureDiam = $(dfName+":instrument:parameters:EApX") // canonical location
1390        elseif(WaveExists($(dfName+":instrument:parameters:derived_parameters:EApXmm")))
1391                Wave wSourceApertureDiam = $(dfName+":instrument:parameters:derived_parameters:EApXmm")
1392        else
1393                print "Can't find Source Aperture Diameter in " + fname
1394        endif   
1395        value = wSourceApertureDiam[0]
1396        KillWaves wSourceApertureDiam
1397       
1398        return(value)
1399end
1400
1401//Source Aperture diameter (millimeters)
1402Function WriteSourceApDiamToHeader(fname,source)
1403        String fname
1404        Variable source
1405       
1406        // your writer here
1407        Wave wsource
1408        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1409        Make /N=(1,1) wsource
1410        String groupName = "/instrument/parameters"
1411        String varName = "EApX"
1412        wsource[0] = source //
1413        // your code returning value
1414        variable err
1415        err = hdfWrite(fname, groupName, varName, wsource)
1416        KillWaves wsource
1417        //err not handled here
1418        return(0)
1419End
1420
1421//sample aperture diameter (millimeters)
1422Function getSampleApertureDiam(fname)
1423        String fname
1424        Variable value
1425        // your code returning value
1426        variable err
1427        string dfName = ""
1428        err = hdfRead(fname, dfName)
1429        //err not handled here
1430
1431        if(WaveExists($(dfName+":sample:diameter"))) //canonical location
1432                Wave wSampleApertureDiam = $(dfName+":sample:diameter")
1433        elseif(WaveExists($(dfName+":instrument:parameters:autoSampleAp:diameter"))) //canonical location
1434                Wave wSampleApertureDiam = $(dfName+":instrument:parameters:autoSampleAp:diameter")
1435        elseif (WaveExists($(dfName+":instrument:sample_aperture:geometry:shape:SApXmm")))
1436                Wave wSampleApertureDiam = $(dfName+":instrument:sample_aperture:geometry:shape:SApXmm")
1437        else
1438                print "Can't find Sample Aperture Diameter in " + fname
1439        endif   
1440        value = wSampleApertureDiam[0]
1441        KillWaves wSampleApertureDiam
1442
1443        return(value)   
1444end
1445
1446//Sample Aperture diameter (millimeters)
1447Function WriteSampleApDiamToHeader(fname,source)
1448        String fname
1449        Variable source
1450       
1451        // your writer here
1452        Wave wsource
1453        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1454        Make /N=(1,1) wsource
1455        String groupName = "/sample"
1456        String varName = "diameter"
1457        wsource[0] = source //
1458        // your code returning value
1459        variable err
1460        err = hdfWrite(fname, groupName, varName, wsource)
1461        KillWaves wsource
1462        //err not handled here
1463        return(0)
1464End
1465
1466//source AP to Sample AP distance (meters)
1467Function getSourceToSampleDist(fname)
1468        String fname
1469       
1470        Variable value
1471       
1472        // your code returning value
1473        variable err
1474        string dfName = ""
1475               
1476        err = hdfRead(fname, dfName)
1477        //err not handled here
1478       
1479        if(WaveExists($(dfName+":instrument:parameters:L1"))) //canonical location
1480                Wave wSourceToSampleDist = $(dfName+":instrument:parameters:L1")       
1481        elseif(WaveExists($(dfName+":instrument:parameters:L1mm")))
1482                Wave wSourceToSampleDist = $(dfName+":instrument:parameters:L1mm")
1483        elseif(WaveExists($(dfName+":instrument:parameters:derived_parameters:L1mm")))
1484                Wave wSourceToSampleDist = $(dfName+":instrument:parameters:derived_parameters:L1mm")
1485        else
1486                print "Can't find L1 in " + fname
1487        endif
1488       
1489        value = wSourceToSampleDist[0]/1000
1490        KillWaves wSourceToSampleDist   
1491               
1492        return(value)
1493end
1494
1495//Source aperture to sample aperture distance (meters)
1496Function WriteSrcToSamDistToHeader(fname,SSD)
1497        String fname
1498        Variable SSD
1499       
1500// your writer here
1501        Wave wSSD
1502        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1503        Make /N=(1,1) wSSD
1504        String groupName = "/instrument/parameters"
1505        String varName = "L1"
1506        wSSD[0] = SSD * 1000 //input in metres, converted to mm for saving to file.
1507        // your code returning value
1508        variable err
1509        err = hdfWrite(fname, groupName, varName, wSSD)
1510        KillWaves wSSD
1511        //err not handled here
1512       
1513        return(0)
1514End
1515
1516//wavelength (Angstroms)
1517Function getWavelength(fname)
1518        String fname
1519        Variable value
1520        // your code returning value
1521        variable err
1522        string dfName = ""
1523        err = hdfRead(fname, dfName)
1524        //err not handled here
1525
1526        //      Wave wWavelength = $(dfName+":data:LambdaA")
1527        //change 22/12/09 nha
1528        // all these locations to be deprecated
1529        if(WaveExists($(dfName+":instrument:velocity_selector:Lambda")))  // canonical location
1530                Wave wWavelength = $(dfName+":instrument:velocity_selector:Lambda")
1531        elseif(WaveExists($(dfName+":data:Lambda")))
1532                Wave wWavelength = $(dfName+":data:Lambda")
1533        elseif(WaveExists($(dfName+":data:LambdaA")))
1534                Wave wWavelength = $(dfName+":data:LambdaA")
1535        elseif(WaveExists($(dfName+":instrument:velocity_selector:LambdaA")))
1536                Wave wWavelength = $(dfName+":instrument:velocity_selector:LambdaA")
1537        else
1538                print "Can't find Lambda in " + fname
1539        endif
1540        value = wWavelength[0] 
1541        KillWaves wWavelength
1542       
1543        return(value)
1544end
1545
1546//wavelength (Angstroms)
1547Function WriteWavelengthToHeader(fname,wavelength)
1548        String fname
1549        Variable wavelength
1550       
1551        // your writer here
1552        Wave wWavelength
1553        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1554        Make /N=(1,1) wWavelength
1555        String groupName = "/instrument/velocity_selector"
1556        String varName = "Lambda"
1557        wWavelength[0] = wavelength //
1558        // your code returning value
1559        variable err
1560        err = hdfWrite(fname, groupName, varName, wWavelength)
1561       
1562        //and because Bill Hamilton is not happy with the NeXus naming, we write it to 3 other places
1563        //groupName = "/instrument/parameters"
1564        //err = hdfWrite(fname, groupName, varName, wWavelength)
1565        //velocity_selector group causes Igor crash in some cases
1566        //groupName = "/instrument/velocity_selector"
1567        //err = hdfWrite(fname, groupName, varName, wWavelength)
1568        //
1569        //groupName = "/data"
1570        //varName = "lambda"
1571        //err = hdfWrite(fname, groupName, varName, wWavelength)
1572       
1573        KillWaves wWavelength
1574        //err not handled here
1575
1576        return(0)
1577End
1578
1579
1580
1581//wavelength spread (FWHM)
1582Function getWavelengthSpread(fname)
1583        String fname
1584       
1585        Variable value
1586       
1587        // your code returning value
1588        variable err
1589        string dfName = ""
1590        err = hdfRead(fname, dfName)
1591        //err not handled here
1592       
1593        //velocity_selector group causes Igor crash
1594        if(WaveExists($(dfName+":instrument:velocity_selector:LambdaResFWHM_percent")))  //canonical location
1595                Wave wWavelengthSpread = $(dfName+":instrument:velocity_selector:LambdaResFWHM_percent")
1596        elseif(WaveExists($(dfName+":instrument:parameters:LambdaResFWHM_percent")))
1597                Wave wWavelengthSpread = $(dfName+":instrument:parameters:LambdaResFWHM_percent")
1598        else
1599                print "Can't find Wavelength Spread in " + fname
1600        endif
1601        value = wWavelengthSpread[0]   
1602        KillWaves wWavelengthSpread
1603       
1604        return(value)
1605end
1606
1607//wavelength spread (FWHM)
1608Function WriteWavelengthDistrToHeader(fname,wavelengthSpread)
1609        String fname
1610        Variable wavelengthSpread
1611       
1612        // your writer here
1613        Wave wWavelengthSpread
1614        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1615        Make /N=(1,1) wWavelengthSpread
1616        //velocity_selector group causes Igor crash
1617        String groupName = "/instrument/velocity_selector"
1618        String varName = "LambdaResFWHM_percent"
1619
1620        wWavelengthSpread[0] = wavelengthSpread
1621        // your code returning value
1622        variable err
1623        err = hdfWrite(fname, groupName, varName, wWavelengthSpread)
1624        KillWaves wWavelengthSpread
1625        //err not handled here
1626               
1627        return(0)
1628End
1629
1630// physical x-dimension of a detector pixel, in mm
1631Function getDetectorPixelXSize(fname)
1632        String fname
1633        Variable value
1634       
1635        // your code here returning value
1636        variable err
1637        string dfName = ""
1638        err = hdfRead(fname, dfName)
1639        //err not handled here
1640
1641        Wave wActiveArea = $(dfName+":instrument:detector:active_height")
1642        Wave w_x_bin = $(dfName+":instrument:detector:x_bin")
1643        Variable numPixels = dimsize(w_x_bin, 0)
1644        value = wActiveArea[0]/numPixels
1645        KillWaves wActiveArea
1646        KillWaves w_x_bin
1647       
1648        return(value)
1649end
1650
1651// physical y-dimension of a detector pixel, in mm
1652Function getDetectorPixelYSize(fname)
1653        String fname
1654        Variable value
1655       
1656        // your code here returning value
1657        variable err
1658        string dfName = ""
1659        err = hdfRead(fname, dfName)
1660        //err not handled here
1661
1662        Wave wActiveArea = $(dfName+":instrument:detector:active_width")
1663        Wave w_y_bin = $(dfName+":instrument:detector:y_bin")
1664        Variable numPixels = dimsize(w_y_bin, 0)
1665        value = wActiveArea[0]/numPixels
1666        KillWaves wActiveArea
1667        KillWaves w_y_bin
1668                       
1669        return(value)
1670end
1671
1672//transmission detector count (unused, return 0)
1673//
1674Function getTransDetectorCounts(fname)
1675        String fname
1676       
1677        Variable value
1678       
1679        // your code returning value
1680       
1681        return(0)
1682end
1683
1684
1685//total count time (seconds)
1686Function getCountTime(fname)
1687        String fname
1688        Variable value
1689        // your code returning value
1690        variable err
1691        string dfName = ""
1692        err = hdfRead(fname, dfName)
1693        //err not handled here
1694
1695        Wave wTime1 = $(dfName+":monitor:bm1_time")
1696        value = wTime1[0]       
1697        KillWaves wTime1
1698       
1699        return(value)
1700end
1701
1702
1703Function getPhysicalWidth(fname)
1704        String fname
1705        Variable value
1706        // your code returning value
1707        variable err
1708        string dfName = ""
1709        err = hdfRead(fname, dfName)
1710        //err not handled here
1711
1712        Wave wPhysicalWidth = $(dfName+":instrument:detector:active_width")
1713        value = wPhysicalWidth[0]/10
1714        KillWaves wPhysicalWidth
1715               
1716        return(value)
1717end
1718
1719Function/S getSICSVersion(fname)
1720        String fname
1721        String value
1722        // your code returning value
1723        variable err
1724        string dfName = ""
1725        err = hdfRead(fname, dfName)
1726        //err not handled here
1727
1728        Wave/T wSICSVersion = $(dfName+":sics_release")
1729        value = wSICSVersion[0]
1730        KillWaves wSICSVersion
1731       
1732        return(value)
1733end
1734
1735Function/S getHDFversion(fname)
1736        String fname
1737        String value
1738        // your code returning value
1739        variable err
1740        string dfName = ""
1741        string attribute = "HDF5_Version"
1742        err = hdfReadAttribute(fname, dfName, "/", 1, attribute)
1743//      string attribute ="signal"
1744//      err = hdfReadAttribute(fname, dfName, "/entry/data/hmm_xy", 2, attribute)
1745        //err not handled here
1746
1747        Wave/T wHDF5_Version = $(dfName+":"+attribute)
1748        value = wHDF5_Version[0]       
1749//      KillWaves wHDF5_Version
1750
1751        return(value)
1752end
1753
1754
1755//reads the wavelength from a reduced data file (not very reliable)
1756// - does not work with NSORTed files
1757// - only used in FIT/RPA (which itself is almost NEVER used...)
1758//
1759// DOES NOT NEED TO BE CHANGED IF USING NCNR DATA WRITER
1760Function GetLambdaFromReducedData(tempName)
1761        String tempName
1762       
1763        String junkString
1764        Variable lambdaFromFile, fileVar
1765        lambdaFromFile = 6.0
1766        Open/R/P=catPathName fileVar as tempName
1767        FReadLine fileVar, junkString
1768        FReadLine fileVar, junkString
1769        FReadLine fileVar, junkString
1770        if(strsearch(LowerStr(junkString),"lambda",0) != -1)
1771                FReadLine/N=11 fileVar, junkString
1772                FReadLine/N=10 fileVar, junkString
1773                lambdaFromFile = str2num(junkString)
1774        endif
1775        Close fileVar
1776       
1777        return(lambdaFromFile)
1778End
1779
1780/////   TRANSMISSION RELATED FUNCTIONS    ////////
1781//box coordinate are returned by reference
1782//
1783// box to sum over is bounded (x1,y1) to (x2,y2)
1784//
1785// this function returns the bounding coordinates as stored in
1786// the header
1787//
1788// if not using the NCNR Transmission module, this function default to
1789// returning 0000, and no changes needed
1790Function getXYBoxFromFile(fname,x1,x2,y1,y2)
1791        String fname
1792        Variable &x1,&x2,&y1,&y2
1793       
1794        // return your bounding box coordinates or default values of 0
1795        x1=0
1796        y1=0
1797        x2=0
1798        y2=0
1799
1800        // your code returning value
1801        variable err
1802        string dfName = ""
1803        err = hdfRead(fname, dfName)
1804        //err not handled here
1805
1806       
1807        Wave wX1 = $(dfName+":reduce:x1")
1808        if (waveexists(wX1) == 0)
1809                //Waves don't exists which means an XY box has not been set for this file.
1810                //Hence return 0 bounding boxes (default)
1811        else
1812                x1 = wX1[0]
1813                Wave wX2 = $(dfName+":reduce:x2")
1814                x2 = wX2[0]
1815                Wave wY1 = $(dfName+":reduce:y1")
1816                y1 = wY1[0]
1817                Wave wY2 = $(dfName+":reduce:y2")
1818                y2 = wY2[0]
1819        endif
1820       
1821        KillWaves/Z wX1, wX2, wY1, wY2
1822        return(0)
1823
1824End
1825
1826// Writes the coordinates of the box to the header after user input
1827//
1828// box to sum over bounded (x1,y1) to (x2,y2)
1829//
1830// if not using the NCNR Transmission module, this function is null
1831Function WriteXYBoxToHeader(fname,x1,x2,y1,y2)
1832        String fname
1833        Variable x1,x2,y1,y2
1834
1835        String groupName = "/reduce"
1836        variable err
1837               
1838        Wave wX1
1839        Make/O/N=(1,1) wX1
1840        Wave wX2
1841        Make/O/N=(1,1) wX2
1842        Wave wY1
1843        Make/O/N=(1,1) wY1
1844        Wave wY2
1845        Make/O/N=(1,1) wY2
1846               
1847        wX1[0] = x1
1848        wX2[0] = x2
1849        wY1[0] = y1
1850        wY2[0] = y2     
1851       
1852        String varName = "x1"   
1853        err = hdfWrite(fname, groupName, varName,wX1)
1854        varName = "x2"
1855        err = hdfWrite(fname, groupName, varName,wX2)
1856        varName = "y1"
1857        err = hdfWrite(fname, groupName, varName,wY1)
1858        varName = "y2"
1859        err = hdfWrite(fname, groupName, varName,wY2)   
1860       
1861        KillWaves wX1,wX2,wY1,wY2
1862       
1863        //err not handled here
1864        return(0)       
1865
1866End
1867
1868// for transmission calculation, writes an NCNR-specific alphanumeric identifier
1869// (suffix of the data file)
1870//
1871//AJJ June 3 2010 - Note!! For ANSTO data the "suffix" is just the filename.
1872Function WriteAssocFileSuffixToHeader(fname,assoc_fname)
1873        String fname,assoc_fname
1874               
1875// your writer here
1876        Wave/T wAssoc_fname
1877        //nha ??? Should make this wave in our own DataFolder to avoid clashing names.
1878        Make/T /N=(1,1) wAssoc_fname
1879       
1880        String varName =""
1881        String groupName = "/reduce"
1882        if(isTransFile(fname))
1883                varName = "empty_beam_file_name"
1884        elseif(isScatFile(fname))
1885                varName = "transmission_file_name"
1886        endif
1887
1888        wAssoc_fname[0] = assoc_fname
1889        // your code returning value
1890        variable err
1891        err = hdfWrite(fname, groupName, varName, wAssoc_fname)
1892        KillWaves wAssoc_fname
1893        //err not handled here
1894        return(0)
1895end
1896
1897Function/S GetFileAssoc(fname)
1898        String fname
1899       
1900        String assoc_fname
1901        String groupName = ":reduce:"
1902       
1903        String varName = ""
1904        if(isTransFile(fname))
1905                varName = "empty_beam_file_name"
1906        elseif(isScatFile(fname))
1907                varName = "transmission_file_name"
1908        endif
1909       
1910        variable err
1911        string dfName = ""
1912        err = hdfRead(fname, dfName)
1913        //err not handled here
1914
1915        Wave/T wAssoc_fname = $(dfName+groupName+varName)
1916        if (waveexists(wAssoc_fname) == 1)
1917                assoc_fname =wAssoc_fname[0]   
1918        else
1919                assoc_fname = ""
1920        endif
1921        KillWaves/Z wAssoc_fname
1922       
1923        return(assoc_fname)
1924end
1925
1926Function hdfReadAttribute(fname, dfName, nodeName, nodeType, attributeStr)
1927// this is a copy of hdfRead, and could be incorporated back into hdfRead.
1928       
1929        String fname, &dfName, nodeName, attributeStr
1930        variable nodeType
1931        String nxentryName
1932        variable err=0,fileID   
1933        String cDF = getDataFolder(1), temp
1934        String fname_temp = ParseFilePath(0, fname, ":", 1, 0)
1935
1936       
1937        String fileSuffix
1938       
1939        if(strsearch(fname_temp,".nx.hdf",0,2)>=0)
1940                fileSuffix=".nx.hdf"
1941        else
1942                err = 1
1943                abort "unrecognised file suffix. Not .nx.hdf"
1944        endif
1945       
1946        dfName = "root:packages:quokka:"+removeending(fname_temp,fileSuffix)
1947       
1948        //folder must already exist i.e. hdfRead must have already been called
1949        if(!DataFolderExists(dfName))
1950                // possibly call an hdfRead from here
1951                return err
1952        endif
1953       
1954        //test for the name of nxentry
1955        if(DataFolderExists(dfName+":"+removeending(fname_temp,fileSuffix)))
1956                nxentryName = removeending(fname_temp,fileSuffix)
1957        elseif(DataFolderExists(dfName+":"+"entry1"))
1958                nxentryName = "entry1"
1959        else
1960                print "NXentry not found"
1961                return err
1962        endif
1963       
1964        //this is the stupid bit.
1965        // If you're looking for attributes of the root node, then nodename = "/"
1966        // If you're looking for attributes     of the nxentry node, then e.g. nodename ="/entry/instrument"
1967        // /entry is replaced with nxentryName
1968        nodeName = ReplaceString("entry", nodeName, nxentryName)       
1969       
1970        //convert nodeName to data folder string
1971        String dfNodeName = nodeName
1972        dfNodeName = ReplaceString( "/", nodeName, ":")
1973        dfName = dfName + dfNodeName
1974        if(nodeType == 2) //data item (dataset)
1975                //remove the end of dfName so that it points to a folder and not a dataset
1976                variable length = strlen(dfName)
1977                variable position = strsearch(dfName, ":", length, 1) // search backwards to find last :
1978                // to do - truncate string to remove dataset
1979                string truncate = "\"%0." + num2str(position) + "s\""
1980                sprintf dfName, truncate, dfName
1981        endif
1982       
1983        setDataFolder dfName
1984       
1985        try     
1986                HDF5OpenFile /R /Z fileID  as fname
1987                if(!fileID)
1988                        err = 1
1989                        abort "couldn't load HDF5 file"
1990                endif
1991
1992                HDF5LoadData /O /Q /Z /TYPE=(nodeType) /A=attributeStr, fileID, nodeName
1993
1994                if (V_flag!=0)
1995                        print "couldn't load attribute " + attributeStr
1996                endif
1997        catch
1998
1999        endtry
2000        if(fileID)
2001                HDF5CloseFile /Z fileID
2002        endif
2003
2004// add the name of the root node to dfName
2005// in the case of sensitivity files aka DIV files, don't append a root node to dfName
2006        if(DataFolderExists(dfName+":"+removeending(fname_temp,fileSuffix)))
2007                dfName = dfName+":"+removeending(fname_temp,fileSuffix)  //for files written Dec 2009 and after
2008        elseif(DataFolderExists(dfName+":"+"entry1"))
2009                dfName = dfName+":entry1" //for files written before Dec 2009
2010        endif
2011
2012        setDataFolder $cDF
2013        return err
2014end
2015
2016Function hdfRead(fname, dfName)
2017        //Reads hdf5 file into root:packages:quokka:fname
2018        //N. Hauser. 16/12/08
2019        // 29/1/09 - hdf file must have .nx.hdf or .div suffix
2020       
2021        String fname, &dfName
2022        variable err=0,fileID
2023        String cDF = getDataFolder(1), temp
2024        String fname_temp = ParseFilePath(0, fname, ":", 1, 0)
2025               
2026        String fileSuffix
2027        if(strsearch(fname_temp,".nx.hdf",0,2)>=0)
2028                fileSuffix=".nx.hdf"
2029        //elseif(strsearch(fname_temp,".div",0,2)>=0)
2030        //      fileSuffix=".div"
2031        else
2032                err = 1
2033                abort "unrecognised file suffix. Not .nx.hdf"
2034        endif
2035       
2036        dfName = "root:packages:quokka:"+removeending(fname_temp,fileSuffix)
2037       
2038        //if(!DataFolderExists(dfName))
2039        //      return err
2040        //else         
2041                newDataFolder /O root:packages
2042                newDataFolder /O /S root:packages:quokka
2043                newDataFolder /O /S $dfName
2044        //endif
2045       
2046        try     
2047                HDF5OpenFile /R /Z fileID  as fname
2048                if(!fileID)
2049                        err = 1
2050                        abort "couldn't load HDF5 file"
2051                endif
2052                HDF5LoadGroup /O /R /Z  :, fileID, "."
2053        catch
2054
2055        endtry
2056        if(fileID)
2057                HDF5CloseFile /Z fileID
2058        endif
2059
2060        // add the name of the root node to dfName
2061        // in the case of sensitivity files aka DIV files, don't append a root node to dfName
2062        if(DataFolderExists(dfName+":"+removeending(fname_temp,fileSuffix)))
2063                dfName = dfName+":"+removeending(fname_temp,fileSuffix)  //for files written Dec 2009 and after
2064        elseif(DataFolderExists(dfName+":"+"entry1"))
2065                dfName = dfName+":entry1" //for files written before Dec 2009
2066        endif
2067
2068        setDataFolder $cDF
2069        return err
2070end
2071
2072Function hdfWrite(fname, groupName, varName, wav)
2073        //Write Wave 'wav' to hdf5 file 'fname'
2074        //N. Hauser. nha 8/1/09
2075               
2076        String fname, groupName, varName
2077        Wave wav
2078       
2079        variable err=0, fileID,groupID
2080        String cDF = getDataFolder(1), temp
2081        String fname_temp = ParseFilePath(0, fname, ":", 1, 0)
2082        String NXentry_name
2083                       
2084        try     
2085                HDF5OpenFile/Z fileID  as fname  //open file read-write
2086                if(!fileID)
2087                        err = 1
2088                        abort "HDF5 file does not exist"
2089                endif
2090               
2091                //get the NXentry node name
2092                HDF5ListGroup /TYPE=1 fileID, "/"
2093                //remove trailing ; from S_HDF5ListGroup
2094                NXentry_name = S_HDF5ListGroup
2095                NXentry_name = ReplaceString(";",NXentry_name,"")
2096                if(strsearch(NXentry_name,":",0)!=-1) //more than one entry under the root node
2097                        err = 1
2098                        abort "More than one entry under the root node. Ambiguous"
2099                endif
2100                //concatenate NXentry node name and groupName   
2101                groupName = "/" + NXentry_name + groupName
2102                HDF5OpenGroup /Z fileID , groupName, groupID
2103
2104        //      !! At the moment, there is no entry for sample thickness in our data file
2105        //      therefore create new HDF5 group to enable write / patch command
2106        //      comment out the following group creation once thickness appears in revised file
2107       
2108                if(!groupID)
2109                        HDF5CreateGroup /Z fileID, groupName, groupID
2110                        //err = 1
2111                        //abort "HDF5 group does not exist"
2112                else
2113                        // get attributes and save them
2114                        //HDF5ListAttributes /Z fileID, groupName    this is returning null. expect it to return semicolon delimited list of attributes
2115                        //Wave attributes = S_HDF5ListAttributes
2116                endif
2117       
2118                HDF5SaveData /O /Z /IGOR=0  wav, groupID, varName
2119                if (V_flag != 0)
2120                        err = 1
2121                        abort "Cannot save wave to HDF5 dataset" + varName
2122                endif   
2123               
2124                //HDF5SaveData /O /Z /IGOR=0 /A=attributes groupID, varName
2125                //if (V_flag != 0)
2126                ////    err = 1
2127                //      abort "Cannot save attributes to HDF5 dataset"
2128                //endif
2129        catch
2130
2131        endtry
2132       
2133        if(groupID)
2134                HDF5CloseGroup /Z groupID
2135        endif
2136       
2137        if(fileID)
2138                HDF5CloseFile /Z fileID
2139        endif
2140
2141        setDataFolder $cDF
2142        return err
2143end
2144
2145Function DoEdgeCorrection(type)
2146        String type
2147        variable err = 0
2148       
2149        WAVE typeData=$("root:Packages:NIST:"+type+":data")
2150       
2151        //nha workaround for high counts on edges
2152        //copy second column to first column
2153        typeData[][0] = typeData[p][1]
2154        //copy second last column to last column
2155        typeData[][191] = typeData[p][190]
2156        //copy second row to first row
2157        typeData[0][] = typeData[1][q]
2158        //copy second last row to last row
2159        typeData[191][] = typeData[190][q]
2160       
2161        return err     
2162end
2163
2164////// OCT 2009, facility specific bits from ProDiv()
2165//"type" is the data folder that has the corrected, patched, and normalized DIV data array
2166//
2167// the header of this file is rather unimportant. Filling in a title at least would be helpful/
2168//
2169Function Write_DIV_File(type)
2170        String type
2171       
2172        // Your file writing function here. Don't try to duplicate the VAX binary format...
2173        WriteHeaderAndWork(type)
2174       
2175        return(0)
2176End
2177
2178////// OCT 2009, facility specific bits from MonteCarlo functions()
2179//"type" is the data folder that has the data array that is to be (re)written as a full
2180// data file, as if it was a raw data file
2181//
2182// not really necessary
2183//
2184Function/S Write_RawData_File(type,fullpath,dialog)
2185        String type,fullpath
2186        Variable dialog         //=1 will present dialog for name
2187       
2188        // Your file writing function here. Don't try to duplicate the VAX binary format...
2189        Print "Write_RawData_File stub"
2190       
2191        return(fullpath)
2192End
Note: See TracBrowser for help on using the repository browser.