10

I'm serializing my arcpy geometries as geojson so that I can 'hydrate' them back as geometries later and I'm having 2 problems in the cycle.:

PROBLEM 1: Precision

    R0 = arcpy.SearchCursor(self.shpTest, "FID=0").next().getValue("Shape")          
    geojson = R0.__geo_interface__                        
    R1 = arcpy.AsShape(geojson)
    self.assertTrue(R0.equals(R1)) <<< THIS FAILS

If I check the string representation, the coordinates have slightly changed:

    geojson2 = R1.__geo_interface__
    print geojson
    print geojson2  

    {'type': 'Polygon', 'coordinates': [[(442343.5516410945, 4814166.6184399202), (442772.17749834526, 4811610.7383281607), (441565.67508534156, 4811499.6131059099), (440772.50052100699, 4814184.7808806188), (442343.5516410945, 4814166.6184399202)]]}
    {'type': 'Polygon', 'coordinates': [[(442343.55169677734, 4814166.6185302734), (442772.17749023438, 4811610.73828125), (441565.67510986328, 4811499.6130981445), (440772.50048828125, 4814184.7808837891), (442343.55169677734, 4814166.6185302734)]]}

PROBLEM 2: Holes If the polygon has holes, geo_interface generates an error:

    R0_WithHoles = arcpy.SearchCursor(self.shpTest, "FID=0").next().getValue("Shape")          
    geojson = R0.__geo_interface__  <<< generates this ERROR:

    File "C:\Program Files\ArcGIS\Desktop10.0\arcpy\arcpy\arcobjects\geometries.py", line 68, in __geo_interface__
        return {'type': 'Polygon', 'coordinates': [[(pt.X, pt.Y) for pt in part] for part in self]}
    AttributeError: 'NoneType' object has no attribute 'X'

Any ideas on how to solve these problems?

7
  • Yep, just ran across number 2 myself. And doesn't appear to be much love for this topic. Commented Dec 8, 2011 at 0:43
  • This is still broken in arcpy in ArcGIS 10.1 -- It would be nice if ESRI could comment on the subject. Commented Sep 17, 2012 at 23:42
  • I came across the first and second problems. With me, the coordinated do not seem to change (when you print them) but geom1.equals(geom2) fails me only a few times. I am not sure why that happens as well. The 2nd problem was fixed using @valveLondon 's suggestion. If you found out how to fix the .equals please do share. Commented Nov 7, 2012 at 17:46
  • @MichalisAvraam We had the same issue also and got onto ESRI for a solution - turns out it's a known bug (when you create a geom without a projection it truncates the precision) - take a look at this question also.
    – om_henners
    Commented Nov 12, 2012 at 21:12
  • @om_henners I assumed that. But the arcpy.AsShape() function does not let you specify a spatial reference. I have set all environment variables hoping it would do something (output coords, etc...). The solution then is to manually decode the GeoJSON because ESRI doesn't care about accuracy? Commented Nov 12, 2012 at 23:40

1 Answer 1

5

OK - well I thought I had solved it.

replace line ~80 of this file C:\Python26\ArcGIS10.0\Lib\arcpy\arcobjects\geometries.py from this:

return {'type': 'Polygon', 'coordinates': [[(pt.X, pt.Y) for pt in part] for part in self]}

to this (or something that is more concise and elegant and does the same thing):

  obj = {"type": "Polygon"}
    coordinates = []
    for part in self:
        _part = []
        for pt in part:
            if pt is not None:
                print pt
                _part.append([pt.X,pt.Y])
            else:
                print "none"
                coordinates.append(_part)
                _part=[]
        coordinates.append(_part)
    obj["coordinates"]=coordinates
    return obj

Basically they forgot to consider donuts in the shape which are marked by null point values. This spits out good geoJson (separate parts) but the arcpy.AsShape method trashes GeoJSON.

this code:

import arcpy
gj = {
  'type': 'Polygon', 'coordinates': [
   [[-122.803764, 45.509158], [-122.796246, 45.500050], [-122.808193, 45.500109],
      [-122.803764, 45.509158]],
   [[-122.804206, 45.504509], [-122.802882, 45.502522], [-122.801866, 45.504479], 
      [-122.804206, 45.504509]]
   ]
 }

 p = arcpy.AsShape(gj)
 print p.__geo_interface__

outputs this:

    {'type': 'Polygon', 'coordinates': [[[-122.8037109375, 45.50927734375],  
    [-122.79620361328125, 45.5001220703125], [-122.80810546875, 45.5001220703125],
    [-122.8037109375, 45.50927734375]]]}

I give up. ;)

Update The holes problem has been solved at 10.1 with this chunk of python:

return {'type': 'Polygon', 'coordinates': [[((pt.X, pt.Y) if pt else None)
                                                    for pt in part]
                                                        for part in self]}
11
  • shouldn't that return a dictionary instead of a string representing a dictionary? :)
    – blah238
    Commented Mar 13, 2012 at 1:14
  • Yes, you're right, it should. I changed it to spit out a valid GeoJSON dictionary obj. but after checking against the AsShape method I realized the futility of my efforts. Commented Mar 14, 2012 at 22:58
  • I wonder if that has anything to do with the problem described in this thread: forums.arcgis.com/threads/9763-Errors-in-arcpy-s-Polygon-class -- should have been fixed in 10 SP2 and definitely 10.1.
    – blah238
    Commented Mar 14, 2012 at 23:36
  • 2
    ESRI updated C:\Program Files\ArcGIS\Server\arcpy\arcpy\arcobjects\geometries.py at 10.1 but if you're at 10.0 you can fix it yourself. Commented Jun 25, 2012 at 21:35
  • 3
    Yes, I fixed it in 10.1, the update above is the new source in the .py file. I thought it made it into a service pack for 10 but I guess not. Commented Jun 25, 2012 at 22:03

Not the answer you're looking for? Browse other questions tagged or ask your own question.