source: sasview/DataLoader/readers/red2d_reader.py @ 095ab1b

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 095ab1b was 3cd95c8, checked in by Jae Cho <jhjcho@…>, 15 years ago

new sans qmap 2D reader that enables reading TOF sans data as well as some other type of detector data

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