nav[aria-label="Primary Navigation"] { padding: 0; & ul { list-style: none; width: 100%; display: flex; flex-direction: row; justify-content: start; align-items: start; gap: 30px; padding: 0; & li { margin: 0; } & ul li { list-style: none; } } }

5 million Point buffer analysis

GaneshmoorthiM

Hi,

I have 5 million points and I would like to generate buffer and find the points that falls in each buffer region and perform some calculation.

I have tried with Spatial join and getting out of memory exception and finally I am trying with Multiprocessing.

I would like to know, Is there any other way we can do this more efficiently in python or any other method.

Find more posts tagged with

Accepted answers

ChrisSnyder

The tool 'Generate Near Table' (BTW: only available at the ArcInfo license level) would also work for something like this... Basically it builds a distance-based one-to-many association table. A warning that the resulting output table is often many orders of magnitude larger than the number of input features - the explanation being that one house can be within 100 miles of many different restaurants. Be sure to apply a search radius! As a experiment, you might 1st apply a very small search radius (< 1 mile) so as to see the format of the output data you will be dealing with.... before you accidently end up with an output table that is 300 million records.

If you are running out of RAM, consider installing the 64-bit background geoprocessor: Background Geoprocessing (64-bit)—Help | ArcGIS for Desktop

Not sure exactly what you are doing here, but I think initially if I were doing this in a script, I would just loop through the restaurants, and select the houses, gather their stats, and write out some result. Something like:

foodPnts = r"C:\temp\test.gdb\restaurants"
housePnts = r"C:\temp\test.gdb\houses"
arcpy.MakeFeatureLayer_management(housePnts, "fl")
searchRows = arcpy.da.SearchRows(foodPnts, ["SHAPE@", "RESTAURANT_ID"])
for searchRow in searchRows:
  shapeObj, restaurantId = searchRow
  arcpy.SelectLayerByLocation_management("fl", "INTESECT", shapeObj, "100 MILES")
  recordCount = int(arcpy.GetCount_management("fl").getOutput(0))
  print("There are " + str(recordCount) + " houses within 100 miles of restaurant #" + str(restaurantId))
  #using the selected features/records in "fl", you now have a hook so as to get their names, incomes, etc.
  #blah, blah...
del searchRow, searchRows

All comments

ChrisSnyder

Yes, there is a way to improve performance... Like Dan mentioned, don't use ArcGIS (directly at least).

Here is another approach leveraging Python dictionaries. Note that this code only works for planar coordinates (UTM, State Plane, etc.) so if your data is in Geographic coordinates you need to use another distance function! BTW: This code assumes all fields are fully populated with valid values - if not you'll have to handle those exceptions.

Numpy has some nice ways to do all this stuff too - I usually prefer to roll my own though .

import math
def getDist(x1, y1, x2, y2):
   return math.sqrt((x1-x2)**2 + (y1-y2)**2)
