Hobu in ArcUser

The April-June 2005 issue of ArcUser contains an article by Hobu on pages 34-36. It's an introduction to Python and GIS, and vector, grid, and raster processing with ogr.py and gdal.py for users of ArcGIS 9 -- which includes Python and makes a small subset of its functionality scriptable. Good writing, too bad he seems to have retired from blogging ;)

I ventured further into the magaine and found myself comparing the agility of ArcGIS and the GDAL/OGR Python modules while reading the tutorial on pages 46-49 ("Creating Accurate Footprints of Registered Aerial Images"). It's a rather involved and round-about manual process requiring an expensive ArcGIS add-on (Spatial Analyst). I'm sure that one could also manage to open a beer bottle using dentistry tools, but not as readily as with a bottle opener. Consider the following Python script:

# Create footprint shapefile
import ogr

driver = ogr.GetDriverByName('ESRI Shapefile')
footprints_shp = driver.CreateDataSource('footprints.shp')
footprints = footprints_shp.CreateLayer('footprints',
                 geom_type=ogr.wkbPolygon)
fd = ogr.FieldDefn('FILENAME', ogr.OFTString)
fd.SetWidth(30)
footprints.CreateField(fd)

# Loop over a number of georeferenced images
import glob
import gdal

files = glob.glob('*.jpg')
for file in files:
    # Get georeferencing and size of imagery
    dataset = gdal.Open(file)
    g = dataset.GetGeoTransform()
    pixels = dataset.RasterXSize
    lines = dataset.RasterYSize
    minx = g[0]
    maxx = minx + pixels * g[1]
    maxy = g[3]
    miny = maxy + lines * g[5]

    # append to the 'footprints' layer
    wkt = 'POLYGON ((%f %f, %f %f, %f %f, %f %f, %f %f))' \
        % (minx, miny, minx, maxy, maxx, maxy, maxx, miny, minx, miny)
    g = ogr.CreateGeometryFromWkt(wkt)
    f = ogr.Feature(feature_def=footprints.GetLayerDefn())
    f.SetField(0, file)
    f.SetGeometryDirectly(g)
    footprints.CreateFeature(f)
    f.Destroy()

# destroy footprints_shp to flush and close
footprints_shp.Destroy()

No fooling about -- directly to the imagery bounds and write a polygon. Simple and straight-forward, the right tool for the job.