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

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

Incorporate suggested changes for cleaner separation of VSANS, USANS, and SANS reduction methods. Plus fixes.

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