I am parsing an OpenStreetMap gpx file which I exported from here. I'm trying to put the data into arcgis geodatabase feature classes. Some of the features in the gpx are points (saved as wpt), while others are polygons (saved as trk consisting of trkpt's). Using the script below, I am succesful in copying over the points, their geometry, the names and attributes I want. However somehow creating the geometry does not work for the polygons, although the attributes are copied over without trouble (I left the bit of code for the attributes out of the example below). I found this question which prompted me to try to build the polygons with the points ordered the other way round. However this didn't work either. I built in quite some print statements to check where it goes wrong, see printout below. Also not that printing poly.area
results in 0.0
, printing poly.trueCentroid
results in None
, etc. So the polygon simply isn't generated from the array. To me it appears that poly = arcpy.Polygon(ar)
simply does not work.
GPX example:
<wpt lat="55.7962684" lon="-4.6260225">
<name>Ian Brown Motors</name>
<desc>name=Ian Brown Motors shop=car_repair</desc>
<link href="http://osm.org/browse/node/3452158061" />
</wpt>
<trk>
<name>Key Store Brodie's Minimarket</name>
<desc>building=yes name=Key Store Brodie's Minimarket opening_hours=08:00-18:00 shop=supermarket source=survey</desc>
<link href="http://osm.org/browse/way/28476052" />
<trkseg>
<trkpt lat="56.4387108" lon="-4.7134965" />
<trkpt lat="56.4385955" lon="-4.7135666" />
<trkpt lat="56.438564" lon="-4.7133971" />
<trkpt lat="56.4386793" lon="-4.7133269" />
<trkpt lat="56.4387108" lon="-4.7134965" />
</trkseg>
</trk>
Script:
from xml.etree import ElementTree as Et
import os
import arcpy
arcpy.env.overwriteOutput = 1
wgs = arcpy.SpatialReference(r'C:\Program Files\ArcGIS\Desktop10.0\Coordinate '
r'Systems\Geographic Coordinate Systems\World\WGS '
r'1984.prj')
gdb = PATH_TO_GDB
f = PATH_TO_GPX
all_records = []
p = Et.parse(f)
root = p.getroot()
arcpy.CreateFeatureclass_management(gdb, 'shops_pt', 'POINT',
spatial_reference=wgs)
arcpy.CreateFeatureclass_management(gdb, 'shops_pl', 'POLYGON',
spatial_reference=wgs)
shops_pt = os.path.join(gdb, 'shops_pt')
shops_pl = os.path.join(gdb, 'shops_pl')
for fc in [shops_pt, shops_pl]:
arcpy.AddField_management(fc, 'name', 'TEXT', 255)
arcpy.AddField_management(fc, 'type', 'TEXT', 50)
arcpy.AddField_management(fc, 'opening_hours', 'TEXT', 255)
ic_pt = arcpy.InsertCursor(shops_pt)
ic_pl = arcpy.InsertCursor(shops_pl)
for mainElement in root:
if mainElement.tag == 'metadata':
continue
if mainElement.tag == 'wpt':
r = ic_pt.newRow()
y = float(mainElement.attrib['lat'])
x = float(mainElement.attrib['lon'])
p = arcpy.Point(x, y)
r.shape = p
ic_pt.insertRow(r)
elif mainElement.tag == 'trk':
r = ic_pl.newRow()
for subElement in mainElement:
if subElement.tag == 'trkseg':
ar = arcpy.Array()
for trkpt in subElement:
y = float(trkpt.attrib['lat'])
x = float(trkpt.attrib['lon'])
p = arcpy.Point(x, y)
print p
ar.add(p)
print 'array length', len(ar)
print 'ar[0]', ar.getObject(0), 'ar[-1]', ar.getObject(len(ar)-1)
poly = arcpy.Polygon(ar)
r.shape = poly
print 'poly.centroid try 1', poly.centroid
if poly.centroid is None:
ar2 = arcpy.Array()
for p in ar:
ar2.insert(0, p)
for p in ar2:
print p
poly = arcpy.Polygon(ar2)
print 'poly.centroid try 2', poly.centroid
ic_pl.insertRow(r)
Printout sample for 1 polygon:
-4.7134965 56.4387108 NaN NaN
-4.7135666 56.4385955 NaN NaN
-4.7133971 56.438564 NaN NaN
-4.7133269 56.4386793 NaN NaN
-4.7134965 56.4387108 NaN NaN
array length 5
ar[0] -4.7134965 56.4387108 NaN NaN ar[-1] -4.7134965 56.4387108 NaN NaN
poly.centroid try 1 None
-4.7134965 56.4387108 NaN NaN
-4.7133269 56.4386793 NaN NaN
-4.7133971 56.438564 NaN NaN
-4.7135666 56.4385955 NaN NaN
-4.7134965 56.4387108 NaN NaN
poly.centroid try 2 None
Solution: as suggested by @Vince below in the comments, when putting the Array
in a Polygon
, explicitly declare the spatial reference, so it would become poly = arcpy.Polygon(ar, wgs)
.
poly = arcpy.Polygon(ar,wgs)
wgs = arcpy.SpatialReference(4326)