# How can I perform a Join by Location for a large vector and very large point file?

I have no programming skills, but am thinking I may need to at least figure out how to run this one script from inside QGIS.

I have a vector file of global quarter degree grid cells and a point file containing an incredibly large number of points with population values.

I want to create a polygon file that has the sum of the population values (not just the count of the points) for each grid cell.

This is really just a basic Join by Location (Spatial Join), but the files are so large it will not work.

Can anyone help me?

• What data format are the layers (shapefile, sqlite, PostGIS...)? – artwork21 Oct 4 '17 at 16:57
• I believe they are shapefiles. One is a vector grid I made in QGIS and the other is a point file I generated in ArcMap – Arthur_15 Oct 4 '17 at 17:07
• I love that people are viewing this and commenting. gives me hope that someone might have an answer for me :) – Arthur_15 Oct 4 '17 at 19:46
• Maybe convert it into a sqlite/spatialite db or PostgreSQL/PostGIS layer where you can take advantage of spatial indexes or split out your heavy data into smaller chunks/layers. – artwork21 Oct 4 '17 at 19:53

I find for your particular scenario - points versus a regular grid - I go for a programatic approach (e.g. Python).

I do this for a couple of reasons:-

• it's faster than a traditional spatial join, as it's much easier to test for intersections with cells.
• unlike a database lookup or spatial join, it's easier to monitor progress :)

Also, there's no need to consider every 1/4 degree cell if the data is 'sparse'. The world has 1,036,800 of those cells. 70% + of the cells will be empty (sea!) and most of the remaining land cells will be empty, too. That's a lot of wasted effort, however you do this - whether you use shapefiles or a database.

This code uses a CSV file for the points. I used the NASA Meteorite Database (from here) if you want to try this code out. In my case I'm summing up the weight of meteorites which land in each quarter-degree grid cell.

To keep things simple, it's standalone python and doesn't need any special libraries.

import csv

sums = {} # key=(cell_x,cell_y), value=sum

# load points file, this is a CSV (comma-separated) with long, lat and values in each record

with open("/tmp/Meteorite_Landings.csv","r") as fi:
count = 0
for line in reader:
try:
lat = float(line[7]) # 8th column in my case for latitude
lon = float(line[8]) # 9th column in my case for longitude
val = float(line[4]) # value to sum is in 5th column
cell_x = int((lon + 180.0) * 4) # 4 cells per degree
cell_y = int((lat + 90.0) * 4)
key = (cell_x, cell_y)
if not key in sums:
sums[key] = val
else:
sums[key] += val
except:
pass # ignore header line, lines with errors

# write this out to a TAB separated file.
# can load this into QGIS as a delimited layer.

with open("/tmp/foo.csv","w") as fo:
fo.write("GEOM\tSUM\n")
for key in sums:
cell_x, cell_y = key
total = sums[key]
# convert cell number back to lat,lon of SW corner of cell
nw_x = sw_x = -180.0 + (0.25 * cell_x) # convert to longitude of SW corner
se_y = sw_y = -90.0 + (0.25 * cell_y)  # convert to longitude of SW corner
ne_x = se_x = sw_x + 0.25
nw_y = ne_y = sw_y + 0.25
poly = "POLYGON(({} {},{} {},{} {},{} {},{} {}))".format(sw_x,sw_y,nw_x,nw_y,ne_x,ne_y,se_x,se_y,sw_x,sw_y)
fo.write("{}\t{}\n".format(poly, total))

print("Finished!")


The output is CSV which you can then import into QGIS as a delimited layer. Only those cells with a value exist using this approach..

This summed up 30k records against a 1/4 degree grid in under a second, on a 10-year old laptop.

CSV isn't exactly the best format for huge amounts of data, but I'm guessing your data is sparse enough that this doesn't matter.

I'd try a QGIS 3.0 nightly release first - Select/Extract by Location are many orders of magnitude faster in 3.0.

• thanks for these answers. Excited to try it out. I'm skeptical of my ability to use Python, but maybe I will be able to figure it out. – Arthur_15 Oct 5 '17 at 16:05
• do you have a link where I can download 3.0? – Arthur_15 Oct 5 '17 at 16:11
• Is it actually 2.99 and is there an OS version? – Arthur_15 Oct 5 '17 at 16:29
• Does anyone know where I can download QGIS 3.0? – Arthur_15 Oct 6 '17 at 19:04