This is an archival dump of old wiki content --- see scipy.org for current material.
Please see http://scipy-cookbook.readthedocs.org/

This cookbook example contains a module that implements a reader for a LAS (Log ASCII Standard) well log file (LAS 2.0). See the Canadian Well Logging Society page about this format for more information.

   1 """LAS File Reader
   2 
   3 The main class defined here is LASReader, a class that reads a LAS file
   4 and makes the data available as a Python object.
   5 """
   6 
   7 # Copyright (c) 2011, Warren Weckesser
   8 #
   9 # Permission to use, copy, modify, and/or distribute this software for any
  10 # purpose with or without fee is hereby granted, provided that the above
  11 # copyright notice and this permission notice appear in all copies.
  12 #
  13 # THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
  14 # WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
  15 # MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
  16 # ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
  17 # WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
  18 # ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
  19 # OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
  20 
  21 
  22 import re
  23 import keyword
  24 
  25 import numpy as np
  26 
  27 
  28 def isidentifier(s):
  29     if s in keyword.kwlist:
  30         return False
  31     return re.match(r'^[a-z_][a-z0-9_]*$', s, re.I) is not None
  32 
  33 
  34 def _convert_to_value(s):
  35     try:
  36         value = int(s)
  37     except ValueError:
  38         try:
  39             value = float(s)
  40         except ValueError:
  41             value = s
  42     return value
  43 
  44 
  45 class LASError(Exception):
  46     pass
  47 
  48 
  49 class LASItem(object):
  50     """This class is just a namespace, holding the attributes 'name',
  51     'units', 'data', 'value', and 'descr'.  'value' is the numerical
  52     value of 'data', if it has a numerical value (specifically, if
  53     int() or float() don't raise an exception when given the value
  54     of the 'data' attribute).
  55 
  56     A class method, from_line(cls, line), is provided to parse
  57     a line from a LAS file and create a LASItem instance.
  58     """
  59     def __init__(self, name, units='', data='', descr=''):
  60         self.name = name
  61         self.units = units
  62         self.data = data
  63         self.value = _convert_to_value(data)
  64         self.descr = descr
  65 
  66     def __str__(self):
  67         s = ("name='%s', units='%s', data='%s', descr='%s'" %
  68                 (self.name, self.units, self.data, self.descr))
  69         return s
  70 
  71     def __repr__(self):
  72         s = str(self)
  73         return "LASItem(%s)" % s
  74 
  75     @classmethod
  76     def from_line(cls, line):
  77         first, descr = line.rsplit(':', 1)
  78         descr = descr.strip()
  79         name, mid = first.split('.', 1)
  80         name = name.strip()
  81         if mid.startswith(' '):
  82             # No units
  83             units = ''
  84             data = mid
  85         else:
  86             units_data = mid.split(None, 1)
  87             if len(units_data) == 1:
  88                 units = units_data[0]
  89                 data = ''
  90             else:
  91                 units, data = units_data
  92         return LASItem(name=name, units=units, data=data.strip(),
  93                        descr=descr.strip())
  94 
  95 
  96 def _read_wrapped_row(f, n):
  97     """Read a "row" of data from the Ascii section of a "wrapped" LAS file.
  98 
  99     `f` must be a file object opened for reading.
 100     `n` is the number of fields in the row.
 101 
 102     Returns the list of floats read from the file.
 103     """
 104     depth = float(f.readline().strip())
 105     values = [depth]
 106     while len(values) < n:
 107         new_values = [float(s) for s in f.readline().split()]
 108         values.extend(new_values)
 109     return values
 110 
 111 
 112 def _read_wrapped_data(f, dt):
 113     data = []
 114     ncols = len(dt.names)
 115     while True:
 116         try:
 117             row = _read_wrapped_row(f, ncols)
 118         except Exception:
 119             break
 120         data.append(tuple(row))
 121     data = np.array(data, dtype=dt)
 122     return data
 123 
 124 
 125 class LASSection(object):
 126     """Represents a "section" of a LAS file.
 127 
 128     A section is basically a collection of items, where each item has the
 129     attributes 'name', 'units', 'data' and 'descr'.
 130 
 131     Any item in the section whose name is a valid Python identifier is
 132     also attached to the object as an attribute.  For example, if `s` is a
 133     LASSection instance, and the corresponding section in the LAS file
 134     contained this line:
 135 
 136     FD   .K/M3               999.9999        :  Fluid Density
 137 
 138     then the item may be referred to as `s.FD` (in addition to the longer
 139     `s.items['FD']`).
 140 
 141     Attributes
 142     ----------
 143     items : dict
 144         The keys are the item names, and the values are LASItem instances.
 145     names : list
 146         List of item names, in the order they were read from the LAS file.
 147 
 148     """
 149     def __init__(self):
 150         # Note: In Python 2.7, 'items' could be an OrderedDict, and
 151         # then 'names' would not be necessary--one could use items.keys().
 152         self.items = dict()
 153         self.names = []
 154 
 155     def add_item(self, item):
 156         self.items[item.name] = item
 157         self.names.append(item.name)
 158         if isidentifier(item.name) and not hasattr(self, item.name):
 159             setattr(self, item.name, item)
 160 
 161     def display(self):
 162         for name in self.names:
 163             item = self.items[name]
 164             namestr = name
 165             if item.units != '':
 166                 namestr = namestr + (" (%s)" % item.units)
 167             print "%-16s %-30s [%s]" % (namestr, "'" + item.data + "'",
 168                                         item.descr)
 169 
 170 
 171 class LASReader(object):
 172     """The LASReader class holds data from a LAS file.
 173 
 174     This reader only handles LAS 2.0 files (as far as I know).
 175 
 176     Constructor
 177     -----------
 178     LASReader(f, null_subs=None)
 179 
 180     f : file object or string
 181         If f is a file object, it must be opened for reading.
 182         If f is a string, it must be the filename of a LAS file.
 183         In that case, the file will be opened and read.
 184 
 185     Attributes for LAS Sections
 186     ---------------------------
 187     version : LASSection instance
 188         This LASSection holds the items from the '~V' section.
 189 
 190     well : LASSection instance
 191         This LASSection holds the items from the '~W' section.
 192 
 193     curves : LASection instance
 194         This LASSection holds the items from the '~C' section.
 195 
 196     parameters : LASSection instance
 197         This LASSection holds the items from the '~P' section.
 198 
 199     other : str
 200         Holds the contents of the '~O' section as a single string.
 201 
 202     data : numpy 1D structured array
 203         The numerical data from the '~A' section.  The data type
 204         of the array is constructed from the items in the '~C'
 205         section.
 206 
 207     Other attributes
 208     ----------------
 209     data2d : numpy 2D array of floats
 210         The numerical data from the '~A' section, as a 2D array.
 211         This is a view of the same data as in the `data` attribute.
 212 
 213     wrap : bool
 214         True if the LAS file was wrapped. (More specifically, this
 215         attribute is True if the data field of the item with the
 216         name 'WRAP' in the '~V' section has the value 'YES'.)
 217 
 218     vers : str
 219         The LAS version. (More specifically, the value of the data
 220         field of the item with the name 'VERS' in the '~V' section).
 221 
 222     null : float or None
 223         The numerical value of the 'NULL' item in the '~W' section.
 224         The value will be None if the 'NULL' item was missing.
 225 
 226     null_subs : float or None
 227         The value given in the constructor, to be used as the
 228         replacement value of each occurrence of `null_value` in
 229         the log data.  The value will be None (and no substitution
 230         will be done) if the `null_subs` argument is not given to
 231         the constructor.
 232 
 233     start : float, or None
 234         Numerical value of the 'STRT' item from the '~W' section.
 235         The value will be None if 'STRT' was not given in the file.
 236 
 237     start_units : str
 238         Units of the 'STRT' item from the '~W' section.
 239         The value will be None if 'STRT' was not given in the file.
 240 
 241     stop : float
 242         Numerical value of the 'STOP' item from the '~W' section.
 243         The value will be None if 'STOP' was not given in the file.
 244 
 245     stop_units : str
 246         Units of the 'STOP' item from the '~W' section.
 247         The value will be None if 'STOP' was not given in the file.
 248 
 249     step : float
 250         Numerical value of the 'STEP' item from the '~W' section.
 251         The value will be None if 'STEP' was not given in the file.
 252 
 253     step_units : str
 254         Units of the 'STEP' item from the '~W' section.
 255         The value will be None if 'STEP' was not given in the file.
 256 
 257     """
 258 
 259     def __init__(self, f, null_subs=None):
 260         """f can be a filename (str) or a file object.
 261 
 262         If 'null_subs' is not None, its value replaces any values in the data
 263         that matches the NULL value specified in the Version section of the LAS
 264         file.
 265         """
 266         self.null = None
 267         self.null_subs = null_subs
 268         self.start = None
 269         self.start_units = None
 270         self.stop = None
 271         self.stop_units = None
 272         self.step = None
 273         self.step_units = None
 274 
 275         self.version = LASSection()
 276         self.well = LASSection()
 277         self.curves = LASSection()
 278         self.parameters = LASSection()
 279         self.other = ''
 280         self.data = None
 281 
 282         self._read_las(f)
 283 
 284         self.data2d = self.data.view(float).reshape(-1, len(self.curves.items))
 285         if null_subs is not None:
 286             self.data2d[self.data2d == self.null] = null_subs
 287 
 288     def _read_las(self, f):
 289         """Read a LAS file.
 290 
 291         Returns a dictionary with keys 'V', 'W', 'C', 'P', 'O' and 'A',
 292         corresponding to the sections of a LAS file.  The values associated
 293         with keys 'V', 'W', 'C' and 'P' will be lists of Item instances.  The
 294         value associated with the 'O' key is a list of strings.  The value
 295         associated with the 'A' key is a numpy structured array containing the
 296         log data.  The field names of the array are the mnemonics from the
 297         Curve section of the file.
 298         """
 299         opened_here = False
 300         if isinstance(f, basestring):
 301             opened_here = True
 302             f = open(f, 'r')
 303 
 304         self.wrap = False
 305 
 306         line = f.readline()
 307         current_section = None
 308         current_section_label = ''
 309         while not line.startswith('~A'):
 310             if not line.startswith('#'):
 311                 if line.startswith('~'):
 312                     if len(line) < 2:
 313                         raise LASError("Missing section character after '~'.")
 314                     current_section_label = line[1:2]
 315                     other = False
 316                     if current_section_label == 'V':
 317                         current_section = self.version
 318                     elif current_section_label == 'W':
 319                         current_section = self.well
 320                     elif current_section_label == 'C':
 321                         current_section = self.curves
 322                     elif current_section_label == 'P':
 323                         current_section = self.parameters
 324                     elif current_section_label == 'O':
 325                         current_section = self.other
 326                         other = True
 327                     else:
 328                         raise LASError("Unknown section '%s'" % line)
 329                 elif current_section is None:
 330                     raise LASError("Missing first section.")
 331                 else:
 332                     if other:
 333                         # The 'Other' section is just lines of text, so we
 334                         # assemble them into a single string.
 335                         self.other += line
 336                         current_section = self.other
 337                     else:
 338                         # Parse the line into a LASItem and add it to the
 339                         # current section.
 340                         m = LASItem.from_line(line)
 341                         current_section.add_item(m)
 342                         # Check for the required items whose values we'll
 343                         # store as attributes of the LASReader instance.
 344                         if current_section == self.version:
 345                             if m.name == 'WRAP':
 346                                 if m.data.strip() == 'YES':
 347                                     self.wrap = True
 348                             if m.name == 'VERS':
 349                                 self.vers = m.data.strip()
 350                         if current_section == self.well:
 351                             if m.name == 'NULL':
 352                                 self.null = float(m.data)
 353                             elif m.name == 'STRT':
 354                                 self.start = float(m.data)
 355                                 self.start_units = m.units
 356                             elif m.name == 'STOP':
 357                                 self.stop = float(m.data)
 358                                 self.stop_units = m.units
 359                             elif m.name == 'STEP':
 360                                 self.step = float(m.data)
 361                                 self.step_units = m.units
 362             line = f.readline()
 363 
 364         # Finished reading the header--all that is left is the numerical
 365         # data that follows the '~A' line.  We'll construct a structured
 366         # data type, and, if the data is not wrapped, use numpy.loadtext
 367         # to read the data into an array.  For wrapped rows, we use the
 368         # function _read_wrapped() defined elsewhere in this module.
 369         # The data type is determined by the items from the '~Curves' section.
 370         dt = np.dtype([(name, float) for name in self.curves.names])
 371         if self.wrap:
 372             a = _read_wrapped_data(f, dt)
 373         else:
 374             a = np.loadtxt(f, dtype=dt)
 375         self.data = a
 376 
 377         if opened_here:
 378             f.close()
 379 
 380 
 381 if __name__ == "__main__":
 382     import sys
 383 
 384     las = LASReader(sys.argv[1], null_subs=np.nan)
 385     print "wrap? ", las.wrap
 386     print "vers? ", las.vers
 387     print "null =", las.null
 388     print "start =", las.start
 389     print "stop  =", las.stop
 390     print "step  =", las.step
 391     print "Version ---"
 392     las.version.display()
 393     print "Well ---"
 394     las.well.display()
 395     print "Curves ---"
 396     las.curves.display()
 397     print "Parameters ---"
 398     las.parameters.display()
 399     print "Other ---"
 400     print las.other
 401     print "Data ---"
 402     print las.data2d

Source code: las.py

Here's an example of the use of this module:

>>> import numpy as np
>>> from las import LASReader
>>> sample3 = LASReader('sample3.las', null_subs=np.nan)
>>> print sample3.null
-999.25
>>> print sample3.start, sample3.stop, sample3.step
910.0 909.5 -0.125
>>> print sample3.well.PROV.data, sample3.well.UWI.data
ALBERTA 100123401234W500
>>> from matplotlib.pyplot import plot, show
>>> plot(sample3.data['DEPT'], sample3.data['PHIE'])
[<matplotlib.lines.Line2D object at 0x4c2ae90>]
>>> show()

It creates the following plot:

sample3plot.png

The sample LAS file is here:

sample3.las

SciPy: Cookbook/LASReader (last edited 2015-10-24 17:48:24 by anonymous)