Incorrectly masked result with an out-of-bounds polygon mask

[sidpku]

I tried to mask a raster with a polygon that extends across the boundary of the raster. However, all pixels are masked (set as missing). I built this demo to replicate the error. The latest version of Rasters.jl still has this problem.

using Rasters, Plots, Shapefile
import ArchGDAL

# the raster to be masked
R = Raster(rand(100,100), (Rasters.X(1:100),Rasters.Y(1:100)))

# the geometry to be used as a mask
roi = Shapefile.Polygon(
    Shapefile.Rect(10.0, -10.0, 50.0, 50.0),
    [0],
    [Shapefile.Point(10.0,50.0),
     Shapefile.Point(50.0,50.0),
     Shapefile.Point(50.0,-10.0),
     Shapefile.Point(10.0, -10.0),
     Shapefile.Point(10.0,50.0)
    ],
    []
)

# plot the raster and the geometry(polygon)
plot(R,xlims=(0,100), ylims = (-15,100))
plot!(roi, fill = nothing, linewidth = 2)

The following figure displays the Raster and the polygon mask. The Raster lies above the 0 latitude, whereas the polygon spans the equator and stretches to 10 S.

Rmasked = mask(R, with = roi)
plot(Rmasked)

using above codes, we want to get a masked Raster (I added the polygon to show the mask)

However, we got the wrong results (fig.3), where all pixels are masked.

possilbe reason might be relevant to boolmask!(dest::AbstractRaster, geoms;kw...) or the function burn_geometry! nested in the boolmask!.

PS: The CLI output of Raster in the latest version looks great! I'm really impressed by it! This is my first time submitting an issue, so I apologize for any mistakes or my poor English.

[rafaqz]

May be a bug! I will have a look.

But note that hand-rolling Shapefile.jl geometries is not the best test here, they are not made to be manually defined like that. It's just a format to read esri shape files. I'm not sure that is correct syntax and I wrote 1/3 of the code there...

Better to use GeometryBasics.jl or GeoInterface.jl geometries, although probably that is not made clear anywhere.

[rafaqz]

Ok yes there is a serious bug with polygon rasterization in general for this example (its not just the Shapefile and not just mask).

I'm not sure how this got through testing and the other polygons are passing. Will try to fix it over the weekend.