How to clip a raster with a polygon

Specific Domains Geo
1

Hi all. I have a julia GIS problem that I can’t seem to solve.

I have been using Rasters.jl to load rasterized GIS data into julia. This implements a Raster data structure where latitude and longitude points are associated with a scalar value.

I have a Voronoi tesselation over the same(rectangular) region that the raster is defined over. I want to find the sum of the values of the raster that lie within each Voronoi cell. The Voronoi tessellation is stored as a list of polygons. I want to find the portion of each raster which falls within each polygon.

I have tried using the clip() function from Rasters.jl but this uses the minimum bounding rectangle around a polygon, not the geometry of the polygon itself.

Is there any existing library that can perform something like this?

2

I think the GMT.jl rasterzones function can do what you want but you’ll need to convert your polygons into a vector of GMTdatset (see function mat2ds).

3

The way I would approach this is to rasterize the polygons. This will give you a Bool array that you can use to mask your raster

4

Hi @WillHWThompson

Yes crop is about setting the bounding box, not the values. With Rasters.jl we have mask to actually mask a raster values with geometries. See the doc string for usage. It makes sense to crop and then mask.

But also, your problem appears to be zonal statistics, and for that we have zonal that just does the sum directly. So zonal(sum, raster; of=polygons) where polygons is a geometry or vector of geometries/features.

(@alex-s-gardner mask is easier to call than rasterize and can be faster as we don’t care about the reducing function, fill value or any of those details, we just flatten everything to bools)

2 Likes