Design Review: CenterIndexedArrays.jl
Overview
CenterIndexedArrays provides a single concrete array type, CenterIndexedArray, that wraps any AbstractArray and reindexes it symmetrically around zero. The primary use case is image
registration, where displacement kernels are naturally expressed as signed offsets from a center point. The package also defines SymRange, a symmetric AbstractUnitRange, used internally
as the axis type but not exported.
The package is very small — ~110 lines of source across three files — so the findings below are correspondingly narrow.
Likely Design Issues
- SymRange(0) iteration is broken.
src/symrange.jl:18–21:
function iterate(r::SymRange)
r.n == 0 && return nothing
first(r), first(r)
end
SymRange(0) has first = last = 0 and length = 1. It represents the single-element range containing only 0. But the early-return guard makes it iterate as empty. The standard
AbstractUnitRange machinery (last(r) - first(r) + 1) correctly computes length 1, creating a contradiction: the range claims to have one element but iteration yields none.
SymRange(0) is created by Base.reduced_index (the last line of symrange.jl), so it appears as an axis whenever a dimension is reduced (e.g., minimum(A, dims=1)). Any code that iterates
over those axes — including eachindex, CartesianIndices, and for i in axes(result, 1) — will silently see zero iterations instead of one. The test suite doesn't cover iteration over
SymRange(0), so this has gone undetected.
The fix is a single line: remove the r.n == 0 && return nothing guard; the subsequent iterate(r, s) call (with s = first(r) = 0) will immediately hit s == last(r) and return nothing,
correctly yielding exactly one element.
Design Questions
- SymRange is not exported, but users unavoidably encounter it.
Every call to axes(A) or axes(A, d) returns a SymRange. Users who want to type-annotate, pattern-match, or construct one must write CenterIndexedArrays.SymRange. The README already shows
this, and tests pull it in explicitly with using CenterIndexedArrays: SymRange. This is effectively a public type that is treated as private.
Question for the author: Was the intent to keep SymRange hidden so the axis type is an implementation detail, or was it left unexported by omission? If users are expected to interact
with axes directly (which the README implies they are), exporting SymRange — or at minimum marking it public (Julia 1.11+) — would clarify the contract.
- Slicing with SymRange vs. UnitRange returns different types.
From the tests:
@test @inferred(A[:,SymRange(1)]) == CenterIndexedArray(dat[:,2:4]) # CenterIndexedArray
@test @inferred(A[:,-2:0]) == OffsetArray(dat[:,1:3], -1:1, 1:3) # OffsetArray
Slicing with a SymRange preserves the CenterIndexedArray wrapper; slicing with a UnitRange (even one derived from valid center indices) drops it and returns an OffsetArray. There is a
TODO comment in the source acknowledging this incompleteness (src/CenterIndexedArrays.jl:54):
This is incomplete: ideally we wouldn't need SymAx in the first slot
Question for the author: Is the asymmetry intentional — i.e., a UnitRange slice represents a non-centered sub-view that is semantically different from a centered one? Or is this a known
limitation waiting for the right approach? The current behavior can surprise users who take a slice to pass to a function and find the type has changed.
- Broadcasting returns CenterIndexedArray, but vec and cat return plain arrays.
A .+ 1 returns a CenterIndexedArray (the tests verify this). But vec(A) returns a plain Vector, and cat is commented out in the tests. This is a reasonable distinction — broadcasting is
shape-preserving while vec changes dimensionality — but it is implicit.
Question for the author: Is the intended rule "operations that preserve shape and size return CenterIndexedArray, others return plain arrays"? If so, what about reshape? Documenting this
policy (even briefly in the README) would help users predict what they'll get.
- The interpolation interface is implicit.
CenterIndexedArray{T,N,I<:AbstractInterpolation} has a specialized getindex that accepts Number arguments (allowing fractional indices, the whole point of interpolation). This is a
meaningful behavioral difference. There is no documentation of this capability in the README or docstring, and no dedicated constructor spelling out "here is how you wrap an
interpolation."
Question for the author: Should the interpolation use case be documented as a first-class feature? And is the implicit type-parameter dispatch (A<:AbstractInterpolation) the right
mechanism, or would a wrapper type (CenterInterpolation?) be clearer for users discovering the API?
Observations
- Deprecations are present in deprecated.jl (the old CenterIndexedArray(::Type{T}, dims) constructors). These are handled by the freshen-package workflow and are noted here for
completeness.
- No public annotations anywhere. Julia 1.11 introduced public as a way to mark non-exported names as part of the public API. SymRange is a candidate.
- _halfsize error message uses bare error("Must have all-odd sizes") rather than a typed exception. This is fine for a small utility, but it means callers can't catch this condition
specifically.
- Base.throw_boundserror is overridden (src/CenterIndexedArrays.jl:83) to add @_noinline_meta. This is a micro-optimization that is reasonable but undocumented. Worth keeping an eye on
if Base changes its internal API.
- The copy keyword in iterate(r, s) = copy(s+1), s+1 is interesting — s is an Int, so copy is a no-op. This appears to be a leftover defensive pattern. It does no harm but is slightly
surprising to a reader.
Overall Assessment
The package is coherent and appropriately narrow. Its identity is clear, its type hierarchy is shallow (one type, one axis type), and there are no leaking abstractions or off-topic
exports. The main tension is between internal simplicity and user ergonomics: SymRange is effectively public but not declared as such, the slicing behavior is inconsistent, and the
interpolation capability is invisible to new users.
The highest-leverage change would be fixing the SymRange(0) iteration bug, which is a correctness issue with no downside to fixing. The second-highest would be exporting or marking
public the SymRange type, which closes the gap between what the API promises and what users must write to work with it.
Design Review: CenterIndexedArrays.jl
Overview
CenterIndexedArrays provides a single concrete array type, CenterIndexedArray, that wraps any AbstractArray and reindexes it symmetrically around zero. The primary use case is image
registration, where displacement kernels are naturally expressed as signed offsets from a center point. The package also defines SymRange, a symmetric AbstractUnitRange, used internally
as the axis type but not exported.
The package is very small — ~110 lines of source across three files — so the findings below are correspondingly narrow.
Likely Design Issues
src/symrange.jl:18–21:
function iterate(r::SymRange)
r.n == 0 && return nothing
first(r), first(r)
end
SymRange(0) has first = last = 0 and length = 1. It represents the single-element range containing only 0. But the early-return guard makes it iterate as empty. The standard
AbstractUnitRange machinery (last(r) - first(r) + 1) correctly computes length 1, creating a contradiction: the range claims to have one element but iteration yields none.
SymRange(0) is created by Base.reduced_index (the last line of symrange.jl), so it appears as an axis whenever a dimension is reduced (e.g., minimum(A, dims=1)). Any code that iterates
over those axes — including eachindex, CartesianIndices, and for i in axes(result, 1) — will silently see zero iterations instead of one. The test suite doesn't cover iteration over
SymRange(0), so this has gone undetected.
The fix is a single line: remove the r.n == 0 && return nothing guard; the subsequent iterate(r, s) call (with s = first(r) = 0) will immediately hit s == last(r) and return nothing,
correctly yielding exactly one element.
Design Questions
Every call to axes(A) or axes(A, d) returns a SymRange. Users who want to type-annotate, pattern-match, or construct one must write CenterIndexedArrays.SymRange. The README already shows
this, and tests pull it in explicitly with using CenterIndexedArrays: SymRange. This is effectively a public type that is treated as private.
Question for the author: Was the intent to keep SymRange hidden so the axis type is an implementation detail, or was it left unexported by omission? If users are expected to interact
with axes directly (which the README implies they are), exporting SymRange — or at minimum marking it public (Julia 1.11+) — would clarify the contract.
From the tests:
@test @inferred(A[:,SymRange(1)]) == CenterIndexedArray(dat[:,2:4]) # CenterIndexedArray
@test @inferred(A[:,-2:0]) == OffsetArray(dat[:,1:3], -1:1, 1:3) # OffsetArray
Slicing with a SymRange preserves the CenterIndexedArray wrapper; slicing with a UnitRange (even one derived from valid center indices) drops it and returns an OffsetArray. There is a
TODO comment in the source acknowledging this incompleteness (src/CenterIndexedArrays.jl:54):
This is incomplete: ideally we wouldn't need SymAx in the first slot
Question for the author: Is the asymmetry intentional — i.e., a UnitRange slice represents a non-centered sub-view that is semantically different from a centered one? Or is this a known
limitation waiting for the right approach? The current behavior can surprise users who take a slice to pass to a function and find the type has changed.
A .+ 1 returns a CenterIndexedArray (the tests verify this). But vec(A) returns a plain Vector, and cat is commented out in the tests. This is a reasonable distinction — broadcasting is
shape-preserving while vec changes dimensionality — but it is implicit.
Question for the author: Is the intended rule "operations that preserve shape and size return CenterIndexedArray, others return plain arrays"? If so, what about reshape? Documenting this
policy (even briefly in the README) would help users predict what they'll get.
CenterIndexedArray{T,N,I<:AbstractInterpolation} has a specialized getindex that accepts Number arguments (allowing fractional indices, the whole point of interpolation). This is a
meaningful behavioral difference. There is no documentation of this capability in the README or docstring, and no dedicated constructor spelling out "here is how you wrap an
interpolation."
Question for the author: Should the interpolation use case be documented as a first-class feature? And is the implicit type-parameter dispatch (A<:AbstractInterpolation) the right
mechanism, or would a wrapper type (CenterInterpolation?) be clearer for users discovering the API?
Observations
completeness.
specifically.
if Base changes its internal API.
surprising to a reader.
Overall Assessment
The package is coherent and appropriately narrow. Its identity is clear, its type hierarchy is shallow (one type, one axis type), and there are no leaking abstractions or off-topic
exports. The main tension is between internal simplicity and user ergonomics: SymRange is effectively public but not declared as such, the slicing behavior is inconsistent, and the
interpolation capability is invisible to new users.
The highest-leverage change would be fixing the SymRange(0) iteration bug, which is a correctness issue with no downside to fixing. The second-highest would be exporting or marking
public the SymRange type, which closes the gap between what the API promises and what users must write to work with it.