Construct Points on polyline at 3d distance

JasonBarnes

Greetings. I'm having trouble finding information on the construct points tool. I have about 100 separate hiking trails that have z geometry. I have calculated points at 1 mile intervals along the track using the Construct Points tool in the editor toolbar. The problem I'm having is determining whether these points are based on the surface distance or euclidean distance along the line. I can't seem to find any documentation on this in the help files or in the forum. Any insight? I can upload files if needed.

Thanks.

XanderBakker

It would be nice if there there were additional methods for 3D measurements. All methods are 2D on the polyline geometry. There is the segmentAlongLine method that mentions measures, so I was hopeful that converting the geometry to Polyline ZM would allow for extracting a sub-segment and use the lastPoint of the polyline to get the point. In a previous thread I published a snippet to convert a PolylineZ to PolylineZM and assign the 3D length as M value:

https://community.esri.com/message/516789#comment-516789

That turned out to be wrong. So I used the snippet and shaped it a little bit to create a positionAlongLine3D function.

def main():
    import arcpy
    arcpy.env.overwriteOutput = True
    sr = arcpy.SpatialReference(28992)
    lst_coords = [[[0, 0, 0], [1000, 0, 0], [1000, 1000, 0], [0, 1000, 0], [0, 0, 0]],
                  [[2000, 0, 0], [3000, 0, 1000], [3000, 1000, 2000], [2000, 1000, 3000], [2000, 0, 4000]]]
    lines = []
    points = []
    lines_m = []
    points_m = []
    points_ok = []
    for feature in lst_coords:
        polyline = arcpy.Polyline(
            arcpy.Array([arcpy.Point(*coords) for coords in feature]), sr, True, False)
        lines.append(polyline)
        segment = polyline.segmentAlongLine(0, 2500, False)
        points.append(arcpy.PointGeometry(segment.lastPoint, sr, True, False))
        polyline_zm = CreateZAwareFeature(polyline, sr)
        lines_m.append(polyline_zm)
        segment_m = polyline_zm.segmentAlongLine(0, 2500, False)
        points_m.append(arcpy.PointGeometry(segment_m.lastPoint, sr, True, True))
        pnt_ok = positionAlongLine3D(polyline, 2500)
        points_ok.append(pnt_ok)

    arcpy.CopyFeatures_management(lines, r'C:\GeoNet\Query3DpointAndDistance\test.gdb\lines')
    arcpy.CopyFeatures_management(points, r'C:\GeoNet\Query3DpointAndDistance\test.gdb\points')
    arcpy.CopyFeatures_management(lines_m, r'C:\GeoNet\Query3DpointAndDistance\test.gdb\lines_m')
    arcpy.CopyFeatures_management(points_m, r'C:\GeoNet\Query3DpointAndDistance\test.gdb\points_m')
    arcpy.CopyFeatures_management(points_ok, r'C:\GeoNet\Query3DpointAndDistance\test.gdb\points_ok')

def positionAlongLine3D(polyline_z, dist):
    length = 0
    if dist <= 0:
        return arcpy.PointGeometry(polyline_z.firstPoint, polyline_z.spatialReference, True)
    elif dist >= polyline_z.length:
        return arcpy.PointGeometry(polyline_z.lastPoint, polyline_z.spatialReference, True)
    else:
        for part in polyline_z:
            cnt = 0
            for pnt in part:
                if cnt == 0:
                    pnt_prev = arcpy.Point(pnt.X, pnt.Y, pnt.Z)
                length += get3DLength(pnt, pnt_prev)
                if length > dist:
                    d1 = dist - prev_length
                    seg_length = length - prev_length
                    frac = d1 / seg_length
                    segment = arcpy.Polyline(arcpy.Array([pnt_prev, pnt]), polyline_z.spatialReference)
                    pnt_3d = segment.positionAlongLine(frac, True)
                    z = (1-frac) * pnt_prev.Z + frac * pnt.Z
                    pntz = arcpy.Point(pnt_3d.firstPoint.X, pnt_3d.firstPoint.Y, z)
                    return arcpy.PointGeometry(pntz, polyline_z.spatialReference, True)
                    break
                cnt += 1
                pnt_prev = pnt
                prev_length = length

def CreateZAwareFeature(polyline_z, sr):
    lst_poly_zm = []
    length = 0
    for part in polyline_z:
        part_zm = []
        cnt = 0
        for pnt in part:
            if cnt == 0:
                pnt_prev = arcpy.Point(pnt.X, pnt.Y, pnt.Z)
            length += get3DLength(pnt, pnt_prev)
            pnt_zm = arcpy.Point(pnt.X, pnt.Y, pnt.Z, length)
            part_zm.append(pnt_zm)
            pnt_prev = arcpy.Point(pnt.X, pnt.Y, pnt.Z)
            cnt += 1
        lst_poly_zm.append(part_zm)
    polyline_zm = arcpy.Polyline(arcpy.Array(lst_poly_zm), sr, True, True)
    return polyline_zm

def get3DLength(pnt1, pnt2):
    return math.sqrt((pnt1.X - pnt2.X)**2 + (pnt1.Y - pnt2.Y)**2 + (pnt1.Z - pnt2.Z)**2)

if __name__ == '__main__':
    main()

So the script creates two lines, the first is 2D (no Z) and the second has the same shape but the sides have an angle of 45° in Z direction, which mean that they are root(2) longer as shown in the field Length3D:

If you look at the lines in 2D they look the same, but the (3D) length isn't:

When you want to get the point at a distance of 2500m from the start in 3D the result should be:

... or looking at it in 3D:

The alternative would be to use dynamic segmentation on the M-aware polylines.