4

Below is a code snippet that I have extracted from something that I am working on. It should create two polyline features: the first is expected to be tiny, while the second is expected to be small.

The small polyline is created as expected, but the tiny line is created as a feature of zero length when I think it should be about 4 metres long.

I have experimented with setFalseOriginAndUnits on the SpatialReference object to see if I can get a greater than zero length but have been unsuccessful. I thought that this was worth trying on the basis of my learnings at Controlling XY Resolution on file geodatabase feature class output from Append?

How can I get this tiny feature to be created with a tiny but not zero length?

enter image description here

import arcpy

if arcpy.Exists("C:/Temp/test.gdb"):
    arcpy.Delete_management("C:/Temp/test.gdb")

arcpy.CreateFileGDB_management("C:/Temp","test")
sr = arcpy.SpatialReference("Geocentric Datum of Australia 1994")
sr.setFalseOriginAndUnits(-400.0,-400.0,10000000.0)
arcpy.CreateFeatureDataset_management("C:/Temp/test.gdb","testFD",sr)
arcpy.CreateFeatureclass_management("C:/Temp/test.gdb/testFD","testLines","POLYLINE")
cursor = arcpy.da.InsertCursor("C:/Temp/test.gdb/testFD/testLines",["SHAPE@"])
# Tiny line but NOT zero length
point = arcpy.Point()
array = arcpy.Array()
point.X = 153.06436725
point.Y = -27.49210825
array.add(point)
point.X = 153.06437825
point.Y = -27.49208025
array.add(point)
polyline = arcpy.Polyline(array)
cursor.insertRow([polyline])
# Small line
point = arcpy.Point()
array = arcpy.Array()
point.X = 153.0640
point.Y = -27.4920
array.add(point)
point.X = 153.0650
point.Y = -27.4930
array.add(point)
polyline = arcpy.Polyline(array)
cursor.insertRow([polyline])
del cursor

with arcpy.da.SearchCursor("C:/Temp/test.gdb/testFD/testLines",["Shape_Length"]) as cursor:
    for row in cursor:
        print "Length of line is {0}".format(row[0])

The output that results from the above code is:

>>> 
Length of line is 0.0
Length of line is 0.00141421356242
>>> 
  • You don't apply the sr to the Polyline constructor (second parameter), leaving the default spatial reference with a (relatively) large scale units. – Vince Nov 10 '16 at 1:42
  • @Vince, should that matter? The fc will inherit the sr of the feature dataset. – Fezter Nov 10 '16 at 1:50
  • Thanks @Vince - that was it - feel free to write up a quick answer or I can do one later. – PolyGeo Nov 10 '16 at 1:51
7

In the ArcGIS world view, there is a difference between a Point and a PointGeometry feature -- points are double precision objects with 2, 3, or 4 dimensions, and the feature is an object with a SpatialReference, with location encoded using the false origins and units for the dimensions of that spatial reference. Polyline and Polygon features are constructed from an array of points, but also encoded with the spatial reference properties.

It is perhaps unfortunate that the PointGeometry, Polyline, and Polygon constructors make the spatial reference parameter optional, because without an explicit spatial reference, a Miranda one is applied. Yes, once the feature goes into the feature class, the spatial reference is transformed to that of the feature class, but that means extra work, and possible precision loss and/or vertex collapse if your data isn't in meters.

The solution is to always explicitly apply the spatial reference associated with the feature class at geometry construction.

I've encountered situations using the 'C' API where the coordinate reference used to create the feature class is subtly (and sometimes not so subtly) modified by the layer creation process, so in Python I try to fetch the actual (as-built) SpatialReference from a Describe, and apply that to the geometry:

sr = arcpy.SpatialReference("Geocentric Datum of Australia 1994")
sr.setFalseOriginAndUnits(-400.0,-400.0,10000000.0)
arcpy.CreateFeatureDataset_management("C:/Temp/test.gdb","testFD",sr)
arcpy.CreateFeatureclass_management("C:/Temp/test.gdb/testFD","testLines","POLYLINE")

sr = arcpy.Describe("C:/Temp/test.gdb/testFD").spatialReference
cursor = arcpy.da.InsertCursor("C:/Temp/test.gdb/testFD/testLines",["SHAPE@"])
point = arcpy.Point()
array = arcpy.Array()
point.X = 153.06436725
point.Y = -27.49210825
array.add(point)
point.X = 153.06437825
point.Y = -27.49208025
array.add(point)
polyline = arcpy.Polyline(array,sr)  // Don't forget the second parameter!
// ...

Be careful to only use feature datasets when required by geodatabase functionality constraints, since it adds overhead to nominal operation.

| improve this answer | |

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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