Use Calculate Field GP tool to update feature coordinates

Question

ArcMap 10.7.1 (Oracle SDE.ST_GEOMETRY):

I want to use the Calculate Field GP tool to update the coordinates of a polyline FC's shape column.

---

I have a script that works in a standalone Python IDE (PyScripter):

Source: Set M-values to cumulative length of line (via ArcPy)

import arcpy
connection = "Database Connections\my_conn.sde"
feature_class = connection + "\my_owner.my_fc"
spatial_reference = arcpy.Describe(feature_class).spatialReference

with arcpy.da.Editor(connection) as edit_session:
    with arcpy.da.UpdateCursor(feature_class, "SHAPE@") as cursor:
        for row in cursor:
            geometry = row[0].densify("ANGLE", 10000, 0.174533)
            parts = arcpy.Array()
            for part in geometry:
                points = arcpy.Array()
                for point in part:
                    point.M = geometry.measureOnLine(point)
                    points.append(point)
                parts.append(points)
            row[0] = arcpy.Polyline(parts, spatial_reference)
            cursor.updateRow(row)

---

I want to use that same logic, but in the Calculate Field GP tool:

Expression:
-----------
new_shape( !SHAPE! )

Code Block:
-----------
def new_shape(geometry):
    geometry = geometry.densify("ANGLE", 10000, 0.174533)
    parts = arcpy.Array()
    for part in geometry:
        points = arcpy.Array()
        for point in part:
            point.M = geometry.measureOnLine(point)
            points.append(point)
        parts.append(points)
    return arcpy.Polyline(parts)

---

When I run that Calculate Field code, I get an error:

"geoprocessing describe geometry object' object is not iterable"

ERROR 000539: Error running expression: new_shape( GPVARIANTOBJECT0 ) 
Traceback (most recent call last):
  File "<expression>", line 1, in <module>
  File "<string>", line 4, in new_shape
TypeError: 'geoprocessing describe geometry object' object is not iterable

Failed to execute (CalculateField).

---

What am I doing wrong?

The geometry was iterable when using the full blown script in PyScripter. I don't understand why it isn't iterable in the Calculate Field script too.

Answer 1

Some years ago ESRI changed the behaviour of field calculator. In the past you could directly manipulate and update the Shape field as you are attempting within a field calculator. This functionality was removed for what ever reason and you can only do it through arcpy.

So basically you can't do what you are trying. If your intention was to run the tool in model builder you could try your first code sample within the model only calculate value tool.

Answer 2

Solved:

For an unknown reason, for part in geometry: and for point in part: don't seem to work in the Field Calculator.

But if we replace for part in geometry: with:

for i in range(geom.partCount):
    part = geom.getPart(i)

And replace for point in part: with:

for j in range(part.count):
    point = part.getObject(j)

...then that works as expected. It loops through the parts and the points without issue.

Expression:
-----------
new_shape( !SHAPE! )

Code Block:
-----------
def new_shape(geom):
    spatial_reference = geom.spatialReference
    geom = geom.densify("ANGLE", 10000, 0.174533)
    parts = arcpy.Array()
    for i in range(geom.partCount):
        part = geom.getPart(i)                      #I think "part" is an ArcPy array: https://desktop.arcgis.com/en/arcmap/10.7/analyze/arcpy-classes/array.htm
        points = arcpy.Array()
        for j in range(part.count):                 #"part.count" is the number of points within the part array.
                                                    #Alternatively, this would work in the Calculate Field GP tool too: "for j in range(len(part)):"   ...and this as well: "for j, point in enumerate(part):"
                                                    #...but "for point in part:"  doesn't work in the Calculate Field GP tool. It runs without errors, but it doesn't actually loop through the points. So the script returns an empty shape, which isn't what we want.

                                                    #Note: this doesn't work in the Calculate Field GP tool: "point = part[j]"  ...although that might work in newer versions of ArcPy. See a newer version of the ArcPy Array docs: "The getObject method is equivalent to indexing an object; that is, obj.getObject(0) is equivalent to obj[0]."
            point = part.getObject(j)               #But this does work: ".getObject(j)"
            point.M = geom.measureOnLine(point)
            points.append(point)
        parts.append(points)
    return arcpy.Polyline(parts, spatial_reference)

Screenshot of using the field calculator on the shape column.

Answer 3

Updated Version:

(In this example, I'm using the Field Calculator tool in the Attribute Table --> Select By Attributes. But it would work in the Calculate Field geoprocessing tool too.)

---

Use the Select By Attribute tool in the Attribute Window to select lines that have blank M-values:

asset_id in
(select
    asset_id
from
    my_owner.my_fc
where
    abs(sde.st_maxm(shape) - sde.st_length(shape))   > 1
)

We need to use a subquery because there’s a bug in ArcMap – the WHERE clause errors-out when we run it on its own: Use SDE.ST_GEOMETRY functions in Select by Attributes window.

---

The Python field calculator script:

Pre-Logic Script Code:
----------------------
def new_shape(geom):
    orig_geom = geom
    spatial_reference = geom.spatialReference
    geom = geom.densify("ANGLE", 10000, 0.174533)
    parts = arcpy.Array()
    for i in range(geom.partCount):
        part = geom.getPart(i)
        points = arcpy.Array()
        for j in range(part.count):
            point = part.getObject(j)
            point.M = geom.measureOnLine(point)
            points.append(point)
        parts.append(points)
    new_geom = arcpy.Polyline(parts, spatial_reference)

    #Check the new geometry against the original geometry. If the new geometry is invalid, return the old geometry instead.
    if new_geom   and   new_geom.pointCount>=orig_geom.pointCount   and   abs(orig_geom.length-new_geom.length) < 1   and   new_geom.partCount==orig_geom.partCount:
        return new_geom
    else:
        return orig_geom

SHAPE = 
----------------------
new_shape( !SHAPE! )