source: sans/Dev/branches/nxcansas_writer/NCNR_User_Procedures/Reduction/SANS/Write_SANS_NXcanSAS.ipf @ 1194

Last change on this file since 1194 was 1194, checked in by krzywon, 3 years ago

Output NSORT data as NXcanSAS. Mostly tested, but a couple TODOs and FIXMEs left.

File size: 17.8 KB
Line 
1#pragma rtGlobals=3             // Use modern global access method.
2#pragma version=5.0
3#pragma IgorVersion=6.1
4
5#include <HDF5 Browser>
6
7//************************
8// Vers 1.15 20171003
9//
10//************************
11
12
13///////////////////////////////////////////////////////////////////////////
14// - WriteNxCanSAS1D - Method for writing 1D NXcanSAS data
15// Creates an HDF5 file, with reduced 1D data and stores all meta data
16// If dialog and fullpath are left blank (0 and ""), fake data will be used
17
18Function WriteNxCanSAS1D(type,fullpath,dialog)
19        // Define input variables
20        String type // data location, in memory, relative to root:Packages:NIST:
21        String fullpath // file path and name where data will be saved
22        Variable dialog // if 1, prompt user for file path, otherwise, use fullpath
23       
24        // Define local function variables
25        Variable fileID
26        String destStr="", parentBase, nxcansasBase
27        String/G base = "root:NXcanSAS_file"
28       
29        // Define local waves
30//      Wave/T vals,attr,attrVals
31       
32        // Define folder for data heirarchy
33        NewDataFolder/O/S root:NXcanSAS_file
34       
35        // Check fullpath and dialog
36        fileID = NXcanSAS_OpenOrCreate(dialog,fullpath,base)
37       
38        Variable sasentry = NumVarOrDefault("root:Packages:NIST:gSASEntryNumber", 1)
39        sPrintf parentBase,"%s:sasentry%d",base,sasentry // Igor memory base path for all
40        sPrintf nxcansasBase,"/sasentry%d/",sasentry // HDF5 base path for all
41
42        destStr = "root:Packages:NIST:"+type
43        //*****these waves MUST EXIST, or IGOR Pro will crash, with a type 2 error****
44        WAVE intw = $(destStr + ":integersRead")
45        WAVE rw = $(destStr + ":realsRead")
46        WAVE/T textw=$(destStr + ":textRead")
47        WAVE qvals =$(destStr + ":qval")
48        WAVE inten=$(destStr + ":aveint")
49        WAVE sig=$(destStr + ":sigave")
50        WAVE qbar = $(destStr + ":QBar")
51        WAVE sigmaq = $(destStr + ":SigmaQ")
52        WAVE fsubs = $(destStr + ":fSubS")
53
54        ///////////////////////////////////////////////////////////////////////////
55        // Write all data
56       
57        // Define common attribute waves
58        Make/T/O/N=1 empty = {""}
59        Make/T/O/N=1 units = {"units"}
60        Make/T/O/N=1 inv_cm = {"1/cm"}
61        Make/T/O/N=1 inv_angstrom = {"1/A"}
62       
63        // Run Name and title
64        NewDataFolder/O/S $(parentBase)
65        Make/O/T/N=1 $(parentBase + ":title") = {textw[6]}
66        CreateStrNxCansas(fileID,nxcansasBase,"","title",$(parentBase + ":title"),empty,empty)
67        Make/O/T/N=1 $(parentBase + ":run") = {textw[0]}
68        CreateStrNxCansas(fileID,nxcansasBase,"","run",$(parentBase + ":run"),empty,empty)
69       
70        // SASData
71        String dataParent = nxcansasBase + "sasdata/"
72        // Create SASdata entry
73        String dataBase = parentBase + ":sasdata"
74        NewDataFolder/O/S $(dataBase)
75        Make/O/T/N=5 $(dataBase + ":attr") = {"canSAS_class","signal","I_axes","NX_class","Q_indices", "timestamp"}
76        Make/O/T/N=5 $(dataBase + ":attrVals") = {"SASdata","I","Q","NXdata","0",textw[1]}
77        CreateStrNxCansas(fileID,dataParent,"","",empty,$(dataBase + ":attr"),$(dataBase + ":attrVals"))
78        // Create q entry
79        NewDataFolder/O/S $(dataBase + ":q")
80        Make/O/T/N=2 $(dataBase + ":q:attr") = {"units","resolutions"}
81        Make/O/T/N=2 $(dataBase + ":q:attrVals") = {"1/angstrom","Qdev"}
82        CreateVarNxCansas(fileID,dataParent,"sasdata","Q",qvals,$(dataBase + ":q:attr"),$(dataBase + ":q:attrVals"))
83        // Create i entry
84        NewDataFolder/O/S $(dataBase + ":i")
85        Make/O/T/N=2 $(dataBase + ":i:attr") = {"units","uncertainties"}
86        Make/O/T/N=2 $(dataBase + ":i:attrVals") = {"1/cm","Idev"}
87        CreateVarNxCansas(fileID,dataParent,"sasdata","I",inten,$(dataBase + ":i:attr"),$(dataBase + ":i:attrVals"))
88        // Create idev entry
89        CreateVarNxCansas(fileID,dataParent,"sasdata","Idev",sig,units,inv_cm)
90        // Create qdev entry
91        CreateVarNxCansas(fileID,dataParent,"sasdata","Qdev",sigmaq,units,inv_angstrom)
92        CreateVarNxCansas(fileID,dataParent,"sasdata","Qmean",qbar,units,inv_angstrom)
93       
94        // Write all meta data
95        if (CmpStr(type,"NSORT") == 0)
96                Wave process = $(destStr + ":process")
97                // TODO: Write to file
98        Else
99                WriteMetaData(fileID,parentBase,nxcansasBase,rw,textw)
100        EndIf
101       
102        //
103        ///////////////////////////////////////////////////////////////////////////
104       
105        // Close the file
106        if(fileID)
107                HDF5CloseFile /Z fileID
108        endif
109       
110        KillDataFolder/Z $base
111       
112End
113
114//
115///////////////////////////////////////////////////////////////////////////
116
117
118///////////////////////////////////////////////////////////////////////////
119// - WriteNxCanSAS2D - Method for writing 2D NXcanSAS data
120// Creates an HDF5 file, generates reduced 2D data and stores all meta data
121// If dialog and fullpath are left blank (0 and ""), fake data will be used
122
123Function WriteNxCanSAS2D(type,fullpath,dialog)
124        // Define input variables
125        String type // data location, in memory, relative to root:Packages:NIST:
126        String fullpath // file path and name where data will be saved
127        Variable dialog // if 1, prompt user for file path, otherwise, use fullpath
128       
129        // Define local function variables
130        Variable fileID
131        String destStr="",typeStr="", parentBase, nxcansasBase
132        String/G base = "root:NXcanSAS_file"
133       
134        // Define local waves
135//      Wave/T vals,attr,attrVals
136       
137        // Define folder for data heirarchy
138        NewDataFolder/O/S root:NXcanSAS_file
139       
140        // Check fullpath and dialog
141        fileID = NXcanSAS_OpenOrCreate(dialog,fullpath,base)
142
143        Variable sasentry = NumVarOrDefault("root:Packages:NIST:gSASEntryNumber", 1)
144        sPrintf parentBase,"%s:sasentry%d",base,sasentry // Igor memory base path for all
145        sPrintf nxcansasBase,"/sasentry%d/",sasentry // HDF5 base path for all
146       
147        destStr = "root:Packages:NIST:"+type
148
149        //must select the linear_data to export
150        NVAR isLog = $(destStr+":gIsLogScale")
151        if(isLog==1)
152                typeStr = ":linear_data"
153        else
154                typeStr = ":data"
155        endif
156        NVAR pixelsX = root:myGlobals:gNPixelsX
157        NVAR pixelsY = root:myGlobals:gNPixelsY
158        Wave data=$(destStr+typeStr)
159        Wave data_err=$(destStr+":linear_data_error")
160        WAVE intw=$(destStr + ":integersRead")
161        WAVE rw=$(destStr + ":realsRead")
162        WAVE/T textw=$(destStr + ":textRead")
163       
164        ///////////////////////////////////////////////////////////////////////////
165        // Compute Qx, Qy data from pixel space
166       
167        Duplicate/O data,qx_val,qy_val,z_val,qval,qz_val,phi,r_dist
168       
169        Variable xctr,yctr,sdd,lambda,pixSize
170        xctr = rw[16]
171        yctr = rw[17]
172        sdd = rw[18]
173        lambda = rw[26]
174        pixSize = rw[13]/10             //convert mm to cm (x and y are the same size pixels)
175       
176        qx_val = CalcQx(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)          //+1 converts to detector coordinate system
177        qy_val = CalcQy(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)
178       
179        Redimension/N=(pixelsX*pixelsY) qx_val,qy_val,z_val
180
181        Variable L2 = rw[18]
182        Variable BS = rw[21]
183        Variable S1 = rw[23]
184        Variable S2 = rw[24]
185        Variable L1 = rw[25]
186        Variable lambdaWidth = rw[27]   
187        Variable usingLenses = rw[28]           //new 2007
188
189        Variable vz_1 = 3.956e5         //velocity [cm/s] of 1 A neutron
190        Variable g = 981.0                              //gravity acceleration [cm/s^2]
191        Variable m_h    = 252.8                 // m/h [=] s/cm^2
192
193        Variable acc,ssd,lambda0,yg_d,qstar
194               
195        G = 981.  //!   ACCELERATION OF GRAVITY, CM/SEC^2
196        acc = vz_1              //      3.956E5 //!     CONVERT WAVELENGTH TO VELOCITY CM/SEC
197        SDD = L2        *100    //1317
198        SSD = L1        *100    //1627          //cm
199        lambda0 = lambda                //              15
200        YG_d = -0.5*G*SDD*(SSD+SDD)*(LAMBDA0/acc)^2
201        qstar = -2*pi/lambda0*2*yg_d/sdd
202
203        // the gravity center is not the resolution center
204        // gravity center = beam center
205        // resolution center = offset y = dy + (2)*yg_d
206        ///************
207        // do everything to write out the resolution too
208        // un-comment these if you want to write out qz_val and qval too, then use the proper save command
209        qval = CalcQval(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)
210        qz_val = CalcQz(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)
211        //      phi = FindPhi( pixSize*((p+1)-xctr) , pixSize*((q+1)-yctr))             //(dx,dy)
212        //      r_dist = sqrt(  (pixSize*((p+1)-xctr))^2 +  (pixSize*((q+1)-yctr))^2 )          //radial distance from ctr to pt
213        phi = FindPhi( pixSize*((p+1)-xctr) , pixSize*((q+1)-yctr)+(2)*yg_d)            //(dx,dy+yg_d)
214        r_dist = sqrt(  (pixSize*((p+1)-xctr))^2 +  (pixSize*((q+1)-yctr)+(2)*yg_d)^2 )         //radial distance from ctr to pt
215        Redimension/N=(pixelsX*pixelsY) qz_val,qval,phi,r_dist
216        Make/O/N=(2,pixelsX,pixelsY) qxy_vals
217        //everything in 1D now
218        Duplicate/O qval SigmaQX,SigmaQY
219        Make/O/N=(pixelsX,pixelsY) shadow
220        Make/O/N=(2,pixelsX,pixelsY) SigmaQ_combined
221
222        //Two parameters DDET and APOFF are instrument dependent.  Determine
223        //these from the instrument name in the header.
224        //From conversation with JB on 01.06.99 these are the current good values
225        Variable DDet
226        NVAR apOff = root:myGlobals:apOff               //in cm
227        DDet = rw[10]/10                        // header value (X) is in mm, want cm here
228
229        Variable ret1,ret2,ret3,jj
230        Variable nq = 0
231        Variable ii = 0
232       
233        do
234                jj = 0
235                do
236                        nq = ii * pixelsX + jj
237                        get2DResolution(qval[nq],phi[nq],lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,pixSize,usingLenses,r_dist[nq],ret1,ret2,ret3)
238                        qxy_vals[0][ii][jj] = qx_val[nq]
239                        qxy_vals[1][ii][jj] = qy_val[nq]
240                        SigmaQ_combined[0][ii][jj] = ret1       
241                        SigmaQ_combined[1][ii][jj] = ret2
242                        shadow[ii][jj] = ret3   
243                        jj+=1
244                while(jj<pixelsX)
245                ii+=1
246        while(ii<pixelsY)
247        //
248        ///////////////////////////////////////////////////////////////////////////
249
250       
251        ///////////////////////////////////////////////////////////////////////////
252        // Write all data
253       
254        // Define common attribute waves
255        Make/O/T/N=1 empty = {""}
256        Make/O/T/N=1 units = {"units"}
257        Make/O/T/N=1 inv_cm = {"1/cm"}
258        Make/O/T/N=1 inv_angstrom = {"1/A"}
259       
260        // Run Name and title
261        NewDataFolder/O/S $(parentBase)
262        Make/O/T/N=1 $(parentBase + ":title") = {textw[6]}
263        CreateStrNxCansas(fileID,nxcansasBase,"","title",$(parentBase + ":title"),empty,empty)
264        Make/O/T/N=1 $(parentBase + ":run") = {textw[0]}
265        CreateStrNxCansas(fileID,nxcansasBase,"","run",$(parentBase + ":run"),empty,empty)
266       
267        // SASData
268        String dataParent = nxcansasBase + "sasdata/"
269        // Create SASdata entry
270        String dataBase = parentBase + ":sasdata"
271        NewDataFolder/O/S $(dataBase)
272        Make/O/T/N=5 $(dataBase + ":attr") = {"canSAS_class","signal","I_axes","NX_class","Q_indices", "timestamp"}
273        Make/O/T/N=5 $(dataBase + ":attrVals") = {"SASdata","I","Q,Q","NXdata","0,1",textw[1]}
274        CreateStrNxCansas(fileID,dataParent,"","",empty,$(dataBase + ":attr"),$(dataBase + ":attrVals"))
275        // Create i entry
276        NewDataFolder/O/S $(dataBase + ":i")
277        Make/O/T/N=2 $(dataBase + ":i:attr") = {"units","uncertainties"}
278        Make/O/T/N=2 $(dataBase + ":i:attrVals") = {"1/cm","Idev"}
279        CreateVarNxCansas(fileID,dataParent,"sasdata","I",data,$(dataBase + ":i:attr"),$(dataBase + ":i:attrVals"))
280
281        //
282        // TODO: Reinstate Qdev/resolutions when I can fix the reader issue
283        //
284
285        // Create qx and qy entry
286        NewDataFolder/O/S $(dataBase + ":q")
287        Make/O/T/N=2 $(dataBase + ":q:attr") = {"units"}//,"resolutions"}
288        Make/O/T/N=2 $(dataBase + ":q:attrVals") = {"1/angstrom"}//,"Qdev"}
289        CreateVarNxCansas(fileID,dataParent,"sasdata","Q",qxy_vals,$(dataBase + ":q:attr"),$(dataBase + ":q:attrVals"))
290       
291        // Create idev entry
292        CreateVarNxCansas(fileID,dataParent,"sasdata","Idev",data_err,units,inv_cm)
293        // Create qdev entry
294        CreateVarNxCansas(fileID,dataParent,"sasdata","Qdev",SigmaQ_combined,units,inv_angstrom)
295        // Create shadwfactor entry
296        CreateVarNxCansas(fileID,dataParent,"sasdata","ShadowFactor",shadow,empty,empty)
297       
298        // Write all meta data
299        WriteMetaData(fileID,parentBase,nxcansasBase,rw,textw)
300       
301        // Close the file
302        if(fileID)
303                HDF5CloseFile /Z fileID
304        endif
305       
306        KillDataFolder/Z $base
307       
308End
309
310//
311///////////////////////////////////////////////////////////////////////////
312
313///////////////////////////////////////////////////////////////////////////
314// - WriteMetaData - Method used to write non data elements into NXcanSAS
315// format. This is common between 1D and 2D data sets.
316
317Function WriteMetaData(fileID,base,parentBase,rw,textw)
318        String base,parentBase
319        Variable fileID
320        Wave rw
321        Wave/T textw
322       
323        // Define common attribute waves
324        Make/T/O/N=1 empty = {""}
325        Make/T/O/N=1 units = {"units"}
326        Make/T/O/N=1 m = {"m"}
327        Make/T/O/N=1 mm = {"mm"}
328        Make/T/O/N=1 cm = {"cm"}
329        Make/T/O/N=1 pixel = {"pixel"}
330        Make/T/O/N=1 angstrom = {"A"}
331       
332        // SASinstrument
333        String instrParent = parentBase + "sasinstrument/"
334        // Create SASinstrument entry
335        String instrumentBase = base + ":sasinstrument"
336        NewDataFolder/O/S $(instrumentBase)
337        Make/O/T/N=5 $(instrumentBase + ":attr") = {"canSAS_class","NX_class"}
338        Make/O/T/N=5 $(instrumentBase + ":attrVals") = {"SASinstrument","NXinstrument"}
339        CreateStrNxCansas(fileID,instrParent,"","",empty,$(instrumentBase + ":attr"),$(instrumentBase + ":attrVals"))
340       
341        // SASaperture
342        String apertureParent = instrParent + "sasaperture/"
343        // Create SASaperture entry
344        String apertureBase = instrumentBase + ":sasaperture"
345        NewDataFolder/O/S $(apertureBase)
346        Make/O/T/N=5 $(apertureBase + ":attr") = {"canSAS_class","NX_class"}
347        Make/O/T/N=5 $(apertureBase + ":attrVals") = {"SASaperture","NXaperture"}
348        CreateStrNxCansas(fileID,apertureParent,"","",empty,$(apertureBase + ":attr"),$(apertureBase + ":attrVals"))
349        // Create SASaperture shape entry
350        Make/O/T/N=1 $(apertureBase + ":shape") = {"pinhole"}
351        CreateStrNxCansas(fileID,apertureParent,"sasaperture","shape",$(apertureBase + ":shape"),empty,empty)
352        // Create SASaperture x_gap entry
353        Make/O/N=1 $(apertureBase + ":x_gap") = {rw[24]}
354        CreateVarNxCansas(fileID,apertureParent,"sasaperture","x_gap",$(apertureBase + ":x_gap"),units,mm)
355        // Create SASaperture y_gap entry
356        Make/O/N=1 $(apertureBase + ":y_gap") = {rw[24]}
357        CreateVarNxCansas(fileID,apertureParent,"sasaperture","y_gap",$(apertureBase + ":y_gap"),units,mm)
358       
359        // SAScollimation
360        String collimationParent = instrParent + "sascollimation/"
361        // Create SAScollimation entry
362        String collimationBase = instrumentBase + ":sascollimation"
363        NewDataFolder/O/S $(collimationBase)
364        Make/O/T/N=5 $(collimationBase + ":attr") = {"canSAS_class","NX_class"}
365        Make/O/T/N=5 $(collimationBase + ":attrVals") = {"SAScollimation","NXcollimator"}
366        CreateStrNxCansas(fileID,collimationParent,"","",empty,$(collimationBase + ":attr"),$(collimationBase + ":attrVals"))
367        // Create SAScollimation distance entry
368        Make/O/N=1 $(collimationBase + ":distance") = {rw[25]}
369        CreateVarNxCansas(fileID,collimationParent,"sasaperture","distance",$(collimationBase + ":distance"),units,m)
370       
371        // SASdetector
372        String detectorParent = instrParent + "sasdetector/"
373        // Create SASdetector entry
374        String detectorBase = instrumentBase + ":sasdetector"
375        NewDataFolder/O/S $(detectorBase)
376        Make/O/T/N=5 $(detectorBase + ":attr") = {"canSAS_class","NX_class"}
377        Make/O/T/N=5 $(detectorBase + ":attrVals") = {"SASdetector","NXdetector"}
378        CreateStrNxCansas(fileID,detectorParent,"","",empty,$(detectorBase + ":attr"),$(detectorBase + ":attrVals"))
379        // Create SASdetector name entry
380        Make/O/T/N=1 $(detectorBase + ":name") = {textw[9]}
381        CreateStrNxCansas(fileID,detectorParent,"","name",$(detectorBase + ":name"),empty,empty)
382        // Create SASdetector distance entry
383        Make/O/N=1 $(detectorBase + ":SDD") = {rw[18]}
384        CreateVarNxCansas(fileID,detectorParent,"","SDD",$(detectorBase + ":SDD"),units,m)
385        // Create SASdetector beam_center_x entry
386        Make/O/N=1 $(detectorBase + ":beam_center_x") = {rw[16]}
387        CreateVarNxCansas(fileID,detectorParent,"","beam_center_x",$(detectorBase + ":beam_center_x"),units,pixel)
388        // Create SASdetector beam_center_y entry
389        Make/O/N=1 $(detectorBase + ":beam_center_y") = {rw[17]}
390        CreateVarNxCansas(fileID,detectorParent,"","beam_center_y",$(detectorBase + ":beam_center_y"),units,pixel)
391        // Create SASdetector x_pixel_size entry
392        Make/O/N=1 $(detectorBase + ":x_pixel_size") = {rw[10]}
393        CreateVarNxCansas(fileID,detectorParent,"","x_pixel_size",$(detectorBase + ":x_pixel_size"),units,mm)
394        // Create SASdetector y_pixel_size entry
395        Make/O/N=1 $(detectorBase + ":y_pixel_size") = {rw[13]}
396        CreateVarNxCansas(fileID,detectorParent,"","y_pixel_size",$(detectorBase + ":y_pixel_size"),units,mm)
397       
398        // SASsource
399        String sourceParent = instrParent + "sassource/"
400        // Create SASdetector entry
401        String sourceBase = instrumentBase + ":sassource"
402        NewDataFolder/O/S $(sourceBase)
403        Make/O/T/N=5 $(sourceBase + ":attr") = {"canSAS_class","NX_class"}
404        Make/O/T/N=5 $(sourceBase + ":attrVals") = {"SASsource","NXsource"}
405        CreateStrNxCansas(fileID,sourceParent,"","",empty,$(sourceBase + ":attr"),$(sourceBase + ":attrVals"))
406        // Create SASsource radiation entry
407        Make/O/T/N=1 $(sourceBase + ":radiation") = {"Reactor Neutron Source"}
408        CreateStrNxCansas(fileID,sourceParent,"","radiation",$(sourceBase + ":radiation"),empty,empty)
409        // Create SASsource incident_wavelength entry
410        Make/O/N=1 $(sourceBase + ":incident_wavelength") = {rw[26]}
411        CreateVarNxCansas(fileID,sourceParent,"","incident_wavelength",$(sourceBase + ":incident_wavelength"),units,angstrom)
412        // Create SASsource incident_wavelength_spread entry
413        Make/O/N=1 $(sourceBase + ":incident_wavelength_spread") = {rw[27]}
414        CreateVarNxCansas(fileID,sourceParent,"","incident_wavelength_spread",$(sourceBase + ":incident_wavelength_spread"),units,angstrom)
415       
416        // SASsample
417        String sampleParent = parentBase + "sassample/"
418        // Create SASsample entry
419        String sampleBase = base + ":sassample"
420        NewDataFolder/O/S $(sampleBase)
421        Make/O/T/N=5 $(sampleBase + ":attr") = {"canSAS_class","NX_class"}
422        Make/O/T/N=5 $(sampleBase + ":attrVals") = {"SASsample","NXsample"}
423        CreateStrNxCansas(fileID,sampleParent,"","",empty,$(sampleBase + ":attr"),$(sampleBase + ":attrVals"))
424        // Create SASsample name entry
425        Make/O/T/N=1 $(sampleBase + ":name") = {textw[6]}
426        CreateStrNxCansas(fileID,sampleParent,"","name",$(sampleBase + ":name"),empty,empty)
427        // Create SASsample thickness entry
428        Make/O/N=1 $(sampleBase + ":thickness") = {rw[5]}
429        CreateVarNxCansas(fileID,sampleParent,"","thickness",$(sampleBase + ":thickness"),units,cm)
430        // Create SASsample transmission entry
431        Make/O/N=1 $(sampleBase + ":transmission") = {rw[4]}
432        CreateVarNxCansas(fileID,sampleParent,"","transmission",$(sampleBase + ":transmission"),empty,empty)
433End
434       
435//
436///////////////////////////////////////////////////////////////////////////
437
Note: See TracBrowser for help on using the repository browser.