I have a 4.29 GB Shapefile from which I want to extract features as GeoJSON, defined by a bounding box.

This should be as fast as possible, because I want to use it in a web application.

My first attempt at this was to convert to SQL using shp2psql and import it into PostGIS. I created a spatial index (gist) and did a full ANALYZE VACUUM on the table. The VM has 4 GB of memory.

The following query works, but took 623 seconds(!):

    ST_AsGeoJSON(geom)AS geometry,
    geom && ST_GeogFromText(
        'POLYGON((5.117705762386322 52.09243663877818,5.118895322084427 52.09243663877818,5.118895322084427 52.09314276719956,5.117705762386322 52.09314276719956,5.117705762386322 52.09243663877818))'
AND ST_Intersects(
        'POLYGON((5.117705762386322 52.09243663877818,5.118895322084427 52.09243663877818,5.118895322084427 52.09314276719956,5.117705762386322 52.09314276719956,5.117705762386322 52.09243663877818))'

An equivalent using ST_DWithin took even longer.

I decided to try an Amazon RDS PostgreSQL 9.3.1 instance with 68 GB of memory and 26 CPU's. Query time dropped to around 30 seconds. Faster but still not fast enough.

Preprocessing will take forever, even on Amazon RDS. And it will cost a lot of money.

I also tried Fiona and 'ogr2ogr -spat', but they're about as slow as PostGIS.

What other options exist for quickly getting GeoJSON out of a large Shapefile for a certain spatial extent?


  • 10 minutes to 30 seconds, such an improvement for a 4 GB shape. How many seconds would you expect for such big file? – Gery Jan 20 '14 at 15:23
  • Well, I thought that maybe if it'll fit into memory entirely (on Amazon) it'd be near instantaneous.. No swapping etc. But perhaps it's a more difficult operation than I imagine. – Gijs Jan 20 '14 at 16:00
  • Perhaps I should take another approach and make GeoJSON files per area and load/unload these with Leaflet as the user pans? – Gijs Jan 20 '14 at 16:01
  • that's the way, divide and conquer will apply nicely here. – Gery Jan 20 '14 at 16:02
  • Your sample query is duplicative (the && test is embedded inside ST_Intersects already) but that's not the problem. "it's big and slow!" problems usually boil down to: you're pulling a large chunk of data out (query isn't very selective) or; the objects in the table are very large and/or overlapping. And there's no magic bullet solution to those problems (well, the latter problem can be partially solved by chopping put big stuff into smaller stuff first) – Paul Ramsey Jan 22 '14 at 20:35

Since you will go for it, I would like to suggest you some links that may give you some guidance to do it so, here they are:

Definitely you need to adapt them to your needs, but I think it is good to have these links at hand for more people having the same issue you have.

Hope this helps,


What level of resolution do you need the geojson at? One approach could be to simplify the shapefile using ST_SimplifyPreserveTopology() or MapShaper (though it might be too big ) to identify the polygons. You could then return the simplified polygons as json or select them from the original table using the IDs.


Solved! :)

As Paul Ramsey points out in the comments above, the query was clumsy to begin with. On top of that, ST_Intersects turns out to be slow (or not the right tool for this job)

A combination of && and BOX3D works REEEALLY fast. I'm talking <1 second for Z16 and up. That's what I was looking for.

For reference:

SELECT ST_AsGeoJSON(geom) AS geometry FROM shapes WHERE geom && ST_SetSRID('BOX3D(" + bbox + ")' :: box3d, 4326)

Paul, I can't accept your answer since it's a comment.

Thanks a lot!

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.