source: sasview/src/sas/sascalc/dataloader/readers/danse_reader.py @ b9d74f3

Last change on this file since b9d74f3 was b9d74f3, checked in by andyfaff, 8 years ago

MAINT: use raise Exception() not raise Exception

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