How to get feature count of unique values using a point and polygon feature class?

09-19-2019 08:21 AM
New Contributor

I'm using a search and update cursor to extract values from a point feature class to another feature class.  The process (see below) is outlined like this:

  1. Use point feature class 1 (ServiceConnections) to select polygon (Parcels) by location
  2. Use the selected Parcels to select point feature class 2 (CustomerLocations) using a SQL query
  3. Use the PT_ACCOUNT_ID field from CustomerLocations and the PG_PARCEL_ID field from Parcels to update the ACCTID and PARCELID of the original point feature class 1 (ServiceConnections).

This all works fine if the ServiceConnections were all in a single parcel.  But because some of parcels are not drawn in correctly, or haven't been updated, I know there are multiple points in one parcel.

My questions is this, how can I use the ServiceConnections to select Parcels, but identify if the parcel contains multiple ServiceConnections?  Ideally, I'd like to continue the script to process those that are not problematic and skip the ones that have multiple ServiceConnections by populating a field 'MultipleConnectionsInParcel' in the ServiceConnections feature class with a Yes or No.

Here is the script I have so far:

import arcpy

# Set local variables
mxd = arcpy.mapping.MapDocument("C:\Meters.mxd")
df = arcpy.mapping.ListDataFrames(mxd, "*")[0]
serviceConnections = arcpy.mapping.ListLayers(mxd, "Water Service Connections", df)[0]
parcels = arcpy.mapping.ListLayers(mxd, "Parcels", df)[0]
customerLocations = arcpy.mapping.ListLayers(mxd, "CustomerLocations", df)[0]

# Select Parcels using the Water Service Connections
with arcpy.da.SearchCursor(serviceConnections, ["ACCOUNTID"]) as cursor:
    serviceConnectionsCount = 0
    for row in cursor:
        # Count number of water service connections that's about to run
        serviceConnectionsCount += 1
        print("Selected Water Connections Count: " + str(serviceConnectionsCount))
        # Select parcels using the location of the Water Service Connections
        arcpy.SelectLayerByLocation_management(parcels, "INTERSECT", serviceConnections, "", "NEW_SELECTION", "")

        with arcpy.da.SearchCursor(parcels, ["SITEADDRESS", "PG_PARCEL_ID"]) as parCursor:
            parcelCount = 0
            for row in parCursor:
                # Count the number of parcels selected
                parcelCount += 1
                print("Parcel Count: " + str(parcelCount))
                addressValue = row[0]
                parcelidValue = row[0]
                print("The address is: " + addressValue)
                print("The parcel id is: " + parcelidValue)

            where_clause = "USER_SERVICE_ADDRESS = " + "'" + addressValue + "'"
            with arcpy.da.SearchCursor(customerLocations, ["PT_ACCOUNT_ID"], where_clause) as wcpCursor:
                customerPointCount = 0
                for row in wcpCursor:
                    # Count the number of customer locations
                    customerPointCount += 1
                    print("Customer Point Count: " + str(customerPointCount))
                    accountID = row[0]
                    print("The Account ID is: " + accountID)

                # Update the Water Service Connections feature class with the account id, parcel id, and if the parcel
                # contains multiple service connections
                with arcpy.da.UpdateCursor(serviceConnections, ["ACCOUNTID", "PARCELID", "MultipleConnectionsInParcels"]) as updater:
                    for row in updater:
                        row[0] = accountID
                        row[1] = parcelidValue
                        row[2] = ???????
0 Kudos
1 Reply

If I understand correctly. You are looking for polygons (CustomerLocations) that contain more than one point (ServiceConnections)

You could try using contains in the geometry class.

Geometry—ArcPy classes | ArcGIS Desktop 

Essentially how it would work you is have a search cursor iterate through the polygons then nest another loop inside of the polygons loop that would compare the geometry. Then you count the amount of points that are contained within the polygon and if its greater than 1 write to the 'MultipleConnectionsInParcel' field.

When using the search cursor make sure to use "SHAPE@" as one of the field names so it uses the geometry.


Here is a snipit of code that shows how I use the search cursors to compare geometry for two feature classes against each other.

for buffersLoop in arcpy.da.SearchCursor(fatStructures, ["SHAPE@", Struc_Num]):
    for photosLoop in arcpy.da.SearchCursor(photosFC, ["SHAPE@", Photo_Name]):
        polyprint = arcpy.Geometry.contains(buffersLoop[0],photosLoop[0])‍‍‍‍‍‍‍‍‍
0 Kudos