topic Re: 5 million Point buffer analysis in Python Questions

5 million Point buffer analysis

GaneshmoorthiM — Fri, 25 Dec 2015 14:53:25 GMT

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.

Re: 5 million Point buffer analysis

DanPatterson_Retired — Fri, 25 Dec 2015 15:39:14 GMT

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

Re: 5 million Point buffer analysis

GaneshmoorthiM — Sat, 26 Dec 2015 00:09:11 GMT

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.

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 26 Dec 2015 00:53:47 GMT

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

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 26 Dec 2015 01:05:51 GMT

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

Re: 5 million Point buffer analysis

GaneshmoorthiM — Sat, 26 Dec 2015 06:10:54 GMT

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.

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 26 Dec 2015 08:11:36 GMT

first, since you can see your buffer, you could do one of the following

roughly draw a rectangle around your buffer using the selection tool to select the points that are within your buffer zone... selections are used when using most tools in ArcToolbox. That will limit the candidate sites.
you could perform the same selection and save the results to a new file for further analysis
you seem to be interested in euclidean distance and not network distance so that further simplifies the analysis.
you can simply perform a select by spatial location using the previously rectangular area points using them as the target and the buffer circle as the selector. This will give you all the points that are within the buffer circle
if you are really stuck on the idea of doing distance measurements, then determine the points that are within the buffer radius distance using numpy Before I forget ... # 16 ... NumPy vs SciPy... making a point this is extreme overkill and not needed.

In short... don't rely solely on the software to solve the problem...human brain intervention often simplifies problems to acceptable levels...since you can see the points that are within the buffer...select them roughly, then fine-tune the selection using arctoolbox tools...

select by spatial location in the analysis,distance toolset
- Select Layer By Location—Help | ArcGIS for Desktop
the Point distance or Near tool if you have the advanced license
- Point Distance—Help | ArcGIS for Desktop
- Near—Help | ArcGIS for Desktop
or numpy, as per the solution in the link I provided.

So....

limit your candidate sites in a quasi-interactive way or use tools which can handle large point masses, numpy can be used to substitute for tools that you don't have the necessary license for

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 26 Dec 2015 08:28:11 GMT

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.

Re: 5 million Point buffer analysis

GaneshmoorthiM — Sat, 26 Dec 2015 13:47:30 GMT

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.

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 26 Dec 2015 15:51:27 GMT

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

Re: 5 million Point buffer analysis

ChrisSnyder — Sat, 11 Dec 2021 18:18:03 GMT

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

Re: 5 million Point buffer analysis

GaneshmoorthiM — Sat, 11 Dec 2021 18:18:06 GMT

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

Re: 5 million Point buffer analysis

DanPatterson_Retired — Tue, 29 Dec 2015 20:18:45 GMT

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.

Re: 5 million Point buffer analysis

GaneshmoorthiM — Fri, 01 Jan 2016 23:41:46 GMT

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

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 02 Jan 2016 00:28:02 GMT

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

Re: 5 million Point buffer analysis

GaneshmoorthiM — Sat, 02 Jan 2016 00:37:12 GMT

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

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 02 Jan 2016 00:49:05 GMT

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.

Re: 5 million Point buffer analysis

GaneshmoorthiM — Sat, 02 Jan 2016 04:48:52 GMT

Each candidate in the buffer zone has its own statistical value. So, I have to find all candidates in given buffer, get each candidate value and perform additional statistical calculation like sum. It is points within polygon analysis. the only problem which I am facing is performance and out of memory issue.

Re: 5 million Point buffer analysis

DanPatterson_Retired — Sat, 02 Jan 2016 05:19:23 GMT

Then you need to delimit the candidates whose values change.

Step 1 would be to take a snapshot and determine the "value" for each candidate at a point in time. That is your benchmark. Do the long hard work of deriving he benchmarks for the buffers. This is easyto do. You could obtain millions of values in microseconds using arcpy and/or numpy (see my previous posts).
Step 2 get the statistical parameter for all candidates at time step 2. Ditto on the quick thing.
Step 3 determine the candidates that changed... a simple difference... incredibly fast
perform your buffer query on only the candidates in step 3 ... should be fast for all buffers since your candidates of interest are now greatly reduced.
final step, take the changes per buffer and apply them to the benchmark

So if your buffer size and locations remain the same, the problem can be greatly reduced to a smaller set by only querying for the 'changelings' after your benchmark conditions are determined.

Re: 5 million Point buffer analysis

ChrisSnyder — Sat, 11 Dec 2021 18:18:08 GMT

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)))