add BorderArray

.gitignore

@@ -3,3 +3,4 @@
 *.jl.mem
 docs/build
 docs/site
+Manifest.toml

docs/src/function_reference.md

@@ -53,6 +53,7 @@ reflect
 
 ```@docs
 padarray
+BorderArray
 Pad
 Fill
 Inner

src/ImageFiltering.jl

@@ -8,7 +8,10 @@ using Base: Indices, tail, fill_to_length, @pure, depwarn, @propagate_inbounds
 using OffsetArrays: IdentityUnitRange   # using the one in OffsetArrays makes this work with multiple Julia versions
 using Requires
 
-export Kernel, KernelFactors, Pad, Fill, Inner, NA, NoPad, Algorithm,
+export Kernel, KernelFactors,
+    Pad, Fill, Inner, NA, NoPad,
+    BorderArray,
+    Algorithm,
     imfilter, imfilter!,
     mapwindow, mapwindow!,
     imgradients, padarray, centered, kernelfactors, reflect
@@ -65,6 +68,7 @@ using .Kernel: Laplacian, reflect, ando3, ando4, ando5, scharr, bickley, prewitt
 NDimKernel{N,K} = Union{AbstractArray{K,N},ReshapedOneD{K,N},Laplacian{N}}
 
 include("border.jl")
+include("borderarray.jl")
 
 BorderSpec{T} = Union{Pad{0}, Fill{T,0}, Inner{0}}
 BorderSpecNoNa{T} = Union{Pad{0}, Fill{T,0}, Inner{0}}

src/border.jl

@@ -236,8 +236,12 @@ function padfft(indk::AbstractUnitRange, l::Integer)
     range(first(indk), length=nextprod([2,3], l+lk)-l+1)
 end
 
+function error_bad_padding_size(inner, border)
+    ArgumentError("$border lacks the proper padding sizes for an array with $(ndims(inner)) dimensions")
+end
+
 function padindices(img::AbstractArray{<:Any,N}, border::Pad) where N
-    throw(ArgumentError("$border lacks the proper padding sizes for an array with $(ndims(img)) dimensions"))
+    throw(error_bad_padding_size(img, border))
 end
 function padindices(img::AbstractArray{<:Any,N}, border::Pad{N}) where N
     _padindices(border, border.lo, axes(img), border.hi)
@@ -652,13 +656,11 @@ The command `padarray(A,Fill(0,(1,1,1)))` yields
 ---
 
 """
-padarray(img::AbstractArray, border::Pad)  = padarray(eltype(img), img, border)
-function padarray(::Type{T}, img::AbstractArray, border::Pad) where T
-    inds = padindices(img, border)
-    # like img[inds...] except that we can control the element type
-    newinds = map(Base.axes1, inds)
-    dest = similar(img, T, newinds)
-    copydata!(dest, img, inds)
+padarray(img::AbstractArray, border::AbstractBorder)  = padarray(eltype(img), img, border)
+function padarray(::Type{T}, img::AbstractArray, border) where {T}
+    ba = BorderArray(img, border)
+    out = similar(ba, T)
+    copy!(out, ba)
 end
 
 padarray(img, ::Type{P}) where {P} = img[padindices(img, P)...]      # just to throw the nice error
@@ -684,7 +686,8 @@ function copydata!(dest::OffsetArray, img, inds::Tuple{Vararg{OffsetArray}})
     dest
 end
 
-Base.ndims(::Pad{N}) where {N} = N
+Base.ndims(b::AbstractBorder) = ndims(typeof(b))
+Base.ndims(::Type{Pad{N}}) where {N} = N
 
 # Make these separate types because the dispatch almost surely needs to be different
 """
@@ -702,7 +705,6 @@ struct Inner{N} <: AbstractBorder
     hi::Dims{N}
 end
 
-Base.ndims(::Inner{N}) where N = N
 Base.ndims(::Type{Inner{N}}) where N = N
 
 """
@@ -743,10 +745,6 @@ for T in (:Inner, :NA)
     end
 end
 
-padarray(img, border::Inner) = padarray(eltype(img), img, border)
-padarray(::Type{T}, img::AbstractArray{T}, border::Inner) where {T} = copy(img)
-padarray(::Type{T}, img::AbstractArray, border::Inner) where {T} = copyto!(similar(Array{T}, axes(img)), img)
-
 """
 ```julia
     struct Fill{T,N} <: AbstractBorder
@@ -806,6 +804,7 @@ struct Fill{T,N} <: AbstractBorder
     Fill{T,N}(value::T) where {T,N} = new{T,N}(value)
     Fill{T,N}(value::T, lo::Dims{N}, hi::Dims{N}) where {T,N} = new{T,N}(value, lo, hi)
 end
