Skip to content

Fix stride computation for dimensions with shape 0 in ndarray #68

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
cgarling wants to merge 5 commits into
base: main
Choose a base branch
from
Open

Fix stride computation for dimensions with shape 0 in ndarray #68

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
34f9a77
Fix stride computation for dimensions with shape 0 in ndarray
cgarling Aug 6, 2026
65d81f1
update `getindex` docstring
cgarling Aug 15, 2026
1c02d43
shorten test-ndarray.jl test comment
cgarling Aug 15, 2026
96de41c
add `!isempty(data)` check to `getindex`
cgarling Aug 15, 2026
3bda108
simplify zero-length dimension strides test
cgarling Aug 15, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
27 changes: 23 additions & 4 deletions src/ASDF.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -785,9 +785,21 @@ function NDArray(
offset = 0
end
if strides isa Nothing
# Calculate byte strides in C order
# Calculate byte strides in C order. Any dimension of length zero is treated as
# length 1 within this product only (not in `shape` itself); this matches NumPy's
# convention for computing default C-contiguous strides (`PyArray_NewFromDescr`),
# relied on by the reference Python `asdf` package when constructing arrays via
# `np.ndarray(shape, dtype, data, offset, None, order)`. Without the clamp, any
# zero-length dimension collapses the strides of every outer dimension whose
# product includes it down to zero, which then fails the `strides` positivity
# check below even though no data is ever read from a zero-size array. Negative
# entries are left unclamped so the `shape` negativity check below still reports
# them with its own clear error message. STScI Roman L2 `.asdf` products contain
# zero-shape arrays (e.g. `chisq`, `dumo`) with no explicit `strides` key,
# triggering this exact failure prior to the fix.
Comment on lines +788 to +799

@icweaver icweaver Aug 15, 2026 •
edited
Loading

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
# Calculate byte strides in C order. Any dimension of length zero is treated as
# length 1 within this product only (not in `shape` itself); this matches NumPy's
# convention for computing default C-contiguous strides (`PyArray_NewFromDescr`),
# relied on by the reference Python `asdf` package when constructing arrays via
# `np.ndarray(shape, dtype, data, offset, None, order)`. Without the clamp, any
# zero-length dimension collapses the strides of every outer dimension whose
# product includes it down to zero, which then fails the `strides` positivity
# check below even though no data is ever read from a zero-size array. Negative
# entries are left unclamped so the `shape` negativity check below still reports
# them with its own clear error message. STScI Roman L2 `.asdf` products contain
# zero-shape arrays (e.g. `chisq`, `dumo`) with no explicit `strides` key,
# triggering this exact failure prior to the fix.
# Calculate byte strides in C order, treating any zero-length dimension as
# length 1 within this product only (not in `shape` itself). This matches NumPy's
# convention for default C-contiguous strides (`PyArray_NewFromDescr`),
# relied on by the reference Python `asdf` package. Without the clamp, a
# zero-length dimension collapses the stride of every outer dimension whose
# product includes it down to zero, which then fails the `strides` positivity
# check below even though no data is ever read from a zero-size array. STScI Roman L2 `.asdf` products contain such
# zero-shape arrays (e.g. `chisq`, `dumo`) with no explicit `strides` key.

@cgarling cgarling Aug 15, 2026 •
edited
Loading

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would keep "Negative entries are left unclamped so the shape negativity check below still reports them with its own clear error message.", see comment below

sz = sizeof(Type(datatype))
strides = reverse(cumprod([sz; reverse(shape[(begin + 1):end])]))
clamped_shape = [s == 0 ? 1 : s for s in shape[(begin + 1):end]]
strides = reverse(cumprod([sz; reverse(clamped_shape)]))
Comment on lines +801 to +802

@icweaver icweaver Aug 15, 2026 •
edited
Loading

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
clamped_shape = [s == 0 ? 1 : s for s in shape[(begin + 1):end]]
strides = reverse(cumprod([sz; reverse(clamped_shape)]))
strides = reverse(cumprod([sz; reverse(max.(shape[(begin + 1):end], 1))]))

@icweaver icweaver Aug 15, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think this is type-stable and requires fewer changes to the original line here

@cgarling cgarling Aug 15, 2026 •
edited
Loading

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I tried this first but the max will also hit any negative entries in shape which I thought we should avoid, thus the comment above Negative entries are left unclamped so the shape negativity check below still reports them with its own clear error message.. If we use max negatives in shape will no longer error; lines 746-748

        if any(shape .< 0)
            throw(ArgumentError("`shape` cannot have negative elements."))
        end

@icweaver icweaver Aug 17, 2026

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure if I follow. Even when a negative shape entry is clamped within max, it is still passed unmodified to the inner constructor for NDArray afterwards, where the check you pointed to on L746-748 would throw as expected, no? That is what I am seeing on my end at least with the max.(...) version:

julia> let
           shape = Int64[-2, 3] # [2, -3], [-3], etc.
           strides = nothing
           ASDF.NDArray(ASDF.LazyBlockHeaders(), Int64(0), nothing, shape, ASDF.Datatype_float32, ASDF.host_byteorder, Int64(0), strides)
       end
ERROR: ArgumentError: `shape` cannot have negative elements.

end
return NDArray(
lazy_block_headers, source, data, Vector{Int64}(shape), datatype, byteorder, Int64(offset), Vector{Int64}(strides)
Expand Down Expand Up @@ -833,7 +845,7 @@ size(result) == Tuple(reverse(ndarray.shape))
eltype(result) == ASDF.materialized_eltype(ndarray.datatype)
```

For the `ucs4` and `ascii` string datatypes, [`materialized_eltype`](@ref) is a thin `AbstractString` view over the characters ([`UCS4String`](@ref) / [`AsciiString`](@ref); see [`stringify_data`](@ref)). For all other datatypes, `eltype(result) == Type(ndarray.datatype)` and additionally `sizeof(eltype) .* strides(result) == Tuple(reverse(ndarray.strides))`.
For the `ucs4` and `ascii` string datatypes, [`materialized_eltype`](@ref) is a thin `AbstractString` view over the characters ([`UCS4String`](@ref) / [`AsciiString`](@ref); see [`stringify_data`](@ref)). For all other datatypes, `eltype(result) == Type(ndarray.datatype)` and additionally, when the array has at least one element, `sizeof(eltype) .* strides(result) == Tuple(reverse(ndarray.strides))` along every dimension with more than one element. (A dimension of length 1 has no well-defined stride, there is no pair of adjacent elements to space apart along it, and in a zero-element array no dimension does, so `reshape`/`reinterpret` are free to report any value there; such strides are excluded from this check.)
"""
function Base.getindex(ndarray::NDArray)
if ndarray.data !== nothing
Expand Down Expand Up @@ -870,7 +882,14 @@ function Base.getindex(ndarray::NDArray)
# Check array layout
@assert size(data) == Tuple(reverse(ndarray.shape)) # `data` conforms to specified `ndarray.shape`
@assert eltype(data) == Type(ndarray.datatype) # `data` matches type specified by `ndarray.datatype`
if sizeof(eltype(data)) .* Base.strides(data) != Tuple(reverse(ndarray.strides))
# A dimension of length 1 has no meaningful stride (there is no pair of adjacent elements
# to space apart along it), and in an empty array no dimension does — Julia reports
# stride 0 along any dimension whose faster-varying dimensions include a zero length.
# Only compare strides where they are meaningful.
computed_strides = sizeof(eltype(data)) .* Base.strides(data)
expected_strides = Tuple(reverse(ndarray.strides))
data_shape = size(data)
if !isempty(data) && any(data_shape[i] > 1 && computed_strides[i] != expected_strides[i] for i in eachindex(data_shape))
error("`data` has different stride from `ndarray.strides`")
end

Expand Down
28 changes: 28 additions & 0 deletions test/test-ndarray.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -85,6 +85,34 @@ end
)
end

@testset "implicit strides with zero-length dimensions" begin
# 2026/08/06: Regression test for STScI Roman L2 `.asdf` products produced by romanisim,
# which contain `!core/ndarray` nodes (e.g. `chisq`/`dumo`) with a zero-length shape and
# no explicit `strides` key; see the default-strides comment in the `NDArray` outer
# constructor for the NumPy convention involved. Expected values below were verified
# against `np.ndarray(shape, dtype, buffer, offset, strides=None, order='C').strides`.
cases = [
(Int64[0, 0], ASDF.Datatype_float16, Int64[2, 2]),
(Int64[0, 5], ASDF.Datatype_float32, Int64[20, 4]),
(Int64[5, 0], ASDF.Datatype_float32, Int64[4, 4]),
(Int64[0, 0, 3], ASDF.Datatype_float32, Int64[12, 12, 4]),
(Int64[2, 0, 3], ASDF.Datatype_float32, Int64[12, 12, 4]),
(Int64[0, 2, 3], ASDF.Datatype_float32, Int64[24, 12, 4]),
]
for (shape, datatype, expected_strides) in cases
lbh = ASDF.LazyBlockHeaders()
push!(lbh.block_headers, make_block_header(UInt8[]))
nd = make_ndarray(;
lazy_block_headers = lbh, source = Int64(0), data = nothing, shape, datatype, strides = nothing,
)
@test nd.strides == expected_strides
# Materialize from an empty block to catch stride-check false positives at read time.
arr = nd[]
@test size(arr) == Tuple(reverse(shape))
@test eltype(arr) == Type(datatype)
end
end

@testset "getindex" begin
opposite = ASDF.host_byteorder == ASDF.Byteorder_little ? ASDF.Byteorder_big : ASDF.Byteorder_little

Expand Down
Loading