I have thousands of shapefiles covering around the world. I need to calculate coverage percentage. For example, a pixel located in Hawaii was covered by 100 of all shapefiles, then this pixel value of 100 will be output.

I can convert them to rasters, but I cannot stack all rasters in a file for calculation because of their large volume.

So, how can I calculate coverage percentage by using shapefiles or separated rasters?

closed as off-topic by ahmadhanb, Kadir Şahbaz, xunilk, LaughU, whyzar Jul 8 at 13:17

This question appears to be off-topic. The users who voted to close gave this specific reason:

  • "When seeking help to debug/write/improve code always provide the desired behavior, a specific problem/error and the shortest code (as formatted text, not pictures) needed to reproduce it in the question body. Providing a clear problem statement and a code attempt helps others to help you." – ahmadhanb, Kadir Şahbaz, xunilk, LaughU, whyzar
If this question can be reworded to fit the rules in the help center, please edit the question.


A conceptual solution would be to reclass all your rasters with 1 (covered) and 0 (not covered) and then use raster algebra to add all your rasters together. The pixel value will be the total value of all rasters which cover the given area. Take extra care to deal with different raster resolutions.


You wouldn't need to convert all the shapefiles to raster at once, so volume consideration is a non-issue.

Using a raster masking technique is going to result in imprecision at the pixel boundary -- Does a file which include the center of the pixel count as "covered" if not 100% of the pixel is covered? Using a vector technique could attain higher precision, at the cost of a lot of processing.

In raster space, you just need to mask each vector at a time (with consideration as to whether you should take a one pixel (or sqrt(2)*pixeldim) negative buffer before masking), sum the presence/absence flag with a cumulator, and then divide by the number of shapefiles (encoding this as byte (either 0-100 or 0-250) would speed future processing, at the cost of precision)

In vector space, you'd need to:

  • Determine the resolution of your percentage raster (pixel size, type and depth)
  • Create a partition vector to form the basis for subsetting the vector space which is in alignment with the raster (large enough to reduce the iterations, but small enough to save on vector overlay cost)
  • Iterate over the tiles with a clip and union of each shapefile, then compute percentage and convert that to raster, encode the desired pixel value, and mask this into a final raster.

For most approaches to either technique, you're going to need a full GIS, not just GDAL and Python.

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