Ignore:
File:
1 edited

Legend:

Unmodified
Added
Removed

  • src/sas/sascalc/dataloader/readers/IgorReader.py

    rdd11014 r9a5097c  
    1313############################################################################# 
    1414import os 
    15  
     15import numpy as np 
     16import math 
     17#import logging 
    1618from sas.sascalc.dataloader.data_info import Data2D 
    1719from sas.sascalc.dataloader.data_info import Detector 
    1820from sas.sascalc.dataloader.manipulations import reader2D_converter 
    19 import numpy as np 
    2021 
    2122# Look for unit converter 
    … …  
    3940        """ Read file """ 
    4041        if not os.path.isfile(filename): 
    41             raise ValueError("Specified file %s is not a regular " 
    42                              "file" % filename) 
    43          
     42            raise ValueError, \ 
     43            "Specified file %s is not a regular file" % filename 
     44         
     45        # Read file 
     46        f = open(filename, 'r') 
     47        buf = f.read() 
     48         
     49        # Instantiate data object 
    4450        output = Data2D() 
    45  
    4651        output.filename = os.path.basename(filename) 
    4752        detector = Detector() 
    48         if len(output.detector): 
    49             print(str(output.detector[0])) 
     53        if len(output.detector) > 0: 
     54            print str(output.detector[0]) 
    5055        output.detector.append(detector) 
    51  
    52         data_conv_q = data_conv_i = None 
    53          
    54         if has_converter and output.Q_unit != '1/A': 
     56                 
     57        # Get content 
     58        dataStarted = False 
     59         
     60        lines = buf.split('\n') 
     61        itot = 0 
     62        x = [] 
     63        y = [] 
     64         
     65        ncounts = 0 
     66         
     67        xmin = None 
     68        xmax = None 
     69        ymin = None 
     70        ymax = None 
     71         
     72        i_x = 0 
     73        i_y = -1 
     74        i_tot_row = 0 
     75         
     76        isInfo = False 
     77        isCenter = False 
     78        
     79        data_conv_q = None 
     80        data_conv_i = None 
     81         
     82        if has_converter == True and output.Q_unit != '1/A': 
    5583            data_conv_q = Converter('1/A') 
    5684            # Test it 
    5785            data_conv_q(1.0, output.Q_unit) 
    5886             
    59         if has_converter and output.I_unit != '1/cm': 
     87        if has_converter == True and output.I_unit != '1/cm': 
    6088            data_conv_i = Converter('1/cm') 
    6189            # Test it 
    6290            data_conv_i(1.0, output.I_unit) 
    63  
    64         data_row = 0 
    65         wavelength = distance = center_x = center_y = None 
    66         dataStarted = isInfo = isCenter = False 
    67  
    68         with open(filename, 'r') as f: 
    69             for line in f: 
    70                 data_row += 1 
    71                 # Find setup info line 
    72                 if isInfo: 
    73                     isInfo = False 
    74                     line_toks = line.split() 
    75                     # Wavelength in Angstrom 
    76                     try: 
    77                         wavelength = float(line_toks[1]) 
    78                     except ValueError: 
    79                         msg = "IgorReader: can't read this file, missing wavelength" 
    80                         raise ValueError(msg) 
    81                     # Distance in meters 
    82                     try: 
    83                         distance = float(line_toks[3]) 
    84                     except ValueError: 
    85                         msg = "IgorReader: can't read this file, missing distance" 
    86                         raise ValueError(msg) 
    87  
    88                     # Distance in meters 
    89                     try: 
    90                         transmission = float(line_toks[4]) 
    91                     except: 
    92                         msg = "IgorReader: can't read this file, " 
    93                         msg += "missing transmission" 
    94                         raise ValueError(msg) 
    95  
    96                 if line.count("LAMBDA"): 
    97                     isInfo = True 
    98  
    99                 # Find center info line 
    100                 if isCenter: 
    101                     isCenter = False 
    102                     line_toks = line.split() 
    103  
    104                     # Center in bin number: Must subtract 1 because 
    105                     # the index starts from 1 
    106                     center_x = float(line_toks[0]) - 1 
    107                     center_y = float(line_toks[1]) - 1 
    108  
    109                 if line.count("BCENT"): 
    110                     isCenter = True 
    111  
    112                 # Find data start 
    113                 if line.count("***"): 
    114                     # now have to continue to blank line 
    115                     dataStarted = True 
    116  
    117                     # Check that we have all the info 
    118                     if (wavelength is None 
    119                             or distance is None 
    120                             or center_x is None 
    121                             or center_y is None): 
    122                         msg = "IgorReader:Missing information in data file" 
    123                         raise ValueError(msg) 
    124  
    125                 if dataStarted: 
    126                     if len(line.rstrip()): 
    127                         continue 
    128                     else: 
    129                         break 
    130  
    131         # The data is loaded in row major order (last index changing most 
    132         # rapidly). However, the original data is in column major order (first 
    133         # index changing most rapidly). The swap to column major order is done 
    134         # in reader2D_converter at the end of this method. 
    135         data = np.loadtxt(filename, skiprows=data_row) 
    136         size_x = size_y = int(np.rint(np.sqrt(data.size))) 
    137         output.data = np.reshape(data, (size_x, size_y)) 
    138         output.err_data = np.zeros_like(output.data) 
    139  
    140         # Det 640 x 640 mm 
    141         # Q = 4 * pi/lambda * sin(theta/2) 
    142         # Bin size is 0.5 cm 
    143         # Removed +1 from theta = (i_x - center_x + 1)*0.5 / distance 
    144         # / 100.0 and 
    145         # Removed +1 from theta = (i_y - center_y + 1)*0.5 / 
    146         # distance / 100.0 
    147         # ToDo: Need  complete check if the following 
    148         # convert process is consistent with fitting.py. 
    149  
    150         # calculate qx, qy bin centers of each pixel in the image 
    151         theta = (np.arange(size_x) - center_x) * 0.5 / distance / 100. 
    152         qx = 4 * np.pi / wavelength * np.sin(theta/2) 
    153  
    154         theta = (np.arange(size_y) - center_y) * 0.5 / distance / 100. 
    155         qy = 4 * np.pi / wavelength * np.sin(theta/2) 
    156  
    157         if has_converter and output.Q_unit != '1/A': 
    158             qx = data_conv_q(qx, units=output.Q_unit) 
    159             qy = data_conv_q(qx, units=output.Q_unit) 
    160  
    161         xmax = np.max(qx) 
    162         xmin = np.min(qx) 
    163         ymax = np.max(qy) 
    164         ymin = np.min(qy) 
    165  
    166         # calculate edge offset in q. 
     91          
     92        for line in lines: 
     93             
     94            # Find setup info line 
     95            if isInfo: 
     96                isInfo = False 
     97                line_toks = line.split() 
     98                # Wavelength in Angstrom 
     99                try: 
     100                    wavelength = float(line_toks[1]) 
     101                except: 
     102                    msg = "IgorReader: can't read this file, missing wavelength" 
     103                    raise ValueError, msg 
     104                 
     105            #Find # of bins in a row assuming the detector is square. 
     106            if dataStarted == True: 
     107                try: 
     108                    value = float(line) 
     109                except: 
     110                    # Found a non-float entry, skip it 
     111                    continue 
     112                 
     113                # Get total bin number 
     114                 
     115            i_tot_row += 1 
     116        i_tot_row = math.ceil(math.sqrt(i_tot_row)) - 1 
     117        #print "i_tot", i_tot_row 
     118        size_x = i_tot_row  # 192#128 
     119        size_y = i_tot_row  # 192#128 
     120        output.data = np.zeros([size_x, size_y]) 
     121        output.err_data = np.zeros([size_x, size_y]) 
     122      
     123        #Read Header and 2D data 
     124        for line in lines: 
     125            # Find setup info line 
     126            if isInfo: 
     127                isInfo = False 
     128                line_toks = line.split() 
     129                # Wavelength in Angstrom 
     130                try: 
     131                    wavelength = float(line_toks[1]) 
     132                except: 
     133                    msg = "IgorReader: can't read this file, missing wavelength" 
     134                    raise ValueError, msg 
     135                # Distance in meters 
     136                try: 
     137                    distance = float(line_toks[3]) 
     138                except: 
     139                    msg = "IgorReader: can't read this file, missing distance" 
     140                    raise ValueError, msg 
     141                 
     142                # Distance in meters 
     143                try: 
     144                    transmission = float(line_toks[4]) 
     145                except: 
     146                    msg = "IgorReader: can't read this file, " 
     147                    msg += "missing transmission" 
     148                    raise ValueError, msg 
     149                                             
     150            if line.count("LAMBDA") > 0: 
     151                isInfo = True 
     152                 
     153            # Find center info line 
     154            if isCenter: 
     155                isCenter = False 
     156                line_toks = line.split() 
     157                 
     158                # Center in bin number: Must substrate 1 because 
     159                #the index starts from 1 
     160                center_x = float(line_toks[0]) - 1 
     161                center_y = float(line_toks[1]) - 1 
     162 
     163            if line.count("BCENT") > 0: 
     164                isCenter = True 
     165                 
     166            # Find data start 
     167            if line.count("***")>0: 
     168                dataStarted = True 
     169                 
     170                # Check that we have all the info 
     171                if wavelength == None \ 
     172                    or distance == None \ 
     173                    or center_x == None \ 
     174                    or center_y == None: 
     175                    msg = "IgorReader:Missing information in data file" 
     176                    raise ValueError, msg 
     177                 
     178            if dataStarted == True: 
     179                try: 
     180                    value = float(line) 
     181                except: 
     182                    # Found a non-float entry, skip it 
     183                    continue 
     184                 
     185                # Get bin number 
     186                if math.fmod(itot, i_tot_row) == 0: 
     187                    i_x = 0 
     188                    i_y += 1 
     189                else: 
     190                    i_x += 1 
     191                     
     192                output.data[i_y][i_x] = value 
     193                ncounts += 1 
     194                 
     195                # Det 640 x 640 mm 
     196                # Q = 4pi/lambda sin(theta/2) 
     197                # Bin size is 0.5 cm  
     198                #REmoved +1 from theta = (i_x-center_x+1)*0.5 / distance 
     199                # / 100.0 and  
     200                #REmoved +1 from theta = (i_y-center_y+1)*0.5 / 
     201                # distance / 100.0 
     202                #ToDo: Need  complete check if the following 
     203                # covert process is consistent with fitting.py. 
     204                theta = (i_x - center_x) * 0.5 / distance / 100.0 
     205                qx = 4.0 * math.pi / wavelength * math.sin(theta/2.0) 
     206 
     207                if has_converter == True and output.Q_unit != '1/A': 
     208                    qx = data_conv_q(qx, units=output.Q_unit) 
     209 
     210                if xmin == None or qx < xmin: 
     211                    xmin = qx 
     212                if xmax == None or qx > xmax: 
     213                    xmax = qx 
     214                 
     215                theta = (i_y - center_y) * 0.5 / distance / 100.0 
     216                qy = 4.0 * math.pi / wavelength * math.sin(theta / 2.0) 
     217 
     218                if has_converter == True and output.Q_unit != '1/A': 
     219                    qy = data_conv_q(qy, units=output.Q_unit) 
     220                 
     221                if ymin == None or qy < ymin: 
     222                    ymin = qy 
     223                if ymax == None or qy > ymax: 
     224                    ymax = qy 
     225                 
     226                if not qx in x: 
     227                    x.append(qx) 
     228                if not qy in y: 
     229                    y.append(qy) 
     230                 
     231                itot += 1 
     232                   
     233                   
    167234        theta = 0.25 / distance / 100.0 
    168         xstep = 4.0 * np.pi / wavelength * np.sin(theta / 2.0) 
     235        xstep = 4.0 * math.pi / wavelength * math.sin(theta / 2.0) 
    169236         
    170237        theta = 0.25 / distance / 100.0 
    171         ystep = 4.0 * np.pi/ wavelength * np.sin(theta / 2.0) 
     238        ystep = 4.0 * math.pi/ wavelength * math.sin(theta / 2.0) 
    172239         
    173240        # Store all data ###################################### 
    174241        # Store wavelength 
    175         if has_converter and output.source.wavelength_unit != 'A': 
     242        if has_converter == True and output.source.wavelength_unit != 'A': 
    176243            conv = Converter('A') 
    177244            wavelength = conv(wavelength, units=output.source.wavelength_unit) 
    … …  
    179246 
    180247        # Store distance 
    181         if has_converter and detector.distance_unit != 'm': 
     248        if has_converter == True and detector.distance_unit != 'm': 
    182249            conv = Converter('m') 
    183250            distance = conv(distance, units=detector.distance_unit) 
    … …  
    187254        output.sample.transmission = transmission 
    188255         
    189         # Store pixel size (mm) 
     256        # Store pixel size 
    190257        pixel = 5.0 
    191         if has_converter and detector.pixel_size_unit != 'mm': 
     258        if has_converter == True and detector.pixel_size_unit != 'mm': 
    192259            conv = Converter('mm') 
    193260            pixel = conv(pixel, units=detector.pixel_size_unit) 
    … …  
    200267         
    201268        # Store limits of the image (2D array) 
    202         xmin -= xstep / 2.0 
    203         xmax += xstep / 2.0 
    204         ymin -= ystep / 2.0 
    205         ymax += ystep / 2.0 
    206         if has_converter and output.Q_unit != '1/A': 
     269        xmin = xmin - xstep / 2.0 
     270        xmax = xmax + xstep / 2.0 
     271        ymin = ymin - ystep / 2.0 
     272        ymax = ymax + ystep / 2.0 
     273        if has_converter == True and output.Q_unit != '1/A': 
    207274            xmin = data_conv_q(xmin, units=output.Q_unit) 
    208275            xmax = data_conv_q(xmax, units=output.Q_unit) 
    … …  
    215282         
    216283        # Store x and y axis bin centers 
    217         output.x_bins = qx.tolist() 
    218         output.y_bins = qy.tolist() 
     284        output.x_bins = x 
     285        output.y_bins = y 
    219286         
    220287        # Units 
Note: See TracChangeset for help on using the changeset viewer.