Using ArcObjects for coordinate transformations in Python

Some of you may find yourselves in situations where you may need to develop in Python but are limited by the existing geoprocessing methods available.  Luckily, there is a method of accessing the ArcObject COM libraries through Python.  The comtypes package allows you use all of the methods available to .NET developers.  The presentation Using ArcObjects in Python gives an excellent overview of how to accomplish this.

I’m going to walk you through a specific example using comtypes.  In this example, we will perform a coordinate conversion on a point from a Geographic Coordinate System to a Projected Coordinate System.  To get started, download and install comtypes.  Once this is complete it’s time to start coding.  We’ll start off by defining four functions.  These will make more sense later.  First get the COM library path:

def GetLibPath():

   import _winreg

   keyESRI = 
      _winreg.OpenKey(_winreg.HKEY_LOCAL_MACHINE,

                        "SOFTWARE\\ESRI\\ArcGIS")

   return _winreg.QueryValueEx(keyESRI,  
          "InstallDir")[0] + "com\\"

Next we need to create a procedure to get the ArcObjects module:

def GetModule(sModuleName):

   import comtypes

   import comtypes.client as comc

   sLibPath = GetLibPath()

   comc.GetModule(sLibPath + sModuleName)

Create a function that creates a new object:

def NewObj(MyClass, MyInterface):

   import comtypes

   try:

        ptr = comtypes.client.CreateObject(MyClass,
               interface=MyInterface)

        return ptr

   except:

        return None

Finally, create a function for casting objects:

def CType(obj, interface):

   try:

        newobj = obj.QueryInterface(interface)

        return newobj

   except:

        return None

Now we need to define the main function to perform the coordinate transformation:

def TransformCoords(x, y, gcs, pcs, trans, direction):

The x and y are the longitude and latitude of the point in decimal degrees.  In order to perform transformation you will need to look up some ESRI codes.  The list of geographic coordinate system codes are here: Page 1, Page 2, Page 3.  The projected coordinate systems are here: Page 1, Page 2, Page 3, Page 4.  And the transformations are here: Page 1, Page 2, Page 3.   Lastly, you’ll need the direction codes here.

For those familiar with writing ArcObjects in C# I will include the C# equivalent with each piece of code.  Now we need to get the Geometry module and import it:

   GetModule("esriGeometry.olb") 

   import comtypes.gen.esriGeometry as esriGeometry

C#:

using ESRI.ArcGIS.Geometry;

Next we create a new point object and assign the x and y values:

   point = NewObj(esriGeometry.Point,
          esriGeometry.IPoint)

   point.X = x

   point.Y = y

C#:

Point point = new Point();

point.X = x;

point.Y = y;

It’s time to start using the codes we looked up.  Let’s define the spatial reference for the geographic coordinate system:

   sRefEnv = 

       NewObj(esriGeometry.SpatialReferenceEnvironment,

          esriGeometry.ISpatialReferenceFactory)

   incomingCoordSystem =

   CType(sRefEnv.CreateGeographicCoordinateSystem(gcs),

      esriGeometry.ISpatialReference)
   point.SpatialReference = incomingCoordSystem

C#:
SpatialReferenceEnvironment sRefEnv = new  
   SpatialReferenceEnvironment();

ISpatialReference incomingCoordSystem = 

   sRefEnv.CreateGeographicCoordinateSystem(gcs);

point.SpatialReference = incomingCoordSystem;

Define the projected coordinate system:

   outgoingCoordSystem = 

   CType(sRefEnv.CreateProjectedCoordinateSystem(pcs), 

      esriGeometry.ISpatialReference)

C#:

ISpatialReference outgoingCoordSystem = 
sRefEnv.CreateProjectedCoordinateSystem(pcs);

Define the geographic transformation:

   geoTransformation = 

      CType(sRefEnv.CreateGeoTransformation(trans), 

         esriGeometry.IGeoTransformation)

C#:

IGeoTransformation geoTransformation = 
(IGeoTransformation)sRefEnv.CreateGeoTransformation(trans);

Re-project the coordinates and assign the new geometry to the point:

   geometry = CType(point, esriGeometry.IGeometry2)

   geometry.ProjectEx(outgoingCoordSystem, 
          direction,geoTransformation, False, 0, 0)

   point = CType(geometry, esriGeometry.IPoint)

       

C#:

IGeometry2 geometry = point as IGeometry2;

geometry.ProjectEx(outgoingCoordSystem, esriTransformDirection.esriTransformForward, geoTransformation, 
   false, 0, 0);

point = (Point)geometry;

Finally, assign the coordinates to an array and return them:

   points = [point.X, point.Y]

       

   return points

In the end you have an array of the transformed x and y of the original point.  Another method exists that can also be helpful called ctypes.  With ctypes, you can access C libraries much like comtypes.  An excellent example of transformation using ctypes can be found here.