Rasters.jl icon indicating copy to clipboard operation
Rasters.jl copied to clipboard

[WIP] Add a `bylayer` keyword to `zonal`

Open asinghvi17 opened this issue 1 year ago • 5 comments

Successor to #775. This provides a switch that allows you to iterate directly over the values of a rasterstack, instead of iterating over the layers and then the layer values.

It currently defaults to true (the old behaviour) but setting it to false allows iteration over the NamedTuple values of each dimselector in the stack, which is the new behaviour.

Needs tests. This also fails when skipmissing=true, because it's not implemented for RasterStacks. From discussion in #775 we still need to decide whether skipmissing on a rasterstack should skip a cell if it has any missing values, or skip a cell if all its values are missing. Personally, I'm inclined to option 1 :).

asinghvi17 avatar Nov 12 '24 16:11 asinghvi17

Yeah, lets not break things straight away. The main reason to have it as true is so sum and mean just work as-is.

I think we probably need to move away from returning NamedTuple at some point to a similar object that math works on. At that point most workflows wont break when we flip the switch, but until then bylayer=true needs to be the default.

I also think option 1 is the correct skipmissing behaviour. We could also allow passing a Symbol like skipmissing=:somelayer to specify the layer to use.

rafaqz avatar Nov 12 '24 16:11 rafaqz

Here's what the namedtuple warning looks like:

julia> zonal(mean, st; of = [polygon], bylayer=false)
┌ Warning: The function you passed to `zonal` cannot handle a RasterStack's values directly,
│ as they are NamedTuples of the form `(; varname = value, ...)`.
│ 
│ Set `bylayer=true` to apply the function to each layer of the stack separately, 
│ or change the function to accept NamedTuple input.
└ @ Rasters ~/.julia/dev/Rasters/src/methods/zonal.jl:170
ERROR: MethodError: no method matching /(::@NamedTuple{a::Missing, b::Missing}, ::Int64)
The function `/` exists, but no method is defined for this combination of argument types.

Closest candidates are:
  /(::BigInt, ::Union{Int16, Int32, Int64, Int8, UInt16, UInt32, UInt64, UInt8})
   @ Base gmp.jl:560
  /(::BigFloat, ::Union{Int16, Int32, Int64, Int8})
   @ Base mpfr.jl:553
  /(::Missing, ::Number)
   @ Base missing.jl:123
  ...

Stacktrace:
  [1] mean(f::typeof(identity), itr::Base.SkipMissing{RasterStack{…}})
    @ Statistics ~/.julia/packages/Statistics/gbcbG/src/Statistics.jl:69
  [2] mean(itr::Base.SkipMissing{RasterStack{…}})
    @ Statistics ~/.julia/packages/Statistics/gbcbG/src/Statistics.jl:44
  [3] _maybe_skipmissing_call(f::typeof(mean), A::RasterStack{…}, sm::Bool)
    @ Rasters ~/.julia/dev/Rasters/src/methods/zonal.jl:206

asinghvi17 avatar Nov 13 '24 16:11 asinghvi17

It's pretty easy to miss in the giant error stack though. This also doesn't activate in zonal where of is a single geometry, only where it is a vector of geometries.

asinghvi17 avatar Nov 13 '24 16:11 asinghvi17

Was there any objection to this? Or should I fix it up?

asinghvi17 avatar Apr 22 '25 19:04 asinghvi17

No it's great we just forgot about it

rafaqz avatar Apr 22 '25 20:04 rafaqz