+Base.ndims(::Type{Fill{T,N}}) where {T,N} = N
 
 """
 ```julia
@@ -895,22 +894,6 @@ function (p::Fill)(kernel, img, ::FFT)
     Fill(p.value, newinds)
 end
 
-function padarray(::Type{T}, img::AbstractArray, border::Fill) where T
-    throw(ArgumentError("$border lacks the proper padding sizes for an array with $(ndims(img)) dimensions"))
-end
-function padarray(::Type{T}, img::AbstractArray{<:Any,N}, f::Fill{<:Any,N}) where {T,N}
-    paxs = map((l,r,h)->first(r)-l:last(r)+h, f.lo, axes(img), f.hi)
-    A = similar(arraytype(img, T), paxs)
-    try
-        fill!(A, f.value)
-    catch
-        error("Unable to fill! an array of element type $(eltype(A)) with the value $(f.value). Supply an appropriate value to `Fill`, such as `zero(eltype(A))`.")
-    end
-    A[axes(img)...] = img
-    A
-end
-padarray(img::AbstractArray, f::Fill) = padarray(eltype(img), img, f)
-
 # There are other ways to define these, but using `mod` makes it safe
 # for cases where the padding is bigger than length(inds)
 """

src/borderarray.jl

@@ -0,0 +1,186 @@
+using Base: @propagate_inbounds
+
+function compatible_dimensions(arr::AbstractArray, border::Inner)
+    true
+end
+
+function compatible_dimensions(arr,border)
+    ndims(arr) == ndims(border)
+end
+
+function convert_border_eltype(inner, border::Fill)
+    T = eltype(inner)
+    val::T = try
+        convert(T, border.value)
+    catch err
+        msg = "Cannot convert elements of border=$border to eltype(inner)=$T."
+        throw(ArgumentError(msg))
+    end
+    Fill(val, border.lo, border.hi)
+end
+
+function convert_border_eltype(inner, border)
+    border
+end
+
+"""
+    BorderArray(inner::AbstractArray, border::AbstractBorder) <: AbstractArray
+
+Construct a thin wrapper around the array `inner`, with given `border`. No data is copied in the constructor, instead border values are computed on the fly in `getindex` calls. Usful for stencil computations. See also [padarray](@ref).
+
+# Examples
+```julia
+julia> using ImageFiltering
+
+julia> arr = reshape(1:6, (2,3))
+2×3 reshape(::UnitRange{Int64}, 2, 3) with eltype Int64:
+ 1  3  5
+ 2  4  6
+
+julia> BorderArray(arr, Pad((1,1)))
+BorderArray{Int64,2,Base.ReshapedArray{Int64,2,UnitRange{Int64},Tuple{}},Pad{2}} with indices 0:3×0:4:
+ 1  1  3  5  5
+ 1  1  3  5  5
+ 2  2  4  6  6
+ 2  2  4  6  6
+
+julia> BorderArray(arr, Fill(10, (2,1)))
+BorderArray{Int64,2,Base.ReshapedArray{Int64,2,UnitRange{Int64},Tuple{}},Fill{Int64,2}} with indices -1:4×0:4:
+ 10  10  10  10  10
+ 10  10  10  10  10
+ 10   1   3   5  10
+ 10   2   4   6  10
+ 10  10  10  10  10
+ 10  10  10  10  10
+```
+"""
+struct BorderArray{T,N,A,B} <: AbstractArray{T,N}
+    inner::A
+    border::B
+    function BorderArray(inner::AbstractArray{T,N}, border::AbstractBorder) where {T,N}
+        if !compatible_dimensions(inner, border)
+            throw(error_bad_padding_size(inner, border))
+        end
+        border = convert_border_eltype(inner, border)
+        A = typeof(inner)
+        B = typeof(border)
+        new{T,N,A,B}(inner, border)
+    end
+end
+
+# these are adapted from OffsetArrays
+function Base.similar(A::BorderArray, ::Type{T}, dims::Dims) where T
+    B = similar(A.inner, T, dims)
+end
+const OffsetAxis = Union{Integer, UnitRange, Base.OneTo, IdentityUnitRange}
+function Base.similar(A::BorderArray, ::Type{T}, inds::Tuple{OffsetAxis,Vararg{OffsetAxis}}) where T
+    similar(A.inner, T, axes(A))
+end
+
+@inline function Base.axes(o::BorderArray)
+    _outeraxes(o.inner, o.border)
+end
+
+@inline function Base.size(o::BorderArray)
+    map(length, axes(o))
+end
+
+@inline function _outeraxes(arr, border::Inner)
+    axes(arr)
+end
+
+@inline function _outeraxes(arr, border)
+    map(axes(arr), border.lo, border.hi) do r, lo, hi
+        (first(r)-lo):(last(r)+hi)
+    end
+end
+
+@inline @propagate_inbounds function Base.getindex(o::BorderArray{T,N}, inds::Vararg{Int,N}) where {T,N}
+    ci = CartesianIndex(inds)
+    o[ci]
+end
+
+@inline function Base.getindex(o::BorderArray, ci::CartesianIndex)
+    if checkbounds(Bool, o.inner, ci)
+        @inbounds o.inner[ci]
+    else
+        @boundscheck checkbounds(o, ci)
+        getindex_outer_inbounds(o.inner, o.border, ci)
+    end
+end
+
+function getindex_outer_inbounds(arr, b::Fill, index)
+    b.value
+end
+
+function _inner_index(arr, b::Pad, index::CartesianIndex)
+    s = b.style
+    inds = if s == :replicate
+        map(index.I, axes(arr)) do i, r
+            clamp(i, first(r), last(r))
+        end
+    elseif s == :circular
+        map(modrange, index.I, axes(arr))
+    elseif s == :reflect
+        map(index.I, axes(arr)) do i, r
+            if i > last(r)
+                i = 2last(r) - i
+            elseif i < first(r)
+                i = 2first(r) - i
+            end
+            i
+        end
+    elseif s == :symmetric
+        map(index.I, axes(arr)) do i, r
+            if i > last(r)
+                i = 2last(r) + 1 - i
+            elseif i < first(r)
+                i = 2first(r) -1 - i
+            end
+            i
+        end
+    else
+        error("border style $(s) unrecognized")
+    end
+    CartesianIndex(inds)
+end
+
+function getindex_outer_inbounds(arr, b::Pad, index::CartesianIndex)
+    i = _inner_index(arr, b, index)
+    @inbounds arr[i]
+end
+
+function Base.checkbounds(::Type{Bool}, arr::BorderArray, index::CartesianIndex)
+    map(axes(arr), index.I) do r, i
+        checkindex(Bool, r, i)
+    end |> all
+end
+
+function Base.copy!(dst::AbstractArray, src::BorderArray)
+    axes(dst) == axes(src) || throw(DimensionMismatch("axes(dst) == axes(src) must hold."))
+    _copy!(dst, src, src.border)
+end
+function Base.copy!(dst::AbstractArray{T,1} where T, src::BorderArray{T,1,A,B} where B where A where T)
+    # fix ambiguity
+    axes(dst) == axes(src) || throw(DimensionMismatch("axes(dst) == axes(src) must hold."))
+    _copy!(dst, src, src.border)
+end
+
+function _copy!(dst, src, ::Inner)
+    copyto!(dst, src.inner)
+end
+
+function _copy!(dst, src, ::Pad)
+    inds = padindices(src.inner, src.border)
+    copydata!(dst, src.inner, inds)
+end
+
+function _copy!(dst, src, ::Fill)
+    try
+        fill!(dst, src.border.value)
+    catch
+        error("Unable to fill! an array of element type $(eltype(dst)) with the value $(src.border.value). Supply an appropriate value to `Fill`, such as `zero(eltype(A))`.")
+    end
+    dst[axes(src.inner)...] = src.inner
+    dst
+end

