10

I'm having real problems cleaning up and reducing the size of a Shapefile which is published by the UK's Environment Agency.

It shows the extents of the LiDAR open data they publish: each polygon is a survey flight, and there are fields for the resolution, the date of the flight and so on.

Eventually I'd like to end up with a dissolved (merging all features) and simplified Shapefile which is much smaller in size. At the moment:

  • The .shp is 360.6MB
  • The .dbp is 26.9MB
  • There are 121,753 polygons

I think one of the reasons the file is so large is that there are many small 'specks' of data (which aren't important for my purposes):

enter image description here

What I've tried so far:

  1. Dissolving with QGIS: this seemed to make no progress on such a big file so I cancelled it after a while
  2. Dissolving with OGR (ogr2ogr dissolved.shp LIDAR_composite_extents_2015.shp -dialect sqlite -sql "SELECT ST_Union(geometry) AS geometry FROM LIDAR_composite_extents_2015"): after a few hours I would get an error like this: GEOS error: TopologyException: Input geom 0 is invalid: Ring Self-intersection at or near point 221912.50000000093 50580.000000001863 at 221912.50000000093 50580.000000001863 – I guess this will just be the first of many errors
  3. Cleaning the Shapefile, first with QGIS's Check Geometry Validity, then with GRASS's v.clean (I tried bpol), but the cleaned file still fails dissolve (I also tried adding a zero buffer)
  4. Converting the multipart polygons to singleparts, adding a geometry column and removing features smaller than a certain area. This took about 100MB off the filesize, but didn't affect the 'holes' within polygons – to get rid of those I tried to make a difference layer, but the difference operation consistently fails (if I check ignore invalid geometry, it creates lots of features which are visible in the attribute table, but doesn't display them).

I will simplify the file in the end, but don't want to do this before dissolving in case I introduce slivers.

Should I use a different tool within v.clean?

I have not yet tried to split the file into regions and dissolve those, then dissolve the regions.

Improve this question
4
  • I wonder what you want to achieve by dissolving those polygons? As far as I understand, this layer only gives you informations about where LIDAR data was taken when at what resolution etc. And of course the polygons overlap because the same area was overflown at different times.
    – Bernd V.
    Commented Jan 1, 2017 at 22:38
  • I suggest to simplify the geometries before dissolving them because otherwise you won't have neighbouring boundaries between the polygons. Furthermore, I would use a combination of Simplify geometries (first) and Smooth geometry (then) tools. If you are able to use scripts, also have a look here.
    – mgri
    Commented Jan 2, 2017 at 10:07
  • @BerndV. Yes, this data is just the extents of the LiDAR. The end-use is a web application, which needs to illustrate to the user which area is covered with the data.
    – user2950747
    Commented Jan 2, 2017 at 18:54
  • @HowToInQGIS Thanks for the helpful link. When I've used simplify in the past it's introduced slivers between polygons (even when those seem to share the same straight edge) so didn't want to do that before dissolving. I might be wrong on this though.
    – user2950747
    Commented Jan 2, 2017 at 18:55

1 Answer 1

Reset to default
18

The source data seems to be rather hard to handle as vectors as you have noticed. However, this workaroung that goes through an intermediate raster file works well and it is very fast.

1) Use gdal_rasterize http://www.gdal.org/gdal_rasterize.html and for making a black/white raster from the layer

gdal_rasterize -burn 100 -ts 1000 1000 -ot byte LIDAR_composite_extents_2015.shp burned.tif

2) Use gdal_polygonize http://www.gdal.org/gdal_polygonize.html for creating polygons from the raster

gdal_polygonize -f "ESRI Shapefile" burned.tif merged.shp

Select all the polygons which have "0" as the value of the "DN" attribute and delete them. As a result you will have a small 2 MB shapefile about the lidar coverage.

enter image description here

Because of rasterize-polygonize cycle some precision of the original vectors is lost. Increase the output raster size if you need better resolution. Depending on your use case you can perhaps use the raster image as an index without going back to vectors.

Improve this answer
2
  • This is a much more time-efficient method! The loss of data isn't a problem with a slightly higher resolution. I do need to go back to vectors (GeoJSON) at the end, but gdal_polygonize works fine for that. Thanks very much.
    – user2950747
    Commented Jan 2, 2017 at 18:52
  • There are some posts I wish I could upvote more. I had an ex-colleague tasked with doing exactly this -- calculate the extent of Lidar coverage in the UK -- and he used FME, which crashed after 4 weeks. No criticism of FME, just wrong tool for the job.
    – John Powell
    Commented Jan 16, 2017 at 18:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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