# How to clip a raster with a polygon

**URL:** <https://discourse.julialang.org/t/how-to-clip-a-raster-with-a-polygon/105883>\
**Category:** Geo\
**Tags:** question\
**Created:** [November 7, 2023, 2:05am UTC](https://discourse.julialang.org/t/how-to-clip-a-raster-with-a-polygon/105883 "2023-11-07T02:05:13Z")\
**Posts on this page:** 4\
**Page:** 1

<div class="post-metadata">

**Author:** ![WillHWThompson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/willhwthompson/32/44842_2.png) [@WillHWThompson](https://discourse.julialang.org/u/WillHWThompson)\
**Post date:** [November 7, 2023, 2:05am UTC](https://discourse.julialang.org/t/how-to-clip-a-raster-with-a-polygon/105883/1 "2023-11-07T02:05:13Z")

</div>

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?

---

<div class="post-metadata">

**Author:** ![joa-quim](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/joa-quim/32/227_2.png) [@joa-quim](https://discourse.julialang.org/u/joa-quim)\
**Post date:** [November 7, 2023, 3:21am UTC](https://discourse.julialang.org/t/how-to-clip-a-raster-with-a-polygon/105883/2 "2023-11-07T03:21:22Z")

</div>

I think the GMT.jl [rasterzones](https://www.generic-mapping-tools.org/GMTjl_doc/documentation/utilities/rasterzones/index.html#rasterzones) function can do what you want but you’ll need to convert your polygons into a vector of [GMTdatset](https://www.generic-mapping-tools.org/GMTjl_doc/documentation/general/types/index.html#dataset_type) (see function [mat2ds](https://www.generic-mapping-tools.org/GMTjl_doc/documentation/utilities/mat2ds/index.html#mat2ds)).

---

<div class="post-metadata">

**Author:** ![alex-s-gardner](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/alex-s-gardner/32/30210_2.png) [@alex-s-gardner](https://discourse.julialang.org/u/alex-s-gardner)\
**Post date:** [November 7, 2023, 4:03am UTC](https://discourse.julialang.org/t/how-to-clip-a-raster-with-a-polygon/105883/3 "2023-11-07T04:03:00Z")

</div>

The way I would approach this is to [rasterize](https://rafaqz.github.io/Rasters.jl/stable/reference/#Rasters.rasterize) the polygons. This will give you a Bool array that you can use to mask your raster

---

<div class="post-metadata">

**Author:** ![Raf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/raf/32/3383_2.png) [@Raf](https://discourse.julialang.org/u/Raf)\
**Post date:** [November 7, 2023, 8:21am UTC](https://discourse.julialang.org/t/how-to-clip-a-raster-with-a-polygon/105883/4 "2023-11-07T08:21:58Z")

</div>

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)