test/border.jl

@@ -2,6 +2,34 @@ using ImageFiltering, OffsetArrays, Colors, FixedPointNumbers, Random
 using Test
 
 @testset "Border" begin
+    borders = []
+    push!(borders, Inner())
+    push!(borders, Fill(0, (1,2), (3,4)))
+    for _ in 1:10
+        s = rand([:replicate, :circular, :symmetric, :reflect])
+        pad = Pad(s, (rand(0:4),rand(0:4)), (rand(0:4), rand(0:4)))
+        push!(borders, pad)
+    end
+    @testset "BorderArray consistent with padarray $b" for b in borders
+        A = randn(5,5)
+        B = BorderArray(A, b)
+        @test padarray(A, b) == B
+        @inferred B[1,1]
+        for I in CartesianIndices(A)
+            @test A[I] == B[I]
+        end
+    end
+    @testset "BorderArray" begin
+        @test_throws ArgumentError BorderArray(randn(3), Fill(0.))
+        @test_throws ArgumentError BorderArray(randn(3), Fill(nothing, (1,)))
+        @inferred BorderArray(randn(3), Fill(0., (1,)))
+        @inferred BorderArray(randn(3), Fill(0, (1,)))
+        ba1 = BorderArray(randn(3), Fill(zero(Int    ), (1,)))
+        ba2 = BorderArray(randn(3), Fill(zero(Float64), (1,)))
+        @test typeof(ba1) === typeof(ba2)
+        @inferred ba1[1]
+    end
+
     @testset "padarray" begin
         A = reshape(1:25, 5, 5)
         @test @inferred(padarray(A, Fill(0,(2,2),(2,2)))) == OffsetArray(
@@ -190,14 +218,13 @@ using Test
         B[1,1] = 0
         @test B != A
         A = rand(RGB{N0f8}, 3, 5)
-        ret = @test_throws ErrorException padarray(A, Fill(0, (0,0), (0,0)))
-        # FIXME: exact phrase depends on showarg in Interpolations
-        # @test occursin("element type ColorTypes.RGB", ret.value.msg)
-        # This is a temporary substitute:
+        ret = @test_throws ArgumentError padarray(A, Fill(0, (0,0), (0,0)))
         @test occursin("RGB", ret.value.msg)
+        @test occursin("convert", ret.value.msg)
         A = bitrand(3, 5)
-        ret = @test_throws ErrorException padarray(A, Fill(7, (0,0), (0,0)))
-        @test occursin("element type Bool", ret.value.msg)
+        ret = @test_throws ArgumentError padarray(A, Fill(7, (0,0), (0,0)))
+        @test occursin("Bool", ret.value.msg)
+        @test occursin("convert", ret.value.msg)
         @test isa(parent(padarray(A, Fill(false, (1,1), (1,1)))), BitArray)
     end