source: sasview/DataLoader/readers/red2d_reader.py @ 67e258c

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.1.1release-4.1.2release-4.2.2release_4.0.1ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload
Last change on this file since 67e258c was d2539aa, checked in by Jae Cho <jhjcho@…>, 14 years ago

reset the value to zero if nan or inf

  • Property mode set to 100644
File size: 10.9 KB
Line 
1
2#####################################################################
3#This software was developed by the University of Tennessee as part of the
4#Distributed Data Analysis of Neutron Scattering Experiments (DANSE)
5#project funded by the US National Science Foundation.
6#See the license text in license.txt
7#copyright 2008, University of Tennessee
8######################################################################
9
10"""
11    TXT/IGOR 2D Q Map file reader
12"""
13
14
15import os, sys
16import numpy
17import math, logging
18from DataLoader.data_info import Data2D, Detector
19
20# Look for unit converter
21has_converter = True
22try:
23    from data_util.nxsunit import Converter
24except:
25    has_converter = False
26   
27def check_point(x_point):
28    """
29    check point validity
30    """
31    # set zero for non_floats
32    try:
33        return float(x_point)
34    except:
35        return 0
36       
37class Reader:
38    """ Simple data reader for Igor data files """
39    ## File type
40    type_name = "IGOR/DAT 2D Q_map"   
41    ## Wildcards
42    type = ["IGOR/DAT 2D file in Q_map (*.dat)|*.DAT"]
43    ## Extension
44    ext=['.DAT', '.dat']
45
46    def read(self,filename=None):
47        """ Read file """
48        if not os.path.isfile(filename):
49            raise ValueError, \
50            "Specified file %s is not a regular file" % filename
51
52        # Read file
53        f = open(filename,'r')
54        buf = f.read()
55        f.close()     
56        # Instantiate data object
57        output = Data2D()
58        output.filename = os.path.basename(filename)
59        detector = Detector()
60        if len(output.detector)>0: print str(output.detector[0])
61        output.detector.append(detector)
62               
63        # Get content
64        dataStarted = False
65       
66        ## Defaults     
67        lines = buf.split('\n')
68        itot = 0
69        x = []
70        y = []
71       
72        ncounts = 0
73       
74        wavelength   = None
75        distance     = None
76        transmission = None
77       
78        pixel_x = None
79        pixel_y = None
80       
81        i_x    = 0
82        i_y    = -1
83        pixels = 0
84       
85        isInfo   = False
86        isCenter = False
87
88        data_conv_q = None
89        data_conv_i = None
90       
91        # Set units: This is the unit assumed for Q and I in the data file.
92        if has_converter == True and output.Q_unit != '1/A':
93            data_conv_q = Converter('1/A')
94            # Test it
95            data_conv_q(1.0, output.Q_unit)
96           
97        if has_converter == True and output.I_unit != '1/cm':
98            data_conv_i = Converter('1/cm')
99            # Test it
100            data_conv_i(1.0, output.I_unit)           
101       
102             
103        # Remove the last lines before the for loop if the lines are empty
104        # to calculate the exact number of data points
105        count = 0
106        while (len(lines[len(lines)-(count+1)].lstrip().rstrip()) < 1):
107            del lines[len(lines)-(count+1)]
108            count = count + 1
109
110        #Read Header and find the dimensions of 2D data
111        line_num = 0
112        for line in lines:     
113            line_num += 1
114            ## Reading the header applies only to IGOR/NIST 2D q_map data files
115            # Find setup info line
116            if isInfo:
117                isInfo = False
118                line_toks = line.split()
119                # Wavelength in Angstrom
120                try:
121                    wavelength = float(line_toks[1])
122                    # Units
123                    if has_converter==True and output.source.wavelength_unit != 'A':
124                        conv = Converter('A')
125                        wavelength = conv(wavelength, units=output.source.wavelength_unit)
126                except:
127                    #Not required
128                    pass
129                # Distance in mm
130                try:
131                    distance = float(line_toks[3])
132                    # Units
133                    if has_converter==True and detector.distance_unit != 'm':
134                        conv = Converter('m')
135                        distance = conv(distance, units=detector.distance_unit)
136                except:
137                    #Not required
138                    pass
139               
140                # Distance in meters
141                try:
142                    transmission = float(line_toks[4])
143                except:
144                    #Not required
145                    pass
146                                           
147            if line.count("LAMBDA")>0:
148                isInfo = True
149               
150            # Find center info line
151            if isCenter:
152                isCenter = False               
153                line_toks = line.split()
154                # Center in bin number
155                center_x = float(line_toks[0])
156                center_y = float(line_toks[1])
157
158            if line.count("BCENT")>0:
159                isCenter = True
160                       
161            # Find data start
162            if line.count("Data columns") or line.count("ASCII data")>0:
163                dataStarted = True
164                continue
165
166            ## Read and get data.   
167            if dataStarted == True:
168                line_toks = line.split()             
169                if len(line_toks) == 0:
170                    #empty line
171                    continue
172                # the number of columns must be stayed same
173                col_num = len(line_toks)
174                break
175
176       
177        # Make numpy array to remove header lines using index
178        lines_array = numpy.array(lines)
179
180        # index for lines_array
181        lines_index = numpy.arange(len(lines))
182       
183        # get the data lines
184        data_lines = lines_array[lines_index>=(line_num-1)]
185        # Now we get the total number of rows (i.e., # of data points)
186        row_num = len(data_lines)
187        # make it as list again to control the separators
188        data_list = " ".join(data_lines.tolist())
189        # split all data to one big list w/" "separator
190        data_list = data_list.split()
191 
192        # Check if the size is consistent with data, otherwise try the tab(\t) separator
193        # (this may be removed once get the confidence the former working all cases).
194        if len(data_list) != (len(data_lines)) * col_num:
195            data_list = "\t".join(data_lines.tolist())
196            data_list = data_list.split()
197
198        # Change it(string) into float
199        #data_list = map(float,data_list)
200        data_list1 = map(check_point,data_list)
201     
202        # numpy array form
203        data_array = numpy.array(data_list1)
204        # Redimesion based on the row_num and col_num, otherwise raise an error.
205        try:
206            data_point = data_array.reshape(row_num,col_num).transpose()
207        except:
208            raise ValueError, "red2d_reader: Can't read this file: Not a proper file format"
209       
210        ## Get the all data: Let's HARDcoding; Todo find better way
211        # Defaults
212        dqx_data = numpy.zeros(0)
213        dqy_data = numpy.zeros(0)
214        qz_data = numpy.zeros(row_num)
215        mask = numpy.ones(row_num,dtype=bool)
216        # Get from the array
217        qx_data = data_point[0]
218        qy_data = data_point[1]
219        data = data_point[2]
220        if col_num >3: qz_data = data_point[3]
221        if col_num >4: dqx_data = data_point[4]
222        if col_num >5: dqy_data = data_point[5]
223        if col_num >6: mask[data_point[6]<1] = False
224        q_data = numpy.sqrt(qx_data*qx_data+qy_data*qy_data+qz_data*qz_data)
225           
226        # Extra protection(it is needed for some data files):
227        # If all mask elements are False, put all True
228        if not mask.any(): mask[mask==False] = True   
229           
230        # Store limits of the image in q space
231        xmin    = numpy.min(qx_data)
232        xmax    = numpy.max(qx_data)
233        ymin    = numpy.min(qy_data)
234        ymax    = numpy.max(qy_data)
235
236        # units
237        if has_converter == True and output.Q_unit != '1/A':
238            xmin = data_conv_q(xmin, units=output.Q_unit)
239            xmax = data_conv_q(xmax, units=output.Q_unit)
240            ymin = data_conv_q(ymin, units=output.Q_unit)
241            ymax = data_conv_q(ymax, units=output.Q_unit)
242           
243        ## calculate the range of the qx and qy_data
244        x_size = math.fabs(xmax - xmin)
245        y_size = math.fabs(ymax - ymin)
246       
247        # calculate the number of pixels in the each axes
248        npix_y = math.floor(math.sqrt(len(data)))
249        npix_x = math.floor(len(data)/npix_y)
250       
251        # calculate the size of bins     
252        xstep = x_size/(npix_x-1)
253        ystep = y_size/(npix_y-1)
254       
255        # store x and y axis bin centers in q space
256        x_bins  = numpy.arange(xmin,xmax+xstep,xstep)
257        y_bins  = numpy.arange(ymin,ymax+ystep,ystep)
258       
259        # get the limits of q values
260        xmin = xmin - xstep/2
261        xmax = xmax + xstep/2
262        ymin = ymin - ystep/2
263        ymax = ymax + ystep/2
264       
265        #Store data in outputs 
266        #TODO: Check the lengths
267        output.data     = data
268        output.err_data = numpy.sqrt(numpy.abs(data))
269        output.qx_data  = qx_data
270        output.qy_data  = qy_data             
271        output.q_data   = q_data
272        output.mask     = mask
273       
274        output.x_bins = x_bins
275        output.y_bins = y_bins
276               
277        output.xmin = xmin
278        output.xmax = xmax
279        output.ymin = ymin
280        output.ymax = ymax
281       
282        output.source.wavelength = wavelength
283       
284        # Store pixel size in mm
285        detector.pixel_size.x = pixel_x
286        detector.pixel_size.y = pixel_y
287       
288        # Store the sample to detector distance
289        detector.distance = distance
290       
291        # optional data: if all of dq data == 0, do not pass to output
292        if len(dqx_data) == len(qx_data) and dqx_data.any()!=0: 
293            # if no dqx_data, do not pass dqy_data.(1 axis dq is not supported yet).
294            if len(dqy_data) == len(qy_data) and dqy_data.any()!=0:
295                output.dqx_data = dqx_data
296                output.dqy_data = dqy_data
297       
298        # Units of axes
299        if data_conv_q is not None:
300            output.xaxis("\\rm{Q_{x}}", output.Q_unit)
301            output.yaxis("\\rm{Q_{y}}", output.Q_unit)
302        else:
303            output.xaxis("\\rm{Q_{x}}", 'A^{-1}')
304            output.yaxis("\\rm{Q_{y}}", 'A^{-1}')           
305        if data_conv_i is not None:
306            output.zaxis("\\rm{Intensity}", output.I_unit)
307        else:
308            output.zaxis("\\rm{Intensity}","cm^{-1}")
309   
310        # Store loading process information
311        output.meta_data['loader'] = self.type_name
312
313        return output
314   
315if __name__ == "__main__": 
316    reader = Reader()
317    print reader.read("../test/exp18_14_igor_2dqxqy.dat") 
318       
319
Note: See TracBrowser for help on using the repository browser.