I am writing scripts that take in a layers geometry at a specific location and create a new layer based on it. In order to construct the new geometry I need access to the coordinates of the existing geometry, but cant seem to find an easy way to get it in arcpy. what I have now is some mixture of SelectLayerByLocation to select the objects at a location and then operating on selected objects to get a list of points. Is there a better way to get either a list of points from objects that are at a given location, or geometry objects that represent those objects?
I am starting with X and Y coordinates which I use to create a pointGeometry object.
I then use this through arcpy.management.SelectLayerByLocation() to select anything in the targeted layer that exists at the same location as the point.
from there its fairly easy to either get the points of selected objects using arcpy.da.FeatureClassToNumPyArray() with explode_to_points or the geometry objects by getting a list of selected FIDs and using a SearchCursor to get the geometry objects.
Ideally I would be able to get the same without the selecting or deselecting anything and without any strange middle-man steps
A spatial join is the only other thing that comes to mind but I would prefer SelectLayerByLocation(). There's a code sample at the bottom of the Geometry doc page that shows using arcpy.Geometry() as the output for CopyFeatures() and having immediate access to the geometry for everything. You could try that after they're selected instead of the numpy array or search cursor.
If you don't want to deal with selections, you can check if the point intersects the geometry objects of the targeted layer. Depending on the object count, this might be a lot slower.
# construct your query geometry query_geometry = ... # optional: buffer query_geometry = query_geometry.buffer(100) # read the geometries of the target layer layer_geometries = [r for r in arcpy.da.SearchCursor(layer, ["SHAPE@"])] # filter by intersection with query_geometry intersecting_geometries = [g for g in layer_geometries if not g.disjoint(query_geometry)]