1 | """ |
---|
2 | DANSE/SANS file reader |
---|
3 | """ |
---|
4 | ############################################################################ |
---|
5 | #This software was developed by the University of Tennessee as part of the |
---|
6 | #Distributed Data Analysis of Neutron Scattering Experiments (DANSE) |
---|
7 | #project funded by the US National Science Foundation. |
---|
8 | #If you use DANSE applications to do scientific research that leads to |
---|
9 | #publication, we ask that you acknowledge the use of the software with the |
---|
10 | #following sentence: |
---|
11 | #This work benefited from DANSE software developed under NSF award DMR-0520547. |
---|
12 | #copyright 2008, University of Tennessee |
---|
13 | ############################################################################# |
---|
14 | import math |
---|
15 | import os |
---|
16 | import sys |
---|
17 | import numpy as np |
---|
18 | import logging |
---|
19 | from sas.sascalc.dataloader.data_info import Data2D, Detector |
---|
20 | from sas.sascalc.dataloader.manipulations import reader2D_converter |
---|
21 | |
---|
22 | # Look for unit converter |
---|
23 | has_converter = True |
---|
24 | try: |
---|
25 | from sas.sascalc.data_util.nxsunit import Converter |
---|
26 | except: |
---|
27 | has_converter = False |
---|
28 | |
---|
29 | |
---|
30 | class Reader: |
---|
31 | """ |
---|
32 | Example data manipulation |
---|
33 | """ |
---|
34 | ## File type |
---|
35 | type_name = "DANSE" |
---|
36 | ## Wildcards |
---|
37 | type = ["DANSE files (*.sans)|*.sans"] |
---|
38 | ## Extension |
---|
39 | ext = ['.sans', '.SANS'] |
---|
40 | |
---|
41 | def read(self, filename=None): |
---|
42 | """ |
---|
43 | Open and read the data in a file |
---|
44 | @param file: path of the file |
---|
45 | """ |
---|
46 | |
---|
47 | read_it = False |
---|
48 | for item in self.ext: |
---|
49 | if filename.lower().find(item) >= 0: |
---|
50 | read_it = True |
---|
51 | |
---|
52 | if read_it: |
---|
53 | try: |
---|
54 | datafile = open(filename, 'r') |
---|
55 | except: |
---|
56 | raise RuntimeError,"danse_reader cannot open %s" % (filename) |
---|
57 | |
---|
58 | # defaults |
---|
59 | # wavelength in Angstrom |
---|
60 | wavelength = 10.0 |
---|
61 | # Distance in meter |
---|
62 | distance = 11.0 |
---|
63 | # Pixel number of center in x |
---|
64 | center_x = 65 |
---|
65 | # Pixel number of center in y |
---|
66 | center_y = 65 |
---|
67 | # Pixel size [mm] |
---|
68 | pixel = 5.0 |
---|
69 | # Size in x, in pixels |
---|
70 | size_x = 128 |
---|
71 | # Size in y, in pixels |
---|
72 | size_y = 128 |
---|
73 | # Format version |
---|
74 | fversion = 1.0 |
---|
75 | |
---|
76 | output = Data2D() |
---|
77 | output.filename = os.path.basename(filename) |
---|
78 | detector = Detector() |
---|
79 | output.detector.append(detector) |
---|
80 | |
---|
81 | output.data = np.zeros([size_x,size_y]) |
---|
82 | output.err_data = np.zeros([size_x, size_y]) |
---|
83 | |
---|
84 | data_conv_q = None |
---|
85 | data_conv_i = None |
---|
86 | |
---|
87 | if has_converter == True and output.Q_unit != '1/A': |
---|
88 | data_conv_q = Converter('1/A') |
---|
89 | # Test it |
---|
90 | data_conv_q(1.0, output.Q_unit) |
---|
91 | |
---|
92 | if has_converter == True and output.I_unit != '1/cm': |
---|
93 | data_conv_i = Converter('1/cm') |
---|
94 | # Test it |
---|
95 | data_conv_i(1.0, output.I_unit) |
---|
96 | |
---|
97 | read_on = True |
---|
98 | while read_on: |
---|
99 | line = datafile.readline() |
---|
100 | if line.find("DATA:") >= 0: |
---|
101 | read_on = False |
---|
102 | break |
---|
103 | toks = line.split(':') |
---|
104 | if toks[0] == "FORMATVERSION": |
---|
105 | fversion = float(toks[1]) |
---|
106 | if toks[0] == "WAVELENGTH": |
---|
107 | wavelength = float(toks[1]) |
---|
108 | elif toks[0] == "DISTANCE": |
---|
109 | distance = float(toks[1]) |
---|
110 | elif toks[0] == "CENTER_X": |
---|
111 | center_x = float(toks[1]) |
---|
112 | elif toks[0] == "CENTER_Y": |
---|
113 | center_y = float(toks[1]) |
---|
114 | elif toks[0] == "PIXELSIZE": |
---|
115 | pixel = float(toks[1]) |
---|
116 | elif toks[0] == "SIZE_X": |
---|
117 | size_x = int(toks[1]) |
---|
118 | elif toks[0] == "SIZE_Y": |
---|
119 | size_y = int(toks[1]) |
---|
120 | |
---|
121 | # Read the data |
---|
122 | data = [] |
---|
123 | error = [] |
---|
124 | if fversion == 1.0: |
---|
125 | data_str = datafile.readline() |
---|
126 | data = data_str.split(' ') |
---|
127 | else: |
---|
128 | read_on = True |
---|
129 | while read_on: |
---|
130 | data_str = datafile.readline() |
---|
131 | if len(data_str) == 0: |
---|
132 | read_on = False |
---|
133 | else: |
---|
134 | toks = data_str.split() |
---|
135 | try: |
---|
136 | val = float(toks[0]) |
---|
137 | err = float(toks[1]) |
---|
138 | if data_conv_i is not None: |
---|
139 | val = data_conv_i(val, units=output._yunit) |
---|
140 | err = data_conv_i(err, units=output._yunit) |
---|
141 | data.append(val) |
---|
142 | error.append(err) |
---|
143 | except: |
---|
144 | logging.info("Skipping line:%s,%s" %(data_str, |
---|
145 | sys.exc_value)) |
---|
146 | |
---|
147 | # Initialize |
---|
148 | x_vals = [] |
---|
149 | y_vals = [] |
---|
150 | ymin = None |
---|
151 | ymax = None |
---|
152 | xmin = None |
---|
153 | xmax = None |
---|
154 | |
---|
155 | # Qx and Qy vectors |
---|
156 | theta = pixel / distance / 100.0 |
---|
157 | stepq = 4.0 * math.pi / wavelength * math.sin(theta / 2.0) |
---|
158 | for i_x in range(size_x): |
---|
159 | theta = (i_x - center_x + 1) * pixel / distance / 100.0 |
---|
160 | qx = 4.0 * math.pi / wavelength * math.sin(theta / 2.0) |
---|
161 | |
---|
162 | if has_converter == True and output.Q_unit != '1/A': |
---|
163 | qx = data_conv_q(qx, units=output.Q_unit) |
---|
164 | |
---|
165 | x_vals.append(qx) |
---|
166 | if xmin == None or qx < xmin: |
---|
167 | xmin = qx |
---|
168 | if xmax == None or qx > xmax: |
---|
169 | xmax = qx |
---|
170 | |
---|
171 | ymin = None |
---|
172 | ymax = None |
---|
173 | for i_y in range(size_y): |
---|
174 | theta = (i_y - center_y + 1) * pixel / distance / 100.0 |
---|
175 | qy = 4.0 * math.pi / wavelength * math.sin(theta/2.0) |
---|
176 | |
---|
177 | if has_converter == True and output.Q_unit != '1/A': |
---|
178 | qy = data_conv_q(qy, units=output.Q_unit) |
---|
179 | |
---|
180 | y_vals.append(qy) |
---|
181 | if ymin == None or qy < ymin: |
---|
182 | ymin = qy |
---|
183 | if ymax == None or qy > ymax: |
---|
184 | ymax = qy |
---|
185 | |
---|
186 | # Store the data in the 2D array |
---|
187 | i_x = 0 |
---|
188 | i_y = -1 |
---|
189 | |
---|
190 | for i_pt in range(len(data)): |
---|
191 | try: |
---|
192 | value = float(data[i_pt]) |
---|
193 | except: |
---|
194 | # For version 1.0, the data were still |
---|
195 | # stored as strings at this point. |
---|
196 | msg = "Skipping entry (v1.0):%s,%s" % (str(data[i_pt]), |
---|
197 | sys.exc_value) |
---|
198 | logging.info(msg) |
---|
199 | |
---|
200 | # Get bin number |
---|
201 | if math.fmod(i_pt, size_x) == 0: |
---|
202 | i_x = 0 |
---|
203 | i_y += 1 |
---|
204 | else: |
---|
205 | i_x += 1 |
---|
206 | |
---|
207 | output.data[i_y][i_x] = value |
---|
208 | if fversion>1.0: |
---|
209 | output.err_data[i_y][i_x] = error[i_pt] |
---|
210 | |
---|
211 | # Store all data |
---|
212 | # Store wavelength |
---|
213 | if has_converter == True and output.source.wavelength_unit != 'A': |
---|
214 | conv = Converter('A') |
---|
215 | wavelength = conv(wavelength, |
---|
216 | units=output.source.wavelength_unit) |
---|
217 | output.source.wavelength = wavelength |
---|
218 | |
---|
219 | # Store distance |
---|
220 | if has_converter == True and detector.distance_unit != 'm': |
---|
221 | conv = Converter('m') |
---|
222 | distance = conv(distance, units=detector.distance_unit) |
---|
223 | detector.distance = distance |
---|
224 | |
---|
225 | # Store pixel size |
---|
226 | if has_converter == True and detector.pixel_size_unit != 'mm': |
---|
227 | conv = Converter('mm') |
---|
228 | pixel = conv(pixel, units=detector.pixel_size_unit) |
---|
229 | detector.pixel_size.x = pixel |
---|
230 | detector.pixel_size.y = pixel |
---|
231 | |
---|
232 | # Store beam center in distance units |
---|
233 | detector.beam_center.x = center_x * pixel |
---|
234 | detector.beam_center.y = center_y * pixel |
---|
235 | |
---|
236 | # Store limits of the image (2D array) |
---|
237 | xmin = xmin - stepq / 2.0 |
---|
238 | xmax = xmax + stepq / 2.0 |
---|
239 | ymin = ymin - stepq /2.0 |
---|
240 | ymax = ymax + stepq / 2.0 |
---|
241 | |
---|
242 | if has_converter == True and output.Q_unit != '1/A': |
---|
243 | xmin = data_conv_q(xmin, units=output.Q_unit) |
---|
244 | xmax = data_conv_q(xmax, units=output.Q_unit) |
---|
245 | ymin = data_conv_q(ymin, units=output.Q_unit) |
---|
246 | ymax = data_conv_q(ymax, units=output.Q_unit) |
---|
247 | output.xmin = xmin |
---|
248 | output.xmax = xmax |
---|
249 | output.ymin = ymin |
---|
250 | output.ymax = ymax |
---|
251 | |
---|
252 | # Store x and y axis bin centers |
---|
253 | output.x_bins = x_vals |
---|
254 | output.y_bins = y_vals |
---|
255 | |
---|
256 | # Units |
---|
257 | if data_conv_q is not None: |
---|
258 | output.xaxis("\\rm{Q_{x}}", output.Q_unit) |
---|
259 | output.yaxis("\\rm{Q_{y}}", output.Q_unit) |
---|
260 | else: |
---|
261 | output.xaxis("\\rm{Q_{x}}", 'A^{-1}') |
---|
262 | output.yaxis("\\rm{Q_{y}}", 'A^{-1}') |
---|
263 | |
---|
264 | if data_conv_i is not None: |
---|
265 | output.zaxis("\\rm{Intensity}", output.I_unit) |
---|
266 | else: |
---|
267 | output.zaxis("\\rm{Intensity}", "cm^{-1}") |
---|
268 | |
---|
269 | if not fversion >= 1.0: |
---|
270 | msg = "Danse_reader can't read this file %s" % filename |
---|
271 | raise ValueError, msg |
---|
272 | else: |
---|
273 | logging.info("Danse_reader Reading %s \n" % filename) |
---|
274 | |
---|
275 | # Store loading process information |
---|
276 | output.meta_data['loader'] = self.type_name |
---|
277 | output = reader2D_converter(output) |
---|
278 | return output |
---|
279 | |
---|
280 | return None |
---|