How do I use rasterio/python to mask a raster using a shapefile, to set the raster pixels inside the polygons to zero?

Question

I am trying to create a land mask to apply to satellite imagery, that will set the pixels in a raster intersecting with a land mass to 0.

After experimenting with gdal, skimage, pyplot etc. I've found the method given in the rasterio cookbook to be quick and easy. However, it is setting the pixels outside the polygons to 0, whereas I am trying to do the inverse of this.

Keep to using rasterio if possible - you don't have to calculate pixel locations of geospatial coordinates or deal with clipping features that are beyond the extents of the raster turning negative. It's also FAST which is important for the file size of the raw imagery I'm working with.

From here: https://mapbox.s3.amazonaws.com/playground/perrygeo/rasterio-docs/cookbook.html#masking-raster-with-a-polygon-feature

My code is as follows:

import fiona
import rasterio
from rasterio.tools.mask import mask

with fiona.open("/Users/Cate/UK_Mainland.shp", "r") as shapefile:
    geoms = [feature["geometry"] for feature in shapefile]

with rasterio.open("jan_clip.tif") as src:
    out_image, out_transform = mask(src, geoms, crop=True)
    out_meta = src.meta.copy()

out_meta.update({"driver": "GTiff",
                 "height": out_image.shape[1],
                 "width": out_image.shape[2],
                 "transform": out_transform})

with rasterio.open("masked2.tif", "w", **out_meta) as dest:
    dest.write(out_image)

How do I mask the areas that intersect with the polygons rather than those that don't?

Answer 1

rasterio.tools.mask.mask (in more recent versions, it is rasterio.mask.mask) includes an option invert. When invert=True, the mask will be applied to pixels that overlap your shape, rather than areas outside the shape. So you can change the line above to:

out_image, out_transform = mask(src, geoms, invert=True)