Made ogrdifference 200% faster by testing bounding boxes.

ogrkit/utilities/difference.py

@@ -2,9 +2,13 @@
 
 import shutil
 
+from osgeo import gdal
 from osgeo import ogr
 
 from ogrkit.cli import OGRKitUtility
+from ogrkit.utils import get_bounding_box
+
+gdal.UseExceptions()
 
 class OGRDifference(OGRKitUtility):
     description = 'Produce a shapefile by subtracting a set of shapefiles from the input.'
@@ -29,21 +33,28 @@ def main(self):
             dest_layer.CreateField(source_layer.GetLayerDefn().GetFieldDefn(i))
 
         mask_features = []
+        mask_boxes = []
 
         for mask in self.args.masks:
             geo = ogr.Open(mask)
             layer = geo.GetLayer(0)
 
             for feature in layer:
                 mask_features.append(feature)
+                mask_boxes.append(get_bounding_box(feature.GetGeometryRef()))
 
         for feature in source_layer:
             masked_feature = ogr.Feature(feature_def=source_layer.GetLayerDefn())
             masked_feature.SetFrom(feature)
 
             masked_geometry = feature.GetGeometryRef().Clone()
 
-            for mask_feature in mask_features:
+            for (i, mask_feature) in enumerate(mask_features):
+                bounding_box = mask_boxes[i]
+
+                if not masked_geometry.Intersects(bounding_box):
+                    continue
+
                 masked_geometry = masked_geometry.Difference(mask_feature.GetGeometryRef())
 
             masked_feature.SetGeometryDirectly(masked_geometry)

ogrkit/utils.py

@@ -0,0 +1,20 @@
+#!/usr/bin/env python
+
+from osgeo import gdal
+from osgeo import ogr
+
+gdal.UseExceptions()
+
+def get_bounding_box(geom):
+    """
+    Gets the Envelope of an OGRGeometry as a Polygon.
+    """
+    min_x, min_y, max_x, max_y = geom.GetEnvelope()
+
+    wkt = 'POLYGON((%s %s,%s %s,%s %s,%s %s,%s %s))' % \
+        (min_x, min_y, min_x, max_y,
+        max_x, max_y, max_x, min_y,
+        min_x, min_y)
+
+    return ogr.CreateGeometryFromWkt(wkt)
+