Buffer and union of large shapefile?

Question

I want to create a buffer of 5km around the coastline, and save that as a separate shapefile.

I have a shapefile with almost 180,000 features. It's a line shapefile, representing the coastlines of the world (I got the data here). The shapefile is in CRS (EPSG:4326, WG84).

So far I have been testing three approaches and all are either very slow or crash:

1) QGIS 2.14.2 Essen (with functioning GRASS functions) I go to 'Vector' --> 'Geoprocessing Tools' --> 'Buffer(s)...' I select the shapefile, For buffer distance I use: 0.05 (EPSG:4326 which is about 5.5km) I check the 'Dissolve buffer results' options. This option doesn't seem to work, as QGIS crashes. I expect it has to do with a shortage of memory available. "Error: minidump written to (filepath)" [

2) PgAdminIII 1.20.0 with the following state I try to create a new table and fill it with the buffer:

CREATE TABLE public.lines_buffer AS SELECT ST_Buffer(geom,0.05) FROM public.lines

But the connection keeps on getting lost and thus breaking the operation. The furthest I got was 2 hours of query time, but that wasn't enough.

3) Python 2.7 (shapely, fiona, descartes) My third attempt, is still in progress. I am trying to loop through the large shapefile and make a buffer. I'm not sure if my code is very efficient. I still have to "dissolve" after the buffer, which would also take a while. I plan on testing shapely.ops.cascaded_union. EDIT: Even running the first 2 items takes very long...[Spyder IDE: 100% memory usage] My current loop just limits the first 3 items in the fiona.collection Credits go to Tom Macwright

def readshape(pathfolder, filename):
filepath = os.path.join(pathfolder, filename)
with collection(filepath) as Input:
    schema = Input.schema
    with collection('buffer_shape.shp', 'w', 'ESRI Shapefile', schema) as output:
        for i in range(0, 3):
            output.write({
                'properties': {
                    'source': Input[i]['properties']['source']
                },
                'geometry': mapping(shape(Input[i]['geometry']).buffer(0.05))
            })

Does anyone have any tips, which software is the quickest and most stable option?

Answer 1

PostGIS should work for your task but:

1. Reproject your data into some Cartesian coordinate system (like Web Mercator) because computations on the WGS84 ellipsoid (geographic coordinates) in PostGIS take longer time. Also for geography ST_Buffer may not behave as expected if object is sufficiently large that it falls between two UTM zones or crosses the dateline, which is true in your case (PostGIS docs).
2. Create a spatial index for your features:

CREATE INDEX IF NOT EXISTS lines_idx ON public.lines USING GIST (geom);

3. Run your query with buffer distance in meters (if you chose a metric CRS):

CREATE TABLE public.lines_buffer AS SELECT ST_Buffer(geom,5500) FROM public.lines;

Answer 2

I've just tried to use the QGIS Vector Geometry tools-> Fixed distance buffer on GSHHS_I_L1 and it completed in less than a second.

So I suspect you are using a much higher resolution version of the shapefile (There is no need if you are studying the whole world) or you are on a much smaller machine than I am.

Answer 3

I doubt you can get much faster than PostGIS. You can script it to run in batches (in your favourite shell or language), so the number of total features won't be a problem. Just create the table and update it in steps, only handling a few features at a time (using LIMIT and OFFSET).