POST
|
Spatial join might work, but does it work for overlapping polygons? It is critical that the points are accurately counted in my case. I adjusted a few parts of my script and it has sped it up considerably. It's still not fast enough to process 300,000 records in under a few hours, but I thought I'd post this improved version in case it was helpful to someone else out there. I still don't understand why this simple procedure has to take so long! :confused: #-------------------------------------------------------------------------------
# Name: Count Points in Polygon
# Purpose:
#
# Author: hcopeland
#
# Created: 22/08/2011
# Copyright: (c) hcopeland 2011
# Licence: <your licence>
#-------------------------------------------------------------------------------
#!/usr/bin/env python
import arcpy
from arcpy import env
#getcontext().prec = 8
ws = r'D:\Work\Projects\SageGrouseEasementScenarios2011\Data\FinalModelData.gdb'
arcpy.env.workspace = ws
# set overwrite outputs to true
arcpy.env.overwriteOutput = 1
POINT = 'D:\Work\Projects\SageGrouseEasementScenarios2011\Data\Final_Development_Buildouts\Cumulative_buildouts.gdb\Cumulative_Baseline'
POLYGON = 'LekParcelFC'
# Converts unselectable feature class into a selectable feature layer
arcpy.MakeFeatureLayer_management(POINT,"point")
rows = arcpy.SearchCursor(POLYGON)
print 'Calculating Points in Polygons...'
### Loop through each row and count the points in each polygon record
try:
for row in rows:
expression2 = "\"OBJECTID\" = " + str(row.OBJECTID) + ""
#expression = "\"OBJECTID\" = " + str(row.OBJECTID) + ""' AND "PARCELNB" <> ' + "'""'"
print expression2
# Select each record inside of the polygon feature class
arcpy.MakeFeatureLayer_management(POLYGON,"polygon2",expression2)
# Select all the points that are inside of the polygon record
arcpy.SelectLayerByLocation_management("point", "intersect", "polygon2", 0, "NEW_SELECTION")
# Count the points that are in each polygon
GetCount = arcpy.GetCount_management("point")
# Calculate the ASSOC_PTS field with the counted points
arcpy.CalculateField_management("polygon2", "CEFeaturesPreventedBaseline", GetCount, "VB", "")
except:
arcpy.GetMessages(2)
arcpy.AddError("Script Bombed")
print 'Count Points in Polygon Script Finished!'
... View more
10-11-2011
05:21 AM
|
0
|
0
|
419
|
POST
|
I wrote this script to iterate through a polygon fc and count the number of points in each polygon and put that information in a new field. The script works, but it is REALLY slow. If I use Hawth's Tools to do the same thing, it takes 20-30 minutes. With this script, it was set to take 4 days. Any ideas what I've done wrong--I have a feeling it's related to the search cursor. Should I not be using the select by location/select by attribute ESRI functions? Thanks, Holly #-------------------------------------------------------------------------------
# Name: Count Points in Polygon
# Purpose:
#
# Author: hcopeland
#
# Created: 22/08/2011
# Copyright: (c) hcopeland 2011
# Licence: <your licence>
#-------------------------------------------------------------------------------
#!/usr/bin/env python
import arcpy
from arcpy import env
#getcontext().prec = 8
ws = r'D:\Work\Projects\SageGrouseEasementScenarios2011\Data\FinalModelData.gdb'
arcpy.env.workspace = ws
# set overwrite outputs to true
arcpy.env.overwriteOutput = 1
POINT = 'D:\Work\Projects\SageGrouseEasementScenarios2011\Data\Final_Development_Buildouts\Cumulative_buildouts.gdb\Cumulative_Unconstrained'
POLYGON = 'LekParcelFC'
print 'Making Multipart into Single Part....'
#Convert multipart to singlepart
#arcpy.MultipartToSinglepart_management(POLYGON,"LekParcelFCSingle")
# Converts unselectable feature class into a selectable feature layer
arcpy.MakeFeatureLayer_management(POINT,"point")
arcpy.MakeFeatureLayer_management(POLYGON,"polygon")
rows = arcpy.SearchCursor("polygon")
print 'Calculating Points in Polygons...'
# Loop through each row and count the points in each polygon record
try:
for row in rows:
# Select each record inside of the polygon feature class
#expression = "\"OBJECTID\" =" + str(row.OBJECTID) + ' AND "PARCELNB" <> ' + "''"
expression = "\"OBJECTID\" = " + str(row.OBJECTID) + ""
print expression
SelPoly = arcpy.SelectLayerByAttribute_management("polygon", "NEW_SELECTION",expression)
# Select all the point that are inside of the polygon record
SelPts = arcpy.SelectLayerByLocation_management("point", "WITHIN", SelPoly, 0, "NEW_SELECTION")
# Count the points that are in each polygon
GetCount = arcpy.GetCount_management(SelPts)
# Calculate the ASSOC_PTS field with the counted points
arcpy.CalculateField_management("polygon", "CEFeaturesPreventedHigh", GetCount, "VB", "")
except:
arcpy.GetMessages(2)
arcpy.AddError("Script Bombed")
print 'Count Points in Polygon Script Finished!'
... View more
10-06-2011
02:44 PM
|
0
|
3
|
2139
|
POST
|
Thanks, Stacy. I don't think the code you sent was quite what I was looking for (but I could be wrong as I am relatively new to Python). This code from a friend was more along the lines of what I needed. Thanks, though! import arcpy
arcpy.env.workspace = r'C:\WorkSpace\Projects\Code\Special Request\Copeland\src\SageGrouseTestCode.gdb'
parcelIDTable = 'ParcelNBTable'
lekFC = 'LekParcelFC'
try:
parcelRows = arcpy.SearchCursor(parcelIDTable, "", "", "ParcelNBTable_PARCELNB", "")
for parcelRow in parcelRows:
print 'Parcel Number: %s' %parcelRow.ParcelNBTable_PARCELNB
try:
if len(parcelRow.ParcelNBTable_PARCELNB) > 0:
expression = '"PARCELNB" = ' + "'" + parcelRow.ParcelNBTable_PARCELNB + "'"
lekRows = arcpy.SearchCursor(lekFC, expression, "", "", "Complex A")
for lekRow in lekRows:
print '\tLek complex name: %s\tCount %i' %(lekRow.Complex, lekRow.MAX_MostRecentPeakCount)
except:
print arcpy.GetMessages(2)
except:
print arcpy.GetMessages(2)
print 'fin'
... View more
07-26-2011
10:15 AM
|
0
|
0
|
180
|
POST
|
Hi Folks, I've searched the archives and ESRI help for example code that does a relate using cursors and can't seem to find one. Here's my problem. Perhaps someone can share code or give an example of how you would code this? It seems like this should be simple enough... I have a ArcGIS data table with an id code "parcelid". I thought that I would use a SearchCursor to walk through each row and then relate through "parcelid" to a featureclass to grab the records that match that parcelid. I then need to do a series of calculations on fields in the table and featureclass and put the result in a new table. Any help would be greatly appreciated! Thanks, Holly
... View more
07-24-2011
06:59 AM
|
0
|
2
|
815
|
POST
|
This code isn't nearly as elegant as yours, Kim, but perhaps it is useful to someone else, so I'm posting it. Holly CODE: # --------------------------------------------------------------------------- # Select Random Parcels.py # Created on: 2011-5-18 # Creator: Holly Copeland (hcopeland@tnc.org) # Requirements: Spatial Analyst Extension # Description: This script randomly selects data from a feature class and # stops selecting when it reaches a threshold value set by the user # It also iterates multiple times and generates summary statistics from a # separate field at the end # --------------------------------------------------------------------------- # Import arcpy module import arcpy, os, sys, random, array from arcpy.sa import * from arcpy import env # Check out the ArcGIS Spatial Analyst extension license arcpy.CheckOutExtension("Spatial") arcpy.env.workspace = "D:\\Work\\SageGrouseEaseStudyBackup\\Data\WYParcels_Merge.gdb" # #Specify Input Parameters # fc = "WYParcelsPrivate_AgValuesMay2011" scenarioValue = 5000000 #this is the THRESHOLD value iterationCount = 0 iterationNumber = 5 #number of iterations sgSumTotal = [] sgMeanTotal = [] while iterationCount < iterationNumber: #Set parcel accumulation value to 0 accumulatorValue = 0 objectIDList = [] #Set sage grouse accumulation value to 0 sagegrouseValue = 0 sagegrouseList = [] sagegrouseListSum = [] # Create a search cursor # # The variables row and rows are initially set to None, so that they # can be deleted in the finally block regardless of where (or if) # script fails. # row, rows = None, None #Start while statement to accumulate parcels #Value in while statement should be the same as the scenario while accumulatorValue < scenarioValue: objectID = random.randrange(1,27577) #generate random number to select parcel object id expression = '"OBJECTID" = ' + str(objectID) try: rows = arcpy.SearchCursor(fc, expression) #select the objectid from above #Iterate through the rows in the cursor of randomly selected parcels for row in rows: #assuming you are sampling without replacement if row.OBJECTID not in objectIDList: accumulatorValue += row.AdjustedParcelValue sagegrouseValue += row.SageGrouseTest objectIDList.append(row.OBJECTID) sagegrouseList.append(row.SageGrouseTest) except: arcpy.AddError("Script Bombed") #calculate statistics for each iteration sgSum = sum(sagegrouseList) sgMax = max(sagegrouseList) sgMin = min(sagegrouseList) sgMean = sum(sagegrouseList)/len(sagegrouseList) sgSumTotal.append(sgSum) sgMeanTotal.append(sgMean) print "Sage Grouse Statistics for Run #" ,iterationCount print "No of parcels:" ,len(objectIDList) print "Max:" ,sgMax print "Min:" ,sgMin print "Sum:" ,sgSum print "Mean:" ,sgMean iterationCount += 1 sgSumTotalFinal = sum(sgSumTotal) sgMeanTotalFinal = sum(sgMeanTotal)/len(sgSumTotal) print "Final Sage Grouse Statistics (Sum, Mean):" ,sgSumTotalFinal, sgMeanTotalFinal
... View more
05-20-2011
02:02 PM
|
0
|
0
|
746
|
POST
|
Thank you, Kim. This is very helpful!!! I'll post my final code too when I get it finished. Cheers, Holly
... View more
05-20-2011
12:05 PM
|
0
|
0
|
746
|
POST
|
Thank you Kim for your advice. I am thankful to hear that this seems doable to you and straightforward. Do you have any example code or quick specific code hints for using a cursor and accumulating random selections from the attribute until I get to the target value? Thank you -- Holly
... View more
05-18-2011
10:47 AM
|
0
|
0
|
746
|
POST
|
Hi all, I need to write a script to do the following: Randomly select polygons (parcels in my case) in a feature class using an attibute with parcel value. It needs to keep selecting parcels randomly until it reaches a total max parcel value. So, if each parcel is valued between 1-100, and I have a max value of 500, it needs to select enough parcels to reach the target value of 500, then I'll do some operations/calculations on these selected parcels. Finally, I need to create a loop to do this 100 or so times. Is what I propose above doable and fairly straightforward. I'm relatively new at Python programming, so any tips about the logic in the code would be so helpful and appreciated. Thanks, Holly
... View more
05-13-2011
11:42 AM
|
0
|
5
|
2773
|
Online Status |
Offline
|
Date Last Visited |
11-11-2020
02:24 AM
|