Cropping a set of lasfiles to rectangular bounds

crop.py

#!/usr/bin/env python

import os
import subprocess


PROJECT_DIR = os.path.abspath(os.path.dirname(__file__))

def path(*a):
    return os.path.join(PROJECT_DIR, *a)


LASFILE_DIR = path("UNC20130419LASv1-0")
CROP_DIR = path("cropped")
WKT_FILE = path("output.wkt", "ASO_April19_2013_0.8_octree_boundary.csv")


PIPELINE_XML = """
<?xml version="1.0" encoding="utf-8"?>
<Pipeline version="1.0">
  <Writer type="drivers.las.writer">
    <Option name="filename">
      %(outfile)s
    </Option>
    <Filter type="filters.crop">
      <Option name="bounds">
        ([259171.9260, 266197.1302], [4193597.2623, 4201173.8345], [0, 100000])
      </Option>
      <Reader type="drivers.las.reader">
        <Option name="filename">
          %(infile)s
        </Option>
      </Reader>
    </Filter>
  </Writer>
</Pipeline>
"""


def main():
    try:
        os.mkdir(CROP_DIR)
    except OSError:
        pass
    with open(WKT_FILE) as f:
        wkt = f.read()
    for infile in os.listdir(LASFILE_DIR):
        xml = PIPELINE_XML % {
                "infile": os.path.join(LASFILE_DIR, infile),
                "outfile": os.path.join(CROP_DIR, infile),
                "wkt": wkt,
                }
        print "Cropping ", os.path.join(LASFILE_DIR, infile)
        p = subprocess.Popen(["pcpipeline", "--stdin", "-v", "5"], stdin=subprocess.PIPE,
                stdout=subprocess.PIPE, stderr=subprocess.PIPE)
        (stdout, stderr) = p.communicate(xml)
        print stdout
        print stderr


if __name__ == "__main__":
    main()