Make all layers broadcastable

SimpleSDMLayers/Project.toml

@@ -10,6 +10,7 @@ Distances = "b4f34e82-e78d-54a5-968a-f98e89d6e8f7"
 GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326"
 PolygonOps = "647866c9-e3ac-4575-94e7-e3d426903924"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
+Revise = "295af30f-e4ad-537b-8983-00126c2a3abe"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"

SimpleSDMLayers/src/SimpleSDMLayers.jl

@@ -13,14 +13,14 @@ using Tables
 # Basic types for the package
 include("lib/types.jl")
 export SimpleSDMLayer, SimpleSDMResponse, SimpleSDMPredictor
+export RasterCell # eltype for the layers, very useful to build table interface later
 
 # Main functions to match coordinates
 include("lib/coordinateconversion.jl")
 
 # Implements a series of interfaces (AbstractArray, iteration, and indexing)
 include("interfaces/common.jl")
 include("interfaces/iteration.jl")
-export RasterCell # eltype for the layers, very useful to build table interface later
 
 include("interfaces/indexing.jl")
 include("interfaces/broadcast.jl")

SimpleSDMLayers/src/interfaces/broadcast.jl

@@ -1,18 +1,70 @@
-Base.broadcastable(layer::SimpleSDMResponse) = layer
-Base.BroadcastStyle(::Type{<:SimpleSDMResponse}) = Broadcast.Style{SimpleSDMResponse}()
+Base.BroadcastStyle(::Type{T}) where {T <: SimpleSDMLayer} = Broadcast.Style{T}()
+Base.BroadcastStyle(
+    ::Broadcast.Style{T},
+    ::S,
+) where {T <: SimpleSDMLayer, S <: Broadcast.BroadcastStyle} = Broadcast.Style{T}()
+Base.BroadcastStyle(
+    ::S,
+    ::Broadcast.Style{T},
+) where {T <: SimpleSDMLayer, S <: Broadcast.BroadcastStyle} = Broadcast.Style{T}()
 
-Base.broadcastable(layer::SimpleSDMPredictor) = similar(layer)
-Base.BroadcastStyle(::Type{<:SimpleSDMPredictor}) = Broadcast.Style{SimpleSDMPredictor}()
+Base.axes(layer::T) where {T <: SimpleSDMLayer} = (Base.OneTo(length(layer)),)
+Base.broadcastable(layer::T) where {T <: SimpleSDMLayer} = layer
 
-function Base.Broadcast.broadcast(f, L::LT) where {LT <: SimpleSDMLayer}
-    newgrid = Array{Any}(nothing, size(L))
-    N = SimpleSDMResponse(newgrid, L)
-    v = filter(!isnothing, L.grid)
-    fv = f.(v)
-    N.grid[findall(!isnothing, L.grid)] .= fv
+find_layer(bc::Base.Broadcast.Broadcasted) = find_layer(bc.args)
+find_layer(args::Tuple) = find_layer(find_layer(args[1]), Base.tail(args))
+find_layer(x) = x
+find_layer(::Tuple{}) = nothing
+find_layer(layer::T, rest) where {T <: SimpleSDMLayer} = layer
+find_layer(::Any, rest) = find_aac(rest)
 
-    internal_types = unique(typeof.(N.grid))
+function Base.similar(
+    bc::Base.Broadcast.Broadcasted{Broadcast.Style{T}},
+    ::Type{ElType},
+) where {T <: SimpleSDMLayer, ElType}
+    @info "similar"
+    layer = find_layer(bc)
+    return similar(layer, ElType)
+end
 
-    RT = LT <: SimpleSDMResponse ? SimpleSDMResponse : SimpleSDMPredictor
-    return RT(convert(Matrix{Union{internal_types...}}, N.grid), N)
+function Base.broadcasted(::Broadcast.Style{T}, f, layer::T) where {T <: SimpleSDMLayer}
+    ElType = typeof(f(collect(layer)[1].value))
+    dest = similar(layer, ElType)
+    for cell in layer
+        dest[cell.longitude, cell.latitude] = f(cell.value)
+    end
+    return dest
+end
+
+function Base.broadcasted(::Broadcast.Style{T}, f, layer::T, x) where {T <: SimpleSDMLayer}
+    ElType = typeof(f(collect(layer)[1].value, x))
+    dest = similar(layer, ElType)
+    for cell in layer
+        dest[cell.longitude, cell.latitude] = f(cell.value, x)
+    end
+    return dest
+end
+
+function Base.broadcasted(::Broadcast.Style{T}, f, x, layer::T) where {T <: SimpleSDMLayer}
+    ElType = typeof(f(x, collect(layer)[1].value))
+    dest = similar(layer, ElType)
+    for cell in layer
+        dest[cell.longitude, cell.latitude] = f(x, cell.value)
+    end
+    return dest
+end
+
+function Base.broadcasted(
+    ::Broadcast.Style{T},
+    f,
+    layer1::T1,
+    layer2::T2,
+) where {T <: SimpleSDMLayer, T1 <: SimpleSDMLayer, T2 <: SimpleSDMLayer}
+    @assert SimpleSDMLayers._layers_are_compatible(layer1, layer2)
+    ElType = typeof(f(collect(layer1)[1].value, collect(layer2)[1].value))
+    dest = similar(layer1, ElType)
+    for cell in layer1
+        dest[cell.longitude, cell.latitude] = f(cell.value, layer2[cell.longitude, cell.latitude])
+    end
+    return dest
 end

SimpleSDMLayers/src/interfaces/iteration.jl

@@ -1,20 +1,5 @@
-struct RasterCell{L <: Number, T <: Any}
-    longitude::L
-    latitude::L
-    value::T
-end
-
-function RasterCell(layer::T, position) where {T <: SimpleSDMLayer}
-    lon = longitudes(layer)[last(position.I)]
-    lat = latitudes(layer)[first(position.I)]
-    val = layer.grid[position]
-    return RasterCell(lon, lat, val)
-end
-
 Base.IteratorSize(::T) where {T <: SimpleSDMLayer} = Base.HasLength()
-function Base.IteratorEltype(layer::T) where {T <: SimpleSDMLayer}
-    return RasterCell{eltype(latitudes(layer)), SimpleSDMLayers._inner_type(layer)}
-end
+Base.IteratorEltype(layer::T) where {T <: SimpleSDMLayer} = Base.HasEltype()
 
 function Base.iterate(layer::T) where {T <: SimpleSDMLayer}
     position = findfirst(!isnothing, layer.grid)

SimpleSDMLayers/src/lib/overloads.jl

@@ -22,7 +22,7 @@ import Base: findmax, findmin, findall
 Shows a textual representation of the layer.
 """
 function Base.show(io::IO, ::MIME"text/plain", layer::T) where {T <: SimpleSDMLayer}
-    itype = eltype(layer)
+    itype = SimpleSDMLayers._inner_type(layer)
     otype = T <: SimpleSDMPredictor ? "predictor" : "response"
     return print(
         io,
@@ -33,7 +33,7 @@ function Base.show(io::IO, ::MIME"text/plain", layer::T) where {T <: SimpleSDMLa
 end
 
 function Base.show(io::IO, layer::T) where {T <: SimpleSDMLayer}
-    itype = eltype(layer)
+    itype = SimpleSDMLayers._inner_type(layer)
     otype = T <: SimpleSDMPredictor ? "predictor" : "response"
     return print(
         io,
@@ -89,11 +89,13 @@ end
 """
     Base.eltype(layer::SimpleSDMLayer{T}) where {T}
 
-Returns the type of the values stored in the grid, where the `Nothing` type is
-omitted.
+Returns the type of the RasterCell describing each entry in the layer.
 """
-Base.eltype(::SimpleSDMResponse{T}) where {T} = T
-Base.eltype(::SimpleSDMPredictor{T}) where {T} = T
+function Base.eltype(layer::T) where {T <: SimpleSDMLayer}
+    coord_type = eltype(latitudes(layer))
+    data_type = SimpleSDMLayers._inner_type(layer)
+    return RasterCell{coord_type, data_type}
+end
 
 """
     Base.stride(layer::T; dims::Union{Nothing,Integer}=nothing) where {T <: SimpleSDMLayer}
@@ -146,7 +148,7 @@ zero for the type. If not, the same result can always be achieved through the
 use of `copy`, manual update, and `convert`.
 """
 function Base.similar(layer::T) where {T <: SimpleSDMLayer}
-    return similar(layer, eltype(layer))
+    return similar(layer, SimpleSDMLayers._inner_type(layer))
 end
 
 """

SimpleSDMLayers/src/lib/types.jl

@@ -143,3 +143,16 @@ for simplesdm_type in simplesdm_types
         end,
     )
 end
+
+struct RasterCell{L <: Number, T <: Any}
+    longitude::L
+    latitude::L
+    value::T
+end
+
+function RasterCell(layer::T, position) where {T <: SimpleSDMLayer}
+    lon = longitudes(layer)[last(position.I)]
+    lat = latitudes(layer)[first(position.I)]
+    val = layer.grid[position]
+    return RasterCell(lon, lat, val)
+end

SimpleSDMLayers/src/operations/mosaic.jl

@@ -21,7 +21,7 @@ function mosaic(f::TF, layers::Vector{T}) where {TF <: Function, T <: SimpleSDML
     end
 
     # Check the types
-    itypes = eltype.(layers)
+    itypes = SimpleSDMLayers._inner_type.(layers)
     if length(unique(itypes)) > 1
         @warn """
         The numeric types of the layers are not unique, this can cause performance issues.