PostGIS Commands in Arcpy

Most of you know, I love spatial SQL.  I also love PostGIS and ArcGIS – can’t we just all get along?  Well, yes we can.  In this case, Python becomes the great mediator.

Check out this video where I analyze data from an ESRI shapefile using Arcpy, call PostGIS to run a procedure on the data using the psycopg2 package, and return the results.

In this case, I am trying to compute the distance between points and line in meters, but the data is in Lat/Lon.  So, I am bringing PostGIS into the mix because it does a very good job of returning distances between points and lines.


Now, you can do all of this in Arcpy directly by getting the geometry object:

import arcpy, numpy

pts = arcpy.MakeFeatureLayer_management("c:/temp/bkpts.shp")
selpts = arcpy.SelectLayerByAttribute_management(pts,"NEW_SELECTION",'POP100 > 10')

rivers = arcpy.MakeFeatureLayer_management("c:/temp/rivers.shp")

for pt in arcpy.da.SearchCursor(selpts, ["SHAPE@"]):
 ptutm = pt[0].projectAs(arcpy.SpatialReference(26918))
 for riv in arcpy.da.SearchCursor(rivers,["SHAPE@"]):
 rivutm = riv[0].projectAs(arcpy.SpatialReference(26918))
 dist = rivutm.distanceTo(ptutm)
 print dist

but, my point here was to show that you can move along with an Arcpy script, take a short and painless detour to use PostGIS, and then get back to your Arcpy script (assuming that PostGIS does something that you are really interested in using).

want to learn more about how to program with Python, Arcpy, and PostGIS?  Check out my courses on on

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s