I am using ArcGIS Pro 3.3 and ArcPy. I have a point feature class with the field "Elevation" containing elevation in meters. I am trying to determine if, for each point, it is the highest within a radius of 500 m using buffers and a spatial join (that is, the highest elevation of a point within 500 m of the point being examined). The program completes OK, but simply gives back the elevation of the point being examined, rather than that the maximum value. It must be my field mapping parameters, but having fiddled with various combinations, I've had no luck (I either get the elevation of the point being examined or the elevation of the first point within the buffer). I should add that the spatial join does get the correct number of points in each buffer.
import arcpy
from arcpy import env
# Outputs can be overwritten
env.overwriteOutput = True
# Set the workspace environment to the geodatabase containing the feature classes
arcpy.env.workspace = r"C:\Users\WMR\Documents\ArcGIS\scratch.gdb"
# Set the name of the input point layer
point_layer = "XBC_Cerros_2024"
# Set the name of the output point layer with maximum elevation
output_layer = "XBC_with_max_elevation3"
# Set the buffer distance
buffer_distance = "500 Meters"
# Set the name of the elevation field
elevation_field = "Elevation"
# Create a buffer around each point
buffered_points = arcpy.Buffer_analysis(point_layer, "in_memory\\Buffered_Points",
buffer_distance)
# Perform spatial join to find maximum elevation within 500 meters for each point
field_mappings = arcpy.FieldMappings()
field_mappings.addTable(buffered_points)
field_mappings.addTable(point_layer)
# Map the elevation field from the input point layer to the output layer with MAX statistic
elevation_map = arcpy.FieldMap()
elevation_map.addInputField(point_layer, elevation_field)
elevation_output = elevation_map.outputField
elevation_output.name = "MaxElev"
elevation_output.aliasName = "MaxElev"
elevation_output.statisticsType = "Max"
elevation_map.outputField = elevation_output
field_mappings.addFieldMap(elevation_map)
# Perform the spatial join using the field mappings
arcpy.analysis.SpatialJoin(buffered_points, point_layer, output_layer,
"JOIN_ONE_TO_ONE", "KEEP_ALL", field_mappings,"INTERSECT")
# Clean up temporary buffer layer
arcpy.management.Delete(buffered_points)
print("Spatial join completed.")