sitesDict = {r[0],r[1]:r[2:] for r in arcpy.da.SearchCursor(sitesFC, ["SHAPE@X","SHAPE@Y","RESTAURANT_ID"])}
customersDict = {r[0],r[1]:r[2:] for r in arcpy.da.SearchCursor(custFC, ["SHAPE@X","SHAPE@Y","INCOME"])}
for x1, y1 in sitesDict:
  incomeList = []
  for x2, y2, in customersDict:
      if getDist(x1, y1, x2, y2) <= 1000: #or whatever map unit theashold...
         income = customersDict[(x2, y2)][0])
         if income > 0:
               incomeList.append(customersDict[(x2, y2)][0])
               #blah
   incomeListSum = sum(incomeList)
   incomeListLen = len(incomeList)
   print ("Restaurant #" + str(sitesDict[(x1,y1)][0] + " has " + str(incomeListLen) + " customers - average income is " = str(incomeListSum / float(incomeListLen)))

DanPatterson_Retired

Then you will need to prune your data set based on some criteria. We still don't know exactly what you are doing and why you need every site within your buffer.

GaneshmoorthiM

I am executing the script as a batch process. I have tried to minimise the candidates, but still it is huge data.

DanPatterson_Retired

As indicated previously, you will have to narrow your search from all the locations to those just near the point-buffer you are looking at. You won't be able to do all for all at once within an arcmap environment

GaneshmoorthiM

I have to run the analysis on daily and it is taking more time (for 95000 points it took almost 9 hrs). I will try with ArcGIS Pro.

one more question, Is there any other way to we can improve the performance

DanPatterson_Retired

not bad...unless you have to do it daily... If you have lots of memory (ie 16-32 Gig) and ArcGIS Pro, why don't you benchmark it in that since, it should run without alteration in that environment. It would be an interesting comparison. You are using the Python 3.4.x version of the print statement already, which is good, and I can't see anything else that would cause fault in the next python as well

print("some stuff {} plus more stuff {}".format(first_stuff,second_stuff))

takes care of type conversion the { } brackets are sequential. See Python's mini-formatting language for more details or one of my blogs on formatting.

GaneshmoorthiM

Hi Dan, Chris,

The scenario which I explained in my previous post was example scenario.

I am using following code to find the customer within 200-mile radius from our distribution center.

It took almost 2 hr to analyse 250 distribution center points against 0.4 Millon customer locations. Is there anyway we can improve the performance

import arcpy

DC_Points = r"C:\temp\TestGDB1.gdb\DistributionSites"  
CL_Points = r"C:\temp\TestGDB1.gdb\CustomerLoaction"  
fcSumPoint=r"C:\temp\TestGDB1.gdb\fcSumPoint"
arcpy.MakeFeatureLayer_management(CL_Points, "cl_points")  
searchRows = arcpy.da.SearchCursor(DC_Points, ["SHAPE@", "site_id","cov1val"])
print("searchcursor")
fields=("SHAPE@","dc_Id","cl_count","acSum")
inCursor = arcpy.da.InsertCursor(fcSumPoint, fields)
print("insert cursor")
for searchRow in searchRows: 
  shapeObj,dc_id,AcSUm = searchRow
 buffSelectionLyr=arcpy.SelectLayerByLocation_management("cl_points", "INTERSECT", shapeObj, "20 MILES","NEW_SELECTION")  
  recordCount = int(arcpy.GetCount_management("cl_points").getOutput(0)) 
  sumValues = 0 
  sumValues += AcSUm
 print(siteid,recordCount,sumValues)
 infiled=[shapeObj,siteid,recordCount,sumValues]
 inCursor.insertRow(infiled)
  print("There are " + str(recordCount) + " points in #" + str(siteid)) 
  #using the selected features/records in "flyou now have a hook so as to get their names, incomes, etc.  
del searchRow, searchRows, inCursor

ChrisSnyder

If you are running out of RAM, consider installing the 64-bit background geoprocessor: Background Geoprocessing (64-bit)—Help | ArcGIS for Desktop

foodPnts = r"C:\temp\test.gdb\restaurants"
housePnts = r"C:\temp\test.gdb\houses"
arcpy.MakeFeatureLayer_management(housePnts, "fl")
searchRows = arcpy.da.SearchRows(foodPnts, ["SHAPE@", "RESTAURANT_ID"])
for searchRow in searchRows:
  shapeObj, restaurantId = searchRow
  arcpy.SelectLayerByLocation_management("fl", "INTESECT", shapeObj, "100 MILES")
  recordCount = int(arcpy.GetCount_management("fl").getOutput(0))
  print("There are " + str(recordCount) + " houses within 100 miles of restaurant #" + str(restaurantId))
  #using the selected features/records in "fl", you now have a hook so as to get their names, incomes, etc.
  #blah, blah...
del searchRow, searchRows

DanPatterson_Retired

after limiting the candidates either method should work. The numpy approach with calculating the distances can be completed using a np.where statement querying for only those locations that have a distance less than your threshold distance

GaneshmoorthiM

thanks for the reply. I have already tried to limit the candidate using following approaches

-used Fishnet (To find buffer regions)

-Creating bounding box around the circle

-Find an Identical location (to eliminate duplicate points(locations))

Currently, I am trying to find the points using Generate Near Table—Help | ArcGIS for Desktop Because it supports finding the more than one near features.

I am currently using arcpy.da.TableToNumPyArrayTableToNumPyArray—Data Access module | ArcGIS for Desktop to store the results in a separate table.

The use case which i explained was a sample one and we have similar different use cases with reasonable use stories.I will try all your above suggestions to restrict the points and will post the update once I am done.

Thank you.

DanPatterson_Retired

In addition to the previous long post...other aspects to consider

what makes you think that people from 100 miles around are going to your restaurant as their primary destination? Perhaps, the target radius is too big. have you performed any analysis on-site that indicates that "people are coming from miles around"... and if so, how many?
the target audience within a 100 mile buffer may bear no resemblence to those within the immediate neighbourhood of your restaurant. Are the socio-economic conditions uniform within the landscape? Will they all likely have the ability to dine out? are they willing to travel large distances?
You are only considering euclidean distance currently... a good first step, however, the actual travel distance by road can be a totally different story.

So these non-distance characteristics can be further used to limit the candidates for the study, pruning them down by attribute and by spatial location making the problem more manageable.

GaneshmoorthiM

Hi Dan,

I am performing proximity overlay analysis to find the points within circle.

Here is the sample senario, I would like to find the houses that are around my restaurant using 100 mile buffer.

As per below image:

Green color point is my restaurant and brown color points are houses.

I have performed proximity analysis to find the houses around my restaurant in 100 mile radius.

Based on 100 mile radius buffer i will select the house and perform some additional calculations.

I am able to perform this analysis successfully with around 2 Million data using spatial join.

the problem is when I try to execute this analysis with 5 million data i am facing out of memory exception. Is there any other option to perform this proximity analysis more efficiently in large data.

DanPatterson_Retired

have a look at kernel density Kernel Density—Help | ArcGIS for Desktop

point density Point Density—Help | ArcGIS for Desktop

optimized hotspot Optimized Hot Spot Analysis—Help | ArcGIS for Desktop

DanPatterson_Retired

that is the point I am not clear on...you say you have 5 million points that you want to buffer and then you want to do a spatial query with 5 million other(?) buffers?

where are those buffers?
what is the purpose of the buffer?
are you looking for some kind of density map? or a proximity map?
do you have access to the spatial analyst?
can you show an example using a smaller data set of what you are trying to do since the problem needs clarification

GaneshmoorthiM

Hi Dan, thanks for the quick reply.

Below is the sample image, The layer has 5Million points and i have to generate 5 Million buffers. After Generating buffer i have to perform spatial query similar point layer which has same amount of data.

DanPatterson_Retired

how many buffers? the obvious thing would be to spatially tile your data since it will be apparent which points can possibly fall within each buffer. multiprocessing, hadoop ( Esri GitHub | Open Source and Example Projects from the Esri Developer Platform ) will probably not be the way to go, since the problem needs to be formed within a different spatial context

Join the Conversation

Signed in members can post, follow updates, and more. New here? Register a free account.

Quick Links

Community Help
Find useful guides, FAQs, and documents to help you navigate and make the most of Esri Community.