Iterate selection and calculate field not working

3316
44
09-22-2018 12:04 PM
ThomasColson
MVP Frequent Contributor

With the following, I am trying to select each polygon by its name, select points from another feature class that are within the selected polygon, then calculate a field using the same name that was used to select the polygon (by attributes). A spatial intersect, or any other method that results in creating a new FC as output and having to do joins won't work here, as this is (hopefully) to be an unattended nightly update script. It iterates through the polygons fine, but the result is that ALL points get the name of the last polygon. 

import os
import errno
import arcpy
import sys 
import traceback 
from arcpy import env
update_feature_class = r'C:\Users\tcolson\Documents\ArcGIS\Projects\MyProject\MyProject.gdb\GRSM_RESEARCH_LOCATIONS'
county_feature_class = r'C:\Users\tcolson\Documents\ArcGIS\Projects\MyProject\MyProject.gdb\CountyBoundaries'

with arcpy.da.SearchCursor(county_feature_class,['SHAPE@','NAME']) as cursor:
    for row in cursor:
        print('Selecting '+row[1])
        expression = "NAME = '"+row[1]+"'"
        arcpy.SelectLayerByAttribute_management (county_feature_class, "NEW_SELECTION", expression)
        arcpy.SelectLayerByLocation_management(update_feature_class, "INTERSECT", county_feature_class)
        arcpy.CalculateField_management(update_feature_class, "COUNTY","'"+row[1]+"'", "PYTHON3")
44 Replies
DanPatterson_Retired
MVP Esteemed Contributor

Does it work? You indentation is off, or that is a copy-paste-format problem.

or are you looking for alternatives

If you don't want to use the existing arctoolbox tools, you could do it with the polygon geometry class's

'contains' method

Polygon—ArcPy classes | ArcGIS Desktop 

it returns a boolean, so if the feature is within the bounds, you can assign the ID to the feature, then move on to the next polygon

Indicates if the base geometry contains the comparison geometry.  contains is the opposite of within.

Reply
0 Kudos
ThomasColson
MVP Frequent Contributor

No, it doesn't work, it calculates all points, regardless of selection, to the the last "Name" in the loop. How would 

contains (second_geometry, {relation}) work any different? You still have to iterate through a select-by-attribute loop?

Reply
0 Kudos
DanPatterson_Retired
MVP Esteemed Contributor

It is not a select by attribute loop.  You have your points and your polygons.  You loop through the points and when it is found in a polygon you 'pop' it out, grab the polygon id and move on to the next polygon 

Reply
0 Kudos
DanPatterson_Retired
MVP Esteemed Contributor

import numpy as np
import arcpy

def id_geom(in_fc):
    """The main code segment which gets the id and shape information and
    explodes the geometry to individual points
    """

    with arcpy.da.SearchCursor(in_fc, ['OBJECTID', 'SHAPE@']) as cursor:
        a = [[row[0], row[1]] for row in cursor]
    return a


def pip_arc(pnt_fc, poly_fc):
    """Naive point-in-polygon using arcpy
    """

    pnts = id_geom(pnt_fc)
    polys = id_geom(poly_fc)
    out = []
    for pnt in pnts:
        pid, p = pnt
        for poly in polys:
            plid, pol = poly
            if p.within(pol):
                out.append([pid, plid])  #[pid, p, plid, pol])
                break
        continue
    return out

# ----
pnt_fc = r"path_to\your.gdb\pnts_2000"      # points featureclass
poly_fc = "path_to\your.gdb\SamplingGrids"  # polygon featureclass

out = pip_arc(pnt_fc, poly_fc)  # ---- run pip_arc and convert to array
out = np.array(out)
dt = [('PntID', '<i4'), ('PolyID', '<i4')]  # ---- make an output array
z = np.ndarray((len(out), ), dtype=dt)      # to be used with ExtendTable
z['PntID'] = out[:, 0]
z['PolyID'] = out[:, 1]

arcpy.da.ExtendTable(pnt_fc, "OBJECTID", z, "PntID")  # --- make sure Pro is closed‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍
ThomasColson
MVP Frequent Contributor

Traceback (most recent call last):
File "C:\Users\bigfoot\Documents\ArcGIS\Projects\MyProject\PIP_TEST.py", line 33, in <module>
out = pip_arc(pnt_fc, poly_fc) # ---- run pip_arc and convert to array
File "C:\Users\bigfoot\Documents\ArcGIS\Projects\MyProject\PIP_TEST.py", line 23, in pip_arc
if p.within(pol):
AttributeError: 'NoneType' object has no attribute 'within'

Both FC"s are GCS NAD 83

Reply
0 Kudos
DanPatterson_Retired
MVP Esteemed Contributor

null point? or not a point featureclass?  it is a point in polygon demo, assume your are working local, adjust if otherwise

ThomasColson
MVP Frequent Contributor

yeah, checked that. 700+ points, each one of them is constrained by a polyogn (7 of them). Data is local FGDB created by Pro

Reply
0 Kudos
DanPatterson_Retired
MVP Esteemed Contributor

no clue, check again, 2000 pnts in 200 polygons... prior to arcpy.da.ExtendTable

out2
array([[   1,   30],
       [   2,   34],
       [   3,   30],
       [   4,   84],
       [   5,   78],
       [   6,   90],
       [   7,   18],
       [   8,   84],
       [   9,   16],
       [  10,   67],
       ...,
       [1991,   24],
       [1992,   96],
       [1993,   70],
       [1994,   41],
       [1995,   81],
       [1996,    8],
       [1997,   61],
       [1998,   29],
       [1999,   26],
       [2000,   67]])
ThomasColson
MVP Frequent Contributor

Ok I got it now, sort of. Turns out there WAS bad geometry, the OGC check in Pro found it. Now I got past that barrier, how does one use arcpy.da.ExtendTable to append the polyId to an existing field in the point table, as opposed to adding a new field? 

Reply
0 Kudos