From bba156041665355611904dcc64efbbd33520195a Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:04:19 -0500 Subject: [PATCH 1/7] fix(format): lay Fixed/Extensible Array chunk indexes out by max dims Both indexes place each chunk at a linear index computed from the dataset's maximum dimensions (libhdf5's max_down_chunks), and the Extensible Array first swizzles its unlimited dimension to the slowest position. We linearised by the current dimensions, so any dataset whose shape was smaller than its maxshape, or whose unlimited dimension was not the first, read back scrambled without an error: h5py libver="latest" files with maxshape (10, None) or (20, 10), and the libhdf5 test files h5fc_ext*.h5 and test_ld.h5. The linearisation now lives in chunk_grid (shared with the writers), and slots beyond the current extent are ignored as the library does. read_fixed_array_chunks / read_extensible_array_chunks take the dataspace's max dimensions. Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/chunk_grid.rs | 201 ++++++++++++++ crates/clawhdf5-format/src/chunked_read.rs | 2 + .../clawhdf5-format/src/extensible_array.rs | 159 ++++------- crates/clawhdf5-format/src/fixed_array.rs | 101 ++----- crates/clawhdf5-format/src/lib.rs | 1 + crates/clawhdf5/tests/chunk_index_interop.rs | 256 ++++++++++++++++++ 6 files changed, 546 insertions(+), 174 deletions(-) create mode 100644 crates/clawhdf5-format/src/chunk_grid.rs create mode 100644 crates/clawhdf5/tests/chunk_index_interop.rs diff --git a/crates/clawhdf5-format/src/chunk_grid.rs b/crates/clawhdf5-format/src/chunk_grid.rs new file mode 100644 index 0000000..c5a06fd --- /dev/null +++ b/crates/clawhdf5-format/src/chunk_grid.rs @@ -0,0 +1,201 @@ +//! Chunk-index linearisation shared by the Fixed Array and Extensible Array +//! chunk indexes (reader and writer). +//! +//! Both indexes store one element per chunk at a *linear* index, and the +//! library derives that index from the chunk's scaled coordinates +//! (`offset / chunk_dim`) using the dataset's **maximum** dimensions, not its +//! current ones (`H5D__farray_idx_get_addr` / `H5D__earray_idx_get_addr`, +//! via `layout->max_down_chunks`). A dataset whose current shape is smaller +//! than its maxshape therefore has gaps in the index, and laying it out by the +//! current shape puts every chunk after the first row in the wrong place. +//! +//! The Extensible Array adds one more step: its one unlimited dimension has no +//! finite chunk count, so the library *swizzles* the coordinates to make that +//! dimension the slowest-varying one (`H5VM_swizzle_coords`, which moves +//! `coords[unlim_dim]` to the front and shifts the dimensions before it right +//! by one) before linearising with `swizzled_max_down_chunks`. When the +//! unlimited dimension is already dimension 0 no swizzle happens. + +#[cfg(not(feature = "std"))] +extern crate alloc; + +#[cfg(not(feature = "std"))] +use alloc::{vec, vec::Vec}; + +use crate::error::FormatError; + +/// How a chunk index maps linear element indexes to chunk coordinates. +#[derive(Debug, Clone)] +pub(crate) struct ChunkGrid { + /// Spatial chunk dimensions, in dataset order. + chunk_dims: Vec, + /// Chunks per dimension covering the *current* extent, in dataset order. + cur_chunks: Vec, + /// Dataset dimension stored at each linearisation position (slowest + /// first). The identity except for a swizzled Extensible Array. + order: Vec, + /// Linear stride of each linearisation position. + down: Vec, +} + +impl ChunkGrid { + /// Grid for a Fixed Array index: row-major over the chunk counts of the + /// maximum dimensions (`max_dims`, falling back to the current dimensions + /// when the dataspace records none). + pub(crate) fn fixed_array( + cur_dims: &[u64], + max_dims: Option<&[u64]>, + chunk_dims: &[u64], + ) -> Result { + Self::build(cur_dims, max_dims, chunk_dims, None) + } + + /// Grid for an Extensible Array index: like the Fixed Array, but the + /// unlimited dimension (the one whose maximum is `H5S_UNLIMITED`) is moved + /// to the slowest-varying position first. + pub(crate) fn extensible_array( + cur_dims: &[u64], + max_dims: Option<&[u64]>, + chunk_dims: &[u64], + ) -> Result { + let unlim = max_dims.and_then(|m| m.iter().position(|&d| d == u64::MAX)); + Self::build(cur_dims, max_dims, chunk_dims, unlim) + } + + fn build( + cur_dims: &[u64], + max_dims: Option<&[u64]>, + chunk_dims: &[u64], + unlim: Option, + ) -> Result { + let rank = chunk_dims.len(); + if cur_dims.len() != rank || max_dims.is_some_and(|m| m.len() != rank) { + return Err(FormatError::ChunkedReadError( + "chunk index rank does not match the dataspace".into(), + )); + } + if chunk_dims.contains(&0) { + return Err(FormatError::ChunkedReadError( + "chunk dimension is zero".into(), + )); + } + let cur_chunks: Vec = cur_dims + .iter() + .zip(chunk_dims) + .map(|(&d, &c)| d.div_ceil(c)) + .collect(); + // Chunk counts of the maximum extent. An unlimited dimension has no + // finite count; it only ever sits in the slowest position, where its + // count never enters a stride. A (corrupt) maximum smaller than the + // current extent is widened so no allocated chunk becomes unreachable. + let max_chunks: Vec = (0..rank) + .map(|d| { + let max = max_dims.map_or(cur_dims[d], |m| m[d]); + if max == u64::MAX { + u64::MAX + } else { + max.div_ceil(chunk_dims[d]).max(cur_chunks[d]) + } + }) + .collect(); + + let mut order: Vec = (0..rank).collect(); + if let Some(u) = unlim { + order.remove(u); + order.insert(0, u); + } + let mut down = vec![1u64; rank]; + for p in (0..rank.saturating_sub(1)).rev() { + let next = max_chunks[order[p + 1]]; + if next == u64::MAX { + // Only reachable with more than one unlimited dimension, which + // neither index type can describe. + return Err(FormatError::ChunkedReadError( + "array chunk index with more than one unlimited dimension".into(), + )); + } + down[p] = down[p + 1].checked_mul(next).ok_or_else(|| { + FormatError::Overflow("chunk index linear stride overflows u64".into()) + })?; + } + Ok(Self { + chunk_dims: chunk_dims.to_vec(), + cur_chunks, + order, + down, + }) + } + + /// Dataset-space offsets of the chunk stored at linear `index`, or `None` + /// when that chunk lies outside the current extent (the index still has a + /// slot for it; the library ignores such chunks on read). + pub(crate) fn offsets(&self, index: u64) -> Option> { + let rank = self.chunk_dims.len(); + let mut offsets = vec![0u64; rank]; + let mut rem = index; + for p in 0..rank { + let d = self.order[p]; + let scaled = rem / self.down[p]; + rem %= self.down[p]; + if scaled >= self.cur_chunks[d] { + return None; + } + offsets[d] = scaled * self.chunk_dims[d]; + } + Some(offsets) + } + + /// Linear index of the chunk with scaled coordinates `scaled` + /// (`offset / chunk_dim` per dimension, in dataset order). + #[allow(dead_code)] // used by the writer + pub(crate) fn linear_index(&self, scaled: &[u64]) -> u64 { + self.order + .iter() + .zip(&self.down) + .map(|(&d, &stride)| scaled[d] * stride) + .sum() + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn fixed_array_uses_max_dims() { + // shape (4, 6), chunks (2, 3), maxshape (20, 10): 10 x 4 chunk grid. + let g = ChunkGrid::fixed_array(&[4, 6], Some(&[20, 10]), &[2, 3]).unwrap(); + assert_eq!(g.offsets(0), Some(vec![0, 0])); + assert_eq!(g.offsets(1), Some(vec![0, 3])); + assert_eq!(g.offsets(2), None); // column chunk 2 is beyond the extent + assert_eq!(g.offsets(4), Some(vec![2, 0])); + assert_eq!(g.offsets(5), Some(vec![2, 3])); + assert_eq!(g.offsets(8), None); // row chunk 2 is beyond the extent + assert_eq!(g.linear_index(&[1, 1]), 5); + } + + #[test] + fn extensible_array_swizzles_unlimited_dim() { + // maxshape (10, None): dim 1 is unlimited and becomes slowest. + let g = ChunkGrid::extensible_array(&[4, 6], Some(&[10, u64::MAX]), &[2, 3]).unwrap(); + // max chunks of dim 0 = 5, so index = c1 * 5 + c0. + assert_eq!(g.linear_index(&[1, 0]), 1); + assert_eq!(g.linear_index(&[0, 1]), 5); + assert_eq!(g.offsets(5), Some(vec![0, 3])); + assert_eq!(g.offsets(6), Some(vec![2, 3])); + assert_eq!(g.offsets(2), None); + } + + #[test] + fn extensible_array_unlimited_first_is_row_major() { + let g = ChunkGrid::extensible_array(&[4, 6], Some(&[u64::MAX, 30]), &[2, 3]).unwrap(); + // max chunks of dim 1 = 10. + assert_eq!(g.linear_index(&[1, 1]), 11); + assert_eq!(g.offsets(11), Some(vec![2, 3])); + } + + #[test] + fn rejects_two_unlimited_dims_after_the_first() { + assert!(ChunkGrid::fixed_array(&[4, 6], Some(&[u64::MAX, u64::MAX]), &[2, 3]).is_err()); + } +} diff --git a/crates/clawhdf5-format/src/chunked_read.rs b/crates/clawhdf5-format/src/chunked_read.rs index 8dade54..f07ae27 100644 --- a/crates/clawhdf5-format/src/chunked_read.rs +++ b/crates/clawhdf5-format/src/chunked_read.rs @@ -593,6 +593,7 @@ pub fn list_chunks( file_data, &header, &dataspace.dimensions, + dataspace.max_dimensions.as_deref(), spatial_chunk_dims, elem_size as u32, offset_size, @@ -608,6 +609,7 @@ pub fn list_chunks( file_data, &header, &dataspace.dimensions, + dataspace.max_dimensions.as_deref(), spatial_chunk_dims, elem_size as u32, offset_size, diff --git a/crates/clawhdf5-format/src/extensible_array.rs b/crates/clawhdf5-format/src/extensible_array.rs index 4416a76..ba1c822 100644 --- a/crates/clawhdf5-format/src/extensible_array.rs +++ b/crates/clawhdf5-format/src/extensible_array.rs @@ -9,6 +9,7 @@ extern crate alloc; #[cfg(not(feature = "std"))] use alloc::{format, vec, vec::Vec}; +use crate::chunk_grid::ChunkGrid; use crate::chunked_read::ChunkInfo; use crate::error::FormatError; @@ -203,8 +204,7 @@ fn read_element( offset_size: u8, chunk_byte_size: u64, linear_index: usize, - num_chunks_per_dim: &[u64], - chunk_dimensions: &[u32], + grid: &ChunkGrid, ) -> Result<(Option, usize), FormatError> { let os = offset_size as usize; @@ -220,7 +220,10 @@ fn read_element( return Ok((None, os)); } let address = read_offset(data, pos, offset_size)?; - let offsets = index_to_chunk_offsets(linear_index, num_chunks_per_dim, chunk_dimensions); + // A slot beyond the current extent is ignored, as the library does. + let Some(offsets) = grid.offsets(linear_index as u64) else { + return Ok((None, os)); + }; Ok(( Some(ChunkInfo { chunk_size: chunk_byte_size as u32, @@ -261,7 +264,9 @@ fn read_element( data[fm_off + 2], data[fm_off + 3], ]); - let offsets = index_to_chunk_offsets(linear_index, num_chunks_per_dim, chunk_dimensions); + let Some(offsets) = grid.offsets(linear_index as u64) else { + return Ok((None, elem_total)); + }; Ok(( Some(ChunkInfo { chunk_size: chunk_size as u32, @@ -274,27 +279,6 @@ fn read_element( } } -/// Convert a linear chunk index to N-dimensional chunk offsets in dataset space. -fn index_to_chunk_offsets( - index: usize, - num_chunks_per_dim: &[u64], - chunk_dimensions: &[u32], -) -> Vec { - let rank = num_chunks_per_dim.len(); - let mut offsets = vec![0u64; rank]; - let mut remaining = index as u64; - for d in (0..rank).rev() { - let nchunks = num_chunks_per_dim[d]; - if nchunks == 0 { - continue; - } - let chunk_idx = remaining % nchunks; - remaining /= nchunks; - offsets[d] = chunk_idx * chunk_dimensions[d] as u64; - } - offsets -} - /// Collect elements from a data block at the given offset. #[allow(clippy::too_many_arguments)] /// Layout of super block `u`, per the HDF5 spec: the number of data blocks it @@ -339,8 +323,7 @@ fn read_data_block_elements( offset_size: u8, chunk_byte_size: u64, start_index: usize, - num_chunks_per_dim: &[u64], - chunk_dimensions: &[u32], + grid: &ChunkGrid, page_init: &[u8], first_page: usize, ) -> Result, FormatError> { @@ -376,8 +359,7 @@ fn read_data_block_elements( offset_size, chunk_byte_size, first_index + i, - num_chunks_per_dim, - chunk_dimensions, + grid, )?; if let Some(ci) = info { chunks.push(ci); @@ -449,25 +431,19 @@ pub fn read_extensible_array_chunks( file_data: &[u8], header: &ExtensibleArrayHeader, dataset_dims: &[u64], + max_dims: Option<&[u64]>, chunk_dimensions: &[u32], element_size: u32, offset_size: u8, _length_size: u8, ) -> Result, FormatError> { - let rank = chunk_dimensions.len(); let os = offset_size as usize; - let mut num_chunks_per_dim = Vec::with_capacity(rank); - for d in 0..rank { - let ch_dim = chunk_dimensions[d] as u64; - if ch_dim == 0 { - return Err(FormatError::ChunkedReadError( - "chunk dimension is zero".into(), - )); - } - let ds_dim = dataset_dims[d]; - num_chunks_per_dim.push(ds_dim.div_ceil(ch_dim)); - } + // Linear indexes follow the maximum dimensions, with the unlimited + // dimension swizzled to the slowest position (see `chunk_grid`). + let dims_u64: Vec = chunk_dimensions.iter().map(|&d| d as u64).collect(); + let grid = ChunkGrid::extensible_array(dataset_dims, max_dims, &dims_u64)?; + let grid = &grid; let chunk_byte_size: u64 = chunk_dimensions.iter().map(|&d| d as u64).product::() * element_size as u64; @@ -557,8 +533,7 @@ pub fn read_extensible_array_chunks( offset_size, chunk_byte_size, i, - &num_chunks_per_dim, - chunk_dimensions, + grid, )?; if let Some(ci) = info { chunks.push(ci); @@ -594,8 +569,7 @@ pub fn read_extensible_array_chunks( offset_size, chunk_byte_size, global_index, - &num_chunks_per_dim, - chunk_dimensions, + grid, &[], 0, )?); @@ -625,8 +599,7 @@ pub fn read_extensible_array_chunks( offset_size, chunk_byte_size, global_index, - &num_chunks_per_dim, - chunk_dimensions, + grid, )?); } global_index = @@ -653,8 +626,7 @@ fn read_super_block( offset_size: u8, chunk_byte_size: u64, start_index: usize, - num_chunks_per_dim: &[u64], - chunk_dimensions: &[u32], + grid: &ChunkGrid, ) -> Result, FormatError> { let os = offset_size as usize; let sb_header_size = 4 + 1 + 1 + os + arr_off_size(header); @@ -710,8 +682,7 @@ fn read_super_block( offset_size, chunk_byte_size, global_idx, - num_chunks_per_dim, - chunk_dimensions, + grid, bitmap, i * npages, )?); @@ -735,35 +706,18 @@ mod tests { } #[test] fn index_to_offsets_1d() { - let num_chunks = vec![5u64]; - let chunk_dims = vec![20u32]; - assert_eq!(index_to_chunk_offsets(0, &num_chunks, &chunk_dims), vec![0]); - assert_eq!( - index_to_chunk_offsets(1, &num_chunks, &chunk_dims), - vec![20] - ); - assert_eq!( - index_to_chunk_offsets(4, &num_chunks, &chunk_dims), - vec![80] - ); + let g = ChunkGrid::fixed_array(&[100], None, &[20]).unwrap(); + assert_eq!(g.offsets(0).unwrap(), vec![0]); + assert_eq!(g.offsets(1).unwrap(), vec![20]); + assert_eq!(g.offsets(4).unwrap(), vec![80]); } #[test] fn index_to_offsets_2d() { - let num_chunks = vec![3u64, 2]; - let chunk_dims = vec![4u32, 3]; - assert_eq!( - index_to_chunk_offsets(0, &num_chunks, &chunk_dims), - vec![0, 0] - ); - assert_eq!( - index_to_chunk_offsets(1, &num_chunks, &chunk_dims), - vec![0, 3] - ); - assert_eq!( - index_to_chunk_offsets(2, &num_chunks, &chunk_dims), - vec![4, 0] - ); + let g = ChunkGrid::fixed_array(&[10, 6], None, &[4, 3]).unwrap(); + assert_eq!(g.offsets(0).unwrap(), vec![0, 0]); + assert_eq!(g.offsets(1).unwrap(), vec![0, 3]); + assert_eq!(g.offsets(2).unwrap(), vec![4, 0]); } #[test] @@ -830,7 +784,7 @@ mod tests { index_block_address: (usize::MAX - 4) as u64, }; let buf = vec![0u8; 64]; - let r = read_extensible_array_chunks(&buf, &header, &[100], &[20], 8, 8, 8); + let r = read_extensible_array_chunks(&buf, &header, &[100], None, &[20], 8, 8, 8); assert!(r.is_err()); } @@ -913,9 +867,17 @@ mod tests { let header = ExtensibleArrayHeader::parse(&file_data, aehd_offset, os, ls).unwrap(); let ds_dims = vec![40u64]; // 2 chunks × 20 elements let chunk_dims = vec![20u32]; - let chunks = - read_extensible_array_chunks(&file_data, &header, &ds_dims, &chunk_dims, 8, os, ls) - .unwrap(); + let chunks = read_extensible_array_chunks( + &file_data, + &header, + &ds_dims, + None, + &chunk_dims, + 8, + os, + ls, + ) + .unwrap(); assert_eq!(chunks.len(), 2); assert_eq!(chunks[0].address, base_addr); @@ -1023,9 +985,17 @@ mod tests { let header = ExtensibleArrayHeader::parse(&file_data, aehd_offset, os, ls).unwrap(); let ds_dims = vec![40u64]; let chunk_dims = vec![10u32]; - let chunks = - read_extensible_array_chunks(&file_data, &header, &ds_dims, &chunk_dims, 8, os, ls) - .unwrap(); + let chunks = read_extensible_array_chunks( + &file_data, + &header, + &ds_dims, + None, + &chunk_dims, + 8, + os, + ls, + ) + .unwrap(); assert_eq!(chunks.len(), 4); for (i, c) in chunks.iter().enumerate() { @@ -1047,10 +1017,8 @@ mod tests { #[test] fn read_element_unallocated() { let data = vec![0xFFu8; 16]; - let num_chunks = vec![5u64]; - let chunk_dims = vec![10u32]; - let (info, consumed) = - read_element(&data, 0, 0, 8, 8, 80, 0, &num_chunks, &chunk_dims).unwrap(); + let grid = ChunkGrid::fixed_array(&[50], None, &[10]).unwrap(); + let (info, consumed) = read_element(&data, 0, 0, 8, 8, 80, 0, &grid).unwrap(); assert!(info.is_none()); assert_eq!(consumed, 8); } @@ -1069,20 +1037,9 @@ mod tests { // Filter mask data[12..16].copy_from_slice(&0u32.to_le_bytes()); - let num_chunks = vec![5u64]; - let chunk_dims = vec![10u32]; - let (info, consumed) = read_element( - &data, - 0, - 1, - elem_size as u8, - os, - 80, - 2, - &num_chunks, - &chunk_dims, - ) - .unwrap(); + let grid = ChunkGrid::fixed_array(&[50], None, &[10]).unwrap(); + let (info, consumed) = + read_element(&data, 0, 1, elem_size as u8, os, 80, 2, &grid).unwrap(); let ci = info.unwrap(); assert_eq!(ci.address, 0x2000); assert_eq!(ci.chunk_size, 120); diff --git a/crates/clawhdf5-format/src/fixed_array.rs b/crates/clawhdf5-format/src/fixed_array.rs index 4c79f03..8546224 100644 --- a/crates/clawhdf5-format/src/fixed_array.rs +++ b/crates/clawhdf5-format/src/fixed_array.rs @@ -6,6 +6,7 @@ extern crate alloc; #[cfg(not(feature = "std"))] use alloc::{format, vec, vec::Vec}; +use crate::chunk_grid::ChunkGrid; use crate::chunked_read::ChunkInfo; use crate::error::FormatError; @@ -151,13 +152,13 @@ pub fn read_fixed_array_chunks( file_data: &[u8], header: &FixedArrayHeader, dataset_dims: &[u64], + max_dims: Option<&[u64]>, chunk_dimensions: &[u32], element_size: u32, offset_size: u8, _length_size: u8, ) -> Result, FormatError> { let db_offset = header.data_block_address as usize; - let rank = chunk_dimensions.len(); // Parse data block header: FADB(4) + version(1) + client_id(1) + header_address(offset_size) let db_header_size = 4 + 1 + 1 + offset_size as usize; @@ -198,19 +199,10 @@ pub fn read_fixed_array_chunks( )) }; - // Compute chunk offsets based on index. - // Chunks are stored in row-major order within the dataset space. - let mut num_chunks_per_dim = Vec::with_capacity(rank); - for d_idx in 0..rank { - let ch_dim = chunk_dimensions[d_idx] as u64; - if ch_dim == 0 { - return Err(FormatError::ChunkedReadError( - "chunk dimension is zero".into(), - )); - } - let ds_dim = dataset_dims[d_idx]; - num_chunks_per_dim.push(ds_dim.div_ceil(ch_dim)); - } + // The index is laid out over the chunk grid of the *maximum* dimensions + // (row-major), so a dataset smaller than its maxshape has gaps. + let dims_u64: Vec = chunk_dimensions.iter().map(|&d| d as u64).collect(); + let grid = ChunkGrid::fixed_array(dataset_dims, max_dims, &dims_u64)?; let chunk_byte_size: u64 = chunk_dimensions.iter().map(|&d| d as u64).product::() * element_size as u64; @@ -226,7 +218,11 @@ pub fn read_fixed_array_chunks( header.element_size, chunk_byte_size, )? { - let offsets = index_to_chunk_offsets(i, &num_chunks_per_dim, chunk_dimensions); + // A slot beyond the current extent is ignored, as the + // library does. + let Some(offsets) = grid.offsets(i as u64) else { + return Ok(()); + }; chunks.push(ChunkInfo { chunk_size, filter_mask, @@ -367,27 +363,6 @@ fn parse_fa_element( } } -/// Convert a linear chunk index to N-dimensional chunk offsets in dataset space. -fn index_to_chunk_offsets( - index: usize, - num_chunks_per_dim: &[u64], - chunk_dimensions: &[u32], -) -> Vec { - let rank = num_chunks_per_dim.len(); - let mut offsets = vec![0u64; rank]; - let mut remaining = index as u64; - for d in (0..rank).rev() { - let nchunks = num_chunks_per_dim[d]; - if nchunks == 0 { - continue; - } - let chunk_idx = remaining % nchunks; - remaining /= nchunks; - offsets[d] = chunk_idx * chunk_dimensions[d] as u64; - } - offsets -} - /// Read a variable-length little-endian unsigned integer. fn read_variable_length(data: &[u8], size: usize) -> Result { if size > 8 || data.len() < size { @@ -416,44 +391,21 @@ mod tests { #[test] fn index_to_offsets_1d() { - let num_chunks = vec![5u64]; - let chunk_dims = vec![20u32]; - assert_eq!(index_to_chunk_offsets(0, &num_chunks, &chunk_dims), vec![0]); - assert_eq!( - index_to_chunk_offsets(1, &num_chunks, &chunk_dims), - vec![20] - ); - assert_eq!( - index_to_chunk_offsets(4, &num_chunks, &chunk_dims), - vec![80] - ); + let g = ChunkGrid::fixed_array(&[100], None, &[20]).unwrap(); + assert_eq!(g.offsets(0).unwrap(), vec![0]); + assert_eq!(g.offsets(1).unwrap(), vec![20]); + assert_eq!(g.offsets(4).unwrap(), vec![80]); } #[test] fn index_to_offsets_2d() { // 10x6 dataset with 4x3 chunks => ceil(10/4)=3, ceil(6/3)=2 => 6 chunks - let num_chunks = vec![3u64, 2]; - let chunk_dims = vec![4u32, 3]; - assert_eq!( - index_to_chunk_offsets(0, &num_chunks, &chunk_dims), - vec![0, 0] - ); - assert_eq!( - index_to_chunk_offsets(1, &num_chunks, &chunk_dims), - vec![0, 3] - ); - assert_eq!( - index_to_chunk_offsets(2, &num_chunks, &chunk_dims), - vec![4, 0] - ); - assert_eq!( - index_to_chunk_offsets(3, &num_chunks, &chunk_dims), - vec![4, 3] - ); - assert_eq!( - index_to_chunk_offsets(5, &num_chunks, &chunk_dims), - vec![8, 3] - ); + let g = ChunkGrid::fixed_array(&[10, 6], None, &[4, 3]).unwrap(); + assert_eq!(g.offsets(0).unwrap(), vec![0, 0]); + assert_eq!(g.offsets(1).unwrap(), vec![0, 3]); + assert_eq!(g.offsets(2).unwrap(), vec![4, 0]); + assert_eq!(g.offsets(3).unwrap(), vec![4, 3]); + assert_eq!(g.offsets(5).unwrap(), vec![8, 3]); } #[test] @@ -517,7 +469,7 @@ mod tests { let read = |f: &[u8], fahd: usize| -> Result, FormatError> { let h = FixedArrayHeader::parse(f, fahd, 8, 8)?; - read_fixed_array_chunks(f, &h, &[60], &[20], 8, 8, 8) + read_fixed_array_chunks(f, &h, &[60], None, &[20], 8, 8, 8) }; let (clean, fahd) = build(); @@ -562,7 +514,7 @@ mod tests { let db = 0x100usize; buf[db..db + 4].copy_from_slice(b"FADB"); let header = FixedArrayHeader::parse(&buf, fahd, 8, 8).unwrap(); - let r = read_fixed_array_chunks(&buf, &header, &[100], &[20], 8, 8, 8); + let r = read_fixed_array_chunks(&buf, &header, &[100], None, &[20], 8, 8, 8); assert!(r.is_err()); } @@ -579,7 +531,7 @@ mod tests { stamp_checksum(&mut buf, fahd, fahd + 24); buf[0x80..0x84].copy_from_slice(b"FADB"); let header = FixedArrayHeader::parse(&buf, fahd, 8, 8).unwrap(); - let r = read_fixed_array_chunks(&buf, &header, &[100], &[20], 8, 8, 8); + let r = read_fixed_array_chunks(&buf, &header, &[100], None, &[20], 8, 8, 8); assert!(r.is_err()); } @@ -602,7 +554,7 @@ mod tests { data_block_address: (usize::MAX - 4) as u64, }; let buf = vec![0u8; 64]; - let r = read_fixed_array_chunks(&buf, &header, &[100], &[20], 8, 8, 8); + let r = read_fixed_array_chunks(&buf, &header, &[100], None, &[20], 8, 8, 8); assert!(r.is_err()); } @@ -664,6 +616,7 @@ mod tests { &file_data, &header, &ds_dims, + None, &chunk_dims, 8, offset_size, @@ -740,6 +693,7 @@ mod tests { &file_data, &header, &ds_dims, + None, &chunk_dims, 8, offset_size, @@ -840,6 +794,7 @@ mod tests { &file_data, &header, &ds_dims, + None, &chunk_dims, 8, offset_size, diff --git a/crates/clawhdf5-format/src/lib.rs b/crates/clawhdf5-format/src/lib.rs index 54bd326..4a6c2a1 100644 --- a/crates/clawhdf5-format/src/lib.rs +++ b/crates/clawhdf5-format/src/lib.rs @@ -54,6 +54,7 @@ pub mod btree_v1; pub mod btree_v2; pub mod checksum; pub mod chunk_cache; +mod chunk_grid; pub mod chunk_index; pub mod chunked_read; pub mod chunked_write; diff --git a/crates/clawhdf5/tests/chunk_index_interop.rs b/crates/clawhdf5/tests/chunk_index_interop.rs new file mode 100644 index 0000000..9c030a3 --- /dev/null +++ b/crates/clawhdf5/tests/chunk_index_interop.rs @@ -0,0 +1,256 @@ +//! Fixed Array / Extensible Array chunk-index interop with libhdf5 (via h5py). +//! +//! Both indexes place each chunk at a linear index computed from the +//! dataset's *maximum* dimensions, and the Extensible Array additionally +//! moves its unlimited dimension to the slowest-varying position. Getting +//! either wrong reads (or writes) every chunk after the first row in the +//! wrong place, silently, so these tests compare every value. +//! +//! Skipped when python3 with h5py is unavailable, unless +//! `CLAWHDF5_REQUIRE_INTEROP=1`. + +use std::process::Command; + +use clawhdf5::File; + +fn python() -> String { + std::env::var("CLAWHDF5_PYTHON").unwrap_or_else(|_| "python3".to_string()) +} + +fn interop_required() -> bool { + std::env::var("CLAWHDF5_REQUIRE_INTEROP").is_ok_and(|v| v == "1") +} + +fn python_available() -> bool { + Command::new(python()) + .args(["-c", "import h5py"]) + .output() + .map(|o| o.status.success()) + .unwrap_or(false) +} + +macro_rules! skip_if_no_python { + () => { + if !python_available() { + assert!( + !interop_required(), + "CLAWHDF5_REQUIRE_INTEROP=1 but python3 with h5py is not available" + ); + eprintln!("SKIP: python3 with h5py not available"); + return; + } + }; +} + +fn run_python(script: &str) -> String { + let output = Command::new(python()) + .args(["-c", script]) + .output() + .expect("failed to run python"); + if !output.status.success() { + panic!( + "Python script failed:\nSTDOUT: {}\nSTDERR: {}", + String::from_utf8_lossy(&output.stdout), + String::from_utf8_lossy(&output.stderr) + ); + } + String::from_utf8_lossy(&output.stdout).trim().to_string() +} + +/// Row-major `arange` of `shape`, cropped to `crop` (the current extent). +fn arange_cropped(full: &[usize], crop: &[usize]) -> Vec { + let n: usize = crop.iter().product(); + let mut out = Vec::with_capacity(n); + for flat in 0..n { + let mut rem = flat; + let mut src = 0usize; + let mut stride = 1usize; + let mut coords = vec![0usize; crop.len()]; + for d in (0..crop.len()).rev() { + coords[d] = rem % crop[d]; + rem /= crop[d]; + } + for d in (0..full.len()).rev() { + src += coords[d] * stride; + stride *= full[d]; + } + out.push(src as i32); + } + out +} + +/// One `i4` dataset, filled with `arange` over `full` and then resized to +/// `shape` (equal to `full` unless the case shrinks it). +struct Case { + name: &'static str, + full: Vec, + shape: Vec, + chunks: Vec, + maxshape: &'static str, + extra: &'static str, + index: &'static str, +} + +fn py_tuple(v: &[usize]) -> String { + let parts: Vec = v.iter().map(|x| x.to_string()).collect(); + format!("({},)", parts.join(",")) +} + +/// Have h5py (`libver="latest"`, so Fixed/Extensible Array indexes) write +/// every case to one file, then read each back and compare every value. +fn check_h5py_written(cases: &[Case]) { + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("h5py_chunk_index.h5"); + let path_str = path.display().to_string(); + + let mut script = + format!("import h5py, numpy as np\nf = h5py.File(r'{path_str}', 'w', libver='latest')\n"); + for c in cases { + script += &format!( + "d = f.create_dataset('{name}', data=np.arange({n}, dtype='i4').reshape({full}), \ + chunks={chunks}, maxshape={maxshape}{extra})\n\ + d.resize({shape})\n", + name = c.name, + n = c.full.iter().product::(), + full = py_tuple(&c.full), + chunks = py_tuple(&c.chunks), + maxshape = c.maxshape, + extra = c.extra, + shape = py_tuple(&c.shape), + ); + } + script += "f.close()\n"; + run_python(&script); + + let file = File::open(&path).unwrap(); + for c in cases { + let ds = file.dataset(c.name).unwrap(); + let shape: Vec = ds.shape().unwrap().iter().map(|&d| d as usize).collect(); + assert_eq!(shape, c.shape, "{}: shape", c.name); + let got = ds.read_i32().unwrap(); + let want = arange_cropped(&c.full, &c.shape); + let bad = got.iter().zip(&want).filter(|(a, b)| a != b).count(); + assert_eq!( + got, + want, + "{}: {bad} of {} values differ (index {})", + c.name, + want.len(), + c.index + ); + } +} + +/// h5py-written Extensible Array whose unlimited dimension is not the first, +/// with the current shape smaller than the finite maximum: the library +/// swizzles the unlimited dimension to the slowest position and strides the +/// rest by their maximum chunk counts. +#[test] +fn h5py_extensible_array_partial_extent_reads_correctly() { + skip_if_no_python!(); + check_h5py_written(&[ + // The `ea_fa_partial.h5` repro from the conformance sweep. + Case { + name: "ea_10_none", + full: vec![4, 6], + shape: vec![4, 6], + chunks: vec![2, 3], + maxshape: "(10, None)", + extra: "", + index: "EA, unlimited dim 1", + }, + Case { + name: "ea_none_10", + full: vec![4, 6], + shape: vec![4, 6], + chunks: vec![2, 3], + maxshape: "(None, 10)", + extra: "", + index: "EA, unlimited dim 0", + }, + Case { + name: "ea_3d_mid", + full: vec![3, 4, 5], + shape: vec![3, 4, 5], + chunks: vec![2, 3, 2], + maxshape: "(5, None, 7)", + extra: "", + index: "EA, unlimited dim 1 of 3", + }, + Case { + name: "ea_3d_last_gzip", + full: vec![3, 4, 5], + shape: vec![3, 4, 5], + chunks: vec![2, 3, 2], + maxshape: "(5, 9, None)", + extra: ", compression='gzip'", + index: "EA, unlimited dim 2 of 3, filtered", + }, + // Many chunks: crosses data blocks, super blocks and paging. + Case { + name: "ea_many", + full: vec![3, 1500], + shape: vec![3, 1500], + chunks: vec![1, 1], + maxshape: "(4, None)", + extra: "", + index: "EA, 4500 slots", + }, + // Shrunk after writing: chunks beyond the extent must be ignored. + Case { + name: "ea_shrunk", + full: vec![8, 9], + shape: vec![3, 4], + chunks: vec![2, 3], + maxshape: "(10, None)", + extra: "", + index: "EA, shrunk", + }, + ]); +} + +/// h5py-written Fixed Array with the current shape smaller than a finite +/// maxshape: the index has one slot per chunk of the *maximum* extent. +#[test] +fn h5py_fixed_array_partial_extent_reads_correctly() { + skip_if_no_python!(); + check_h5py_written(&[ + Case { + name: "fa_20_10", + full: vec![4, 6], + shape: vec![4, 6], + chunks: vec![2, 3], + maxshape: "(20, 10)", + extra: "", + index: "FA", + }, + Case { + name: "fa_3d_gzip", + full: vec![3, 4, 5], + shape: vec![3, 4, 5], + chunks: vec![2, 3, 2], + maxshape: "(6, 8, 10)", + extra: ", compression='gzip'", + index: "FA, filtered", + }, + // Paged (> 1024 slots) with most of them beyond the extent. + Case { + name: "fa_paged", + full: vec![30, 50], + shape: vec![30, 50], + chunks: vec![1, 1], + maxshape: "(40, 60)", + extra: "", + index: "FA, 2400 slots, paged", + }, + Case { + name: "fa_shrunk", + full: vec![8, 9], + shape: vec![5, 2], + chunks: vec![2, 3], + maxshape: "(20, 10)", + extra: "", + index: "FA, shrunk", + }, + ]); +} From 4a1876faf2dd865ed20948ce10d2d647b75d237b Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:07:08 -0500 Subject: [PATCH 2/7] fix(format): page Fixed Array data blocks past 1024 chunks The Fixed Array writer always packed every element into one data block behind one checksum. Past 2^10 elements libhdf5 (and our reader) expect a paged block: a page-init bitmap after the prefix, then one checksummed page per 1024 elements. Any dataset with more than 1024 chunks and no unlimited dimension failed with "incorrect metadata checksum" in h5py, h5dump and our own reader. build_fixed_array_at now takes one Option per array slot so later fixes can leave unallocated slots. Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/chunked_write.rs | 172 +++++++++++-------- crates/clawhdf5/tests/chunk_index_interop.rs | 138 ++++++++++++++- 2 files changed, 240 insertions(+), 70 deletions(-) diff --git a/crates/clawhdf5-format/src/chunked_write.rs b/crates/clawhdf5-format/src/chunked_write.rs index 6816546..c0d02e1 100644 --- a/crates/clawhdf5-format/src/chunked_write.rs +++ b/crates/clawhdf5-format/src/chunked_write.rs @@ -499,108 +499,140 @@ fn serialize_v4_fixed_array( buf } +/// log2 of the elements per Fixed Array data block page (the library's +/// default, `H5D_FARRAY_MAX_DBLK_PAGE_NELMTS_BITS`). +const FA_PAGE_BITS: u8 = 10; + +fn push_addr(buf: &mut Vec, addr: u64, offset_size: u8) { + match offset_size { + 4 => buf.extend_from_slice(&(addr as u32).to_le_bytes()), + _ => buf.extend_from_slice(&addr.to_le_bytes()), + } +} + +/// Width of the chunk-size field of a filtered chunk index element. Must +/// match the library's `H5D_FARRAY_FILT_COMPUTE_CHUNK_SIZE_LEN` (the EA and +/// B-tree v2 indexes use the same formula): +/// `1 + ((log2(unfiltered chunk bytes) + 8) / 8)`, capped at 8. +pub(crate) fn filtered_chunk_size_len(slots: &[Option]) -> usize { + let max_raw = slots + .iter() + .flatten() + .map(|c| c.raw_size) + .max() + .unwrap_or(1); + let log2_val = if max_raw <= 1 { + 0 + } else { + 63 - max_raw.leading_zeros() + }; + (1 + ((log2_val + 8) / 8) as usize).min(8) +} + +/// Append one chunk index element: the chunk's address, plus its stored size +/// and filter mask when the dataset is filtered. `None` is an unallocated +/// chunk (undefined address, zero size and mask). +pub(crate) fn push_index_element( + buf: &mut Vec, + slot: Option<&WrittenChunk>, + offset_size: u8, + chunk_size_bytes: Option, +) { + match slot { + Some(c) => { + push_addr(buf, c.address, offset_size); + if let Some(n) = chunk_size_bytes { + buf.extend_from_slice(&c.compressed_size.to_le_bytes()[..n]); + buf.extend_from_slice(&c.filter_mask.to_le_bytes()); + } + } + None => { + buf.extend(core::iter::repeat_n(0xFF, offset_size as usize)); + if let Some(n) = chunk_size_bytes { + buf.extend(core::iter::repeat_n(0x00, n + 4)); + } + } + } +} + /// Build a complete Fixed Array at a known absolute address. +/// +/// `slots` holds one entry per element of the array, i.e. per chunk of the +/// dataset's *maximum* extent in the order [`crate::chunk_grid`] defines; +/// `None` marks a chunk that is not allocated. An array with more elements +/// than fit in one page (`2^FA_PAGE_BITS`) gets a paged data block: a +/// page-init bitmap after the prefix, then one checksummed page per +/// `2^FA_PAGE_BITS` elements, the last one short (`H5FA__dblock_create`). pub fn build_fixed_array_at( - chunks: &[WrittenChunk], + slots: &[Option], offset_size: u8, length_size: u8, has_filters: bool, fa_base_address: u64, ) -> Vec { let os = offset_size as usize; - let num_elements = chunks.len(); - - // For filtered chunks, compute chunk_size encoding width. - // Must match the HDF5 C library's H5D_FARRAY_FILT_COMPUTE_CHUNK_SIZE_LEN macro: - // chunk_size_len = 1 + ((H5VM_log2_gen(chunk.size) + 8) / 8) - // where chunk.size is the unfiltered chunk size in bytes (product of all chunk dims). - let chunk_size_bytes: usize = if has_filters { - let max_raw = chunks.iter().map(|c| c.raw_size).max().unwrap_or(1); - let log2_val = if max_raw <= 1 { - 0 - } else { - 63 - max_raw.leading_zeros() - }; - let len = 1 + ((log2_val + 8) / 8) as usize; - len.min(8) - } else { - 0 - }; - - let elem_size = if has_filters { - os + chunk_size_bytes + 4 - } else { - os - }; + let num_elements = slots.len(); + let chunk_size_bytes = has_filters.then(|| filtered_chunk_size_len(slots)); + let elem_size = os + chunk_size_bytes.map_or(0, |n| n + 4); let client_id: u8 = if has_filters { 1 } else { 0 }; // FAHD total size - let nelmts_field_size = length_size as usize; - let fahd_total_size = 4 + 1 + 1 + 1 + 1 + nelmts_field_size + os + 4; + let fahd_total_size = 4 + 1 + 1 + 1 + 1 + length_size as usize + os + 4; let fadb_address = fa_base_address + fahd_total_size as u64; - // Build FAHD let mut fahd = Vec::with_capacity(fahd_total_size); fahd.extend_from_slice(b"FAHD"); fahd.push(0); // version fahd.push(client_id); fahd.push(elem_size as u8); - - // max_nelmts_bits: use 10 as default (page_size = 1024), matching h5py convention - let max_bits: u8 = 10; - fahd.push(max_bits); - + fahd.push(FA_PAGE_BITS); match length_size { 4 => fahd.extend_from_slice(&(num_elements as u32).to_le_bytes()), - 8 => fahd.extend_from_slice(&(num_elements as u64).to_le_bytes()), _ => fahd.extend_from_slice(&(num_elements as u64).to_le_bytes()), } - - match offset_size { - 4 => fahd.extend_from_slice(&(fadb_address as u32).to_le_bytes()), - 8 => fahd.extend_from_slice(&fadb_address.to_le_bytes()), - _ => fahd.extend_from_slice(&fadb_address.to_le_bytes()), - } - - // Checksum + push_addr(&mut fahd, fadb_address, offset_size); let checksum = jenkins_lookup3(&fahd); fahd.extend_from_slice(&checksum.to_le_bytes()); - assert_eq!(fahd.len(), fahd_total_size); - // Build FADB + // FADB prefix let mut fadb = Vec::new(); fadb.extend_from_slice(b"FADB"); fadb.push(0); // version fadb.push(client_id); + push_addr(&mut fadb, fa_base_address, offset_size); - // header address - match offset_size { - 4 => fadb.extend_from_slice(&(fa_base_address as u32).to_le_bytes()), - 8 => fadb.extend_from_slice(&fa_base_address.to_le_bytes()), - _ => fadb.extend_from_slice(&fa_base_address.to_le_bytes()), - } - - // Element data - for chunk in chunks { - match offset_size { - 4 => fadb.extend_from_slice(&(chunk.address as u32).to_le_bytes()), - 8 => fadb.extend_from_slice(&chunk.address.to_le_bytes()), - _ => fadb.extend_from_slice(&chunk.address.to_le_bytes()), + let page_nelmts = 1usize << FA_PAGE_BITS; + if num_elements <= page_nelmts { + // Unpaged: the elements follow the prefix, one checksum over both. + for slot in slots { + push_index_element(&mut fadb, slot.as_ref(), offset_size, chunk_size_bytes); } - if has_filters { - // Write compressed size using chunk_size_bytes (variable width) - let cs_bytes = chunk.compressed_size.to_le_bytes(); - fadb.extend_from_slice(&cs_bytes[..chunk_size_bytes]); - fadb.extend_from_slice(&chunk.filter_mask.to_le_bytes()); + let fadb_checksum = jenkins_lookup3(&fadb); + fadb.extend_from_slice(&fadb_checksum.to_le_bytes()); + } else { + // Paged: every page is written, so every page-init bit is set + // (MSB-first, as `H5VM_bit_set` packs them). The prefix and bitmap + // share a checksum; each page carries its own. + let npages = num_elements.div_ceil(page_nelmts); + let mut bitmap = vec![0u8; npages.div_ceil(8)]; + for p in 0..npages { + bitmap[p / 8] |= 0x80 >> (p % 8); + } + fadb.extend_from_slice(&bitmap); + let prefix_checksum = jenkins_lookup3(&fadb); + fadb.extend_from_slice(&prefix_checksum.to_le_bytes()); + for page in slots.chunks(page_nelmts) { + let start = fadb.len(); + for slot in page { + push_index_element(&mut fadb, slot.as_ref(), offset_size, chunk_size_bytes); + } + let page_checksum = jenkins_lookup3(&fadb[start..]); + fadb.extend_from_slice(&page_checksum.to_le_bytes()); } } - // FADB checksum - let fadb_checksum = jenkins_lookup3(&fadb); - fadb.extend_from_slice(&fadb_checksum.to_le_bytes()); - let mut combined = fahd; combined.extend_from_slice(&fadb); combined @@ -734,8 +766,9 @@ pub fn build_chunked_data_from_precompressed( ) } else { let fa_address = base_address + data_buf.len() as u64; + let slots: Vec> = written_chunks.iter().cloned().map(Some).collect(); let fa_bytes = build_fixed_array_at( - &written_chunks, + &slots, offset_size, length_size, pre.has_filters, @@ -747,7 +780,7 @@ pub fn build_chunked_data_from_precompressed( fa_address, offset_size, element_size as u32, - 10, // max_nelmts_bits — matches h5py convention + FA_PAGE_BITS, ) }; @@ -1382,7 +1415,8 @@ mod tests { filter_mask: 0, }, ]; - let fa = build_fixed_array_at(&chunks, 8, 8, false, 0x2000); + let slots: Vec<_> = chunks.into_iter().map(Some).collect(); + let fa = build_fixed_array_at(&slots, 8, 8, false, 0x2000); // Should start with FAHD assert_eq!(&fa[0..4], b"FAHD"); // FAHD size = 4+1+1+1+1+8+8+4 = 28 diff --git a/crates/clawhdf5/tests/chunk_index_interop.rs b/crates/clawhdf5/tests/chunk_index_interop.rs index 9c030a3..cce2df6 100644 --- a/crates/clawhdf5/tests/chunk_index_interop.rs +++ b/crates/clawhdf5/tests/chunk_index_interop.rs @@ -11,7 +11,7 @@ use std::process::Command; -use clawhdf5::File; +use clawhdf5::{File, FileBuilder}; fn python() -> String { std::env::var("CLAWHDF5_PYTHON").unwrap_or_else(|_| "python3".to_string()) @@ -254,3 +254,139 @@ fn h5py_fixed_array_partial_extent_reads_correctly() { }, ]); } + +// =========================================================================== +// Files we write, read back by libhdf5 (h5py and h5dump) and by us +// =========================================================================== + +/// One `i4` dataset we write, filled with `arange` over `shape`. +struct WriteCase { + name: String, + shape: Vec, + chunks: Vec, + maxshape: Option>, + deflate: bool, +} + +fn wcase(name: &str, shape: &[u64], chunks: &[u64], maxshape: Option<&[u64]>) -> WriteCase { + WriteCase { + name: name.to_string(), + shape: shape.to_vec(), + chunks: chunks.to_vec(), + maxshape: maxshape.map(<[u64]>::to_vec), + deflate: false, + } +} + +fn h5dump_available() -> bool { + Command::new("h5dump") + .arg("--version") + .output() + .map(|o| o.status.success()) + .unwrap_or(false) +} + +/// Write every case into one file with our writer, then check that our own +/// reader, h5py and h5dump (when installed) all return every value. Only the +/// libhdf5 half is skipped without h5py. +fn check_we_write(cases: &[WriteCase]) { + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("ours_chunk_index.h5"); + let path_str = path.display().to_string(); + + let mut b = FileBuilder::new(); + for c in cases { + let n: u64 = c.shape.iter().product(); + let data: Vec = (0..n as i32).collect(); + let ds = b.create_dataset(&c.name); + ds.with_i32_data(&data) + .with_shape(&c.shape) + .with_chunks(&c.chunks); + if let Some(ms) = &c.maxshape { + ds.with_maxshape(ms); + } + if c.deflate { + ds.with_deflate(4); + } + } + b.write(&path).unwrap(); + + // Our reader. + let file = File::open(&path).unwrap(); + for c in cases { + let got = file.dataset(&c.name).unwrap().read_i32().unwrap(); + let n: u64 = c.shape.iter().product(); + let bad = got + .iter() + .enumerate() + .filter(|&(i, &v)| v != i as i32) + .count(); + assert!( + got.len() == n as usize && bad == 0, + "{}: our reader: {bad} of {n} values wrong", + c.name + ); + } + + // libhdf5 via h5py. + skip_if_no_python!(); + let mut script = + format!("import h5py, numpy as np\nbad = []\nf = h5py.File(r'{path_str}', 'r')\n"); + for c in cases { + let shape: Vec = c.shape.iter().map(u64::to_string).collect(); + let maxshape: Vec = c + .maxshape + .as_ref() + .unwrap_or(&c.shape) + .iter() + .map(|&d| { + if d == u64::MAX { + "None".to_string() + } else { + d.to_string() + } + }) + .collect(); + script += &format!( + "d = f['{name}']\n\ + want = np.arange({n}, dtype='i4').reshape(({shape},))\n\ + got = d[()]\n\ + if d.maxshape != ({maxshape},): bad.append(('{name}', 'maxshape', d.maxshape))\n\ + elif not np.array_equal(got, want): \ + bad.append(('{name}', int((got != want).sum()), 'of', got.size))\n", + name = c.name, + n = c.shape.iter().product::(), + shape = shape.join(","), + maxshape = maxshape.join(","), + ); + } + script += "print(bad if bad else 'OK')\n"; + let out = run_python(&script); + assert_eq!(out, "OK", "h5py disagrees"); + + // libhdf5's own tool, when installed. + if h5dump_available() { + let o = Command::new("h5dump").arg(&path).output().unwrap(); + let stderr = String::from_utf8_lossy(&o.stderr); + assert!( + o.status.success() && !stderr.to_lowercase().contains("error"), + "h5dump failed: {stderr}" + ); + } +} + +/// A Fixed Array with more than 1024 elements must be paged, or libhdf5 +/// rejects the data block's checksum. +#[test] +fn we_write_paged_fixed_array() { + let mut cases: Vec = [1023u64, 1024, 1025, 2048, 5000] + .iter() + .map(|&n| wcase(&format!("fa_{n}"), &[n * 4], &[4], None)) + .collect(); + // Filtered elements are wider; a 2-D grid pages the same way. + let mut filtered = wcase("fa_1500_deflate", &[1500 * 4], &[4], None); + filtered.deflate = true; + cases.push(filtered); + cases.push(wcase("fa_2d_1100", &[110, 40], &[1, 4], None)); + check_we_write(&cases); +} From 44f5f8b5c5bc0b4ba9960c2fa1aa286e2f630240 Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:09:29 -0500 Subject: [PATCH 3/7] fix(format): index every Extensible Array chunk, not just the first 244 The Extensible Array writer only filled the index block's 4 inline elements and the 6 data blocks it addresses directly (240 elements); its super block addresses were always undefined. Chunks from index 244 on were written to the file but never indexed, so they read back as fill values in our reader and in libhdf5, without an error. The writer now lays out data blocks and super blocks for any element count as H5EA__hdr_init sizes them, pages data blocks larger than 1024 elements (page-init bits in the owning super block), leaves blocks with no defined element unallocated, and records real header statistics (max_idx_set is one past the highest defined index). Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/chunked_write.rs | 8 +- crates/clawhdf5-format/src/ea_writer.rs | 513 +++++++++---------- crates/clawhdf5/tests/chunk_index_interop.rs | 21 + 3 files changed, 271 insertions(+), 271 deletions(-) diff --git a/crates/clawhdf5-format/src/chunked_write.rs b/crates/clawhdf5-format/src/chunked_write.rs index c0d02e1..66bc82d 100644 --- a/crates/clawhdf5-format/src/chunked_write.rs +++ b/crates/clawhdf5-format/src/chunked_write.rs @@ -503,7 +503,7 @@ fn serialize_v4_fixed_array( /// default, `H5D_FARRAY_MAX_DBLK_PAGE_NELMTS_BITS`). const FA_PAGE_BITS: u8 = 10; -fn push_addr(buf: &mut Vec, addr: u64, offset_size: u8) { +pub(crate) fn push_addr(buf: &mut Vec, addr: u64, offset_size: u8) { match offset_size { 4 => buf.extend_from_slice(&(addr as u32).to_le_bytes()), _ => buf.extend_from_slice(&addr.to_le_bytes()), @@ -734,8 +734,9 @@ pub fn build_chunked_data_from_precompressed( let layout_message = if use_extensible { let ea_address = base_address + data_buf.len() as u64; + let slots: Vec> = written_chunks.iter().cloned().map(Some).collect(); let ea_bytes = ea_writer::build_extensible_array_at( - &written_chunks, + &slots, offset_size, length_size, pre.has_filters, @@ -1463,7 +1464,8 @@ mod tests { filter_mask: 0, }, ]; - let ea = ea_writer::build_extensible_array_at(&chunks, 8, 8, false, 0x2000); + let slots: Vec<_> = chunks.into_iter().map(Some).collect(); + let ea = ea_writer::build_extensible_array_at(&slots, 8, 8, false, 0x2000); assert_eq!(&ea[0..4], b"EAHD"); // Find EAIB after EAHD: 12 fixed + 6*8 stats + 8 addr + 4 checksum = 72 let aehd_size = 4 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 6 * 8 + 8 + 4; diff --git a/crates/clawhdf5-format/src/ea_writer.rs b/crates/clawhdf5-format/src/ea_writer.rs index d92fca2..e0f0375 100644 --- a/crates/clawhdf5-format/src/ea_writer.rs +++ b/crates/clawhdf5-format/src/ea_writer.rs @@ -7,7 +7,7 @@ extern crate alloc; use alloc::{vec, vec::Vec}; use crate::checksum::jenkins_lookup3; -use crate::chunked_write::WrittenChunk; +use crate::chunked_write::{WrittenChunk, filtered_chunk_size_len, push_addr, push_index_element}; /// Serialize a v4 Extensible Array layout message. pub(crate) fn serialize_v4_extensible_array( @@ -58,11 +58,11 @@ pub(crate) fn serialize_v4_extensible_array( buf.push(4); // EA creation parameters (must match AEHD and HDF5 C library defaults) - buf.push(32); // max_nelmts_bits - buf.push(4); // idx_blk_elmts - buf.push(4); // super_blk_min_data_ptrs - buf.push(16); // data_blk_min_elmts - buf.push(10); // max_dblk_page_nelmts_bits + buf.push(MAX_NELMTS_BITS); + buf.push(IDX_BLK_ELMTS); + buf.push(SUP_BLK_MIN_DATA_PTRS); + buf.push(DATA_BLK_MIN_ELMTS); + buf.push(MAX_DBLK_PAGE_NELMTS_BITS); // EA header address match offset_size { @@ -74,304 +74,281 @@ pub(crate) fn serialize_v4_extensible_array( buf } +// EA creation parameters — the HDF5 library's defaults for chunk indexes +// (`H5D_EARRAY_*`); the layout message above and the header must agree. +const MAX_NELMTS_BITS: u8 = 32; +const IDX_BLK_ELMTS: u8 = 4; +const SUP_BLK_MIN_DATA_PTRS: u8 = 4; +const DATA_BLK_MIN_ELMTS: u8 = 16; +const MAX_DBLK_PAGE_NELMTS_BITS: u8 = 10; + +/// One data block of the array: its first element (relative to the end of +/// the index block's own elements), element count, and address when it is +/// allocated. +struct DataBlock { + start: usize, + nelmts: usize, + addr: Option, +} + /// Build a complete Extensible Array at a known absolute address. /// -/// For simplicity, we put all elements inline in the index block when the -/// number of chunks is small (up to idx_blk_elmts), otherwise use inline + -/// direct data blocks. +/// `slots[i]` is the element at linear index `i` (see `chunk_grid`); `None` +/// marks an unallocated chunk. The first `IDX_BLK_ELMTS` elements live in +/// the index block, the rest in data blocks grouped by super block level +/// exactly as `H5EA__hdr_init` sizes them: level `u` has `2^(u/2)` data +/// blocks of `DATA_BLK_MIN_ELMTS * 2^ceil(u/2)` elements. The data blocks of +/// the first levels are addressed straight from the index block; later +/// levels go through a super block (EASB). Data blocks larger than a page +/// (`2^MAX_DBLK_PAGE_NELMTS_BITS` elements) are paged, with their page-init +/// bits kept in the owning super block. Only blocks holding a defined element +/// are allocated; the rest keep the undefined address, as in a file the +/// library wrote. pub fn build_extensible_array_at( - chunks: &[WrittenChunk], + slots: &[Option], offset_size: u8, length_size: u8, has_filters: bool, ea_base_address: u64, ) -> Vec { let os = offset_size as usize; - let num_elements = chunks.len(); - - // Compute element encoding size (same logic as Fixed Array) - let chunk_size_bytes: usize = if has_filters { - let max_raw = chunks.iter().map(|c| c.raw_size).max().unwrap_or(1); - let log2_val = if max_raw <= 1 { - 0 - } else { - 63 - max_raw.leading_zeros() - }; - let len = 1 + ((log2_val + 8) / 8) as usize; - len.min(8) - } else { - 0 - }; - - let elem_size = if has_filters { - os + chunk_size_bytes + 4 - } else { - os - }; - + let chunk_size_bytes = has_filters.then(|| filtered_chunk_size_len(slots)); + let elem_size = os + chunk_size_bytes.map_or(0, |n| n + 4); let client_id: u8 = if has_filters { 1 } else { 0 }; + let arr_off_size = (MAX_NELMTS_BITS as usize).div_ceil(8); + let page_nelmts = 1usize << MAX_DBLK_PAGE_NELMTS_BITS; + let idx_blk = IDX_BLK_ELMTS as usize; - // EA creation parameters — must match HDF5 C library defaults exactly - let max_nelmts_bits: u8 = 32; - let idx_blk_elmts: u8 = 4; - let min_dblk_nelmts: u8 = 16; - let super_blk_min_nelmts: u8 = 4; - let max_dblk_nelmts_bits: u8 = 10; + // Elements past the last defined one are never realised + // (`max_idx_set` is one past the highest index ever set). + let max_idx_set = slots.iter().rposition(Option::is_some).map_or(0, |i| i + 1); + let slots = &slots[..max_idx_set]; + let defined_in = |start: usize, n: usize| -> bool { + let lo = idx_blk.saturating_add(start).min(slots.len()); + let hi = idx_blk + .saturating_add(start) + .saturating_add(n) + .min(slots.len()); + slots[lo..hi].iter().any(Option::is_some) + }; - // EAHD size: fixed(12) + 6 stats(6*length_size) + addr(offset_size) + checksum(4) + // Super block levels: (ndblks, dblk_nelmts, first element). + let log2_dmin = (DATA_BLK_MIN_ELMTS as u32).trailing_zeros() as usize; + let nsblks = 1 + MAX_NELMTS_BITS as usize - log2_dmin; + let ndblk_addrs = 2 * (SUP_BLK_MIN_DATA_PTRS as usize - 1); + let mut levels: Vec<(usize, usize, usize)> = Vec::with_capacity(nsblks); + let mut start = 0usize; + for u in 0..nsblks { + let ndblks = 1usize << (u / 2); + let nelmts = (DATA_BLK_MIN_ELMTS as usize) << u.div_ceil(2); + levels.push((ndblks, nelmts, start)); + // Saturate: on 32-bit targets the last levels only need to compare + // as "beyond the end". + start = start.saturating_add(ndblks.saturating_mul(nelmts)); + } + // Levels whose data blocks the index block addresses directly. + let mut direct_levels = 0; + let mut n = 0; + while n < ndblk_addrs { + n += levels[direct_levels].0; + direct_levels += 1; + } + let nsblk_addrs = nsblks - direct_levels; + + let dblk_size = |nelmts: usize| -> usize { + let prefix = 4 + 1 + 1 + os + arr_off_size + 4; + if nelmts > page_nelmts { + prefix + (nelmts / page_nelmts) * (page_nelmts * elem_size + 4) + } else { + prefix + nelmts * elem_size + } + }; + let sblk_bitmap_len = |ndblks: usize, nelmts: usize| -> usize { + if nelmts > page_nelmts { + ndblks * (nelmts / page_nelmts).div_ceil(8) + } else { + 0 + } + }; + + // Plan addresses: header, index block, the direct data blocks, then each + // allocated super block followed by its allocated data blocks. let aehd_size = 4 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 1 + 6 * length_size as usize + os + 4; let aeib_address = ea_base_address + aehd_size as u64; + let aeib_size = 4 + 1 + 1 + os + idx_blk * elem_size + ndblk_addrs * os + nsblk_addrs * os + 4; + let mut cursor = aeib_address + aeib_size as u64; - // Determine how many elements go inline vs data blocks - let n_inline = (idx_blk_elmts as usize).min(num_elements); - let remaining_after_inline = num_elements.saturating_sub(n_inline); + let mut ndata_blks = 0u64; + let mut data_blk_size = 0u64; + let mut nsuper_blks = 0u64; + let mut super_blk_size = 0u64; + let mut realized = idx_blk as u64; - // Compute super block layout per HDF5 spec - let sblk_min = super_blk_min_nelmts as usize; - let log2_dblk_min = if min_dblk_nelmts <= 1 { - 0 - } else { - (min_dblk_nelmts as u32).trailing_zeros() as usize + let mut plan_dblk = |cursor: &mut u64, start: usize, nelmts: usize| -> DataBlock { + let addr = defined_in(start, nelmts).then(|| { + let a = *cursor; + let size = dblk_size(nelmts) as u64; + *cursor += size; + ndata_blks += 1; + data_blk_size += size; + realized += nelmts as u64; + a + }); + DataBlock { + start, + nelmts, + addr, + } }; - let nsblks = (max_nelmts_bits as usize).saturating_sub(log2_dblk_min) + 1; - // Direct data block addresses (from super blocks 0..sblk_min-1) - let mut dblk_sizes: Vec = Vec::new(); - for sblk_idx in 0..sblk_min.min(nsblks) { - let ndblks = 1usize << (sblk_idx / 2); - let dblk_nelmts = (min_dblk_nelmts as usize) * (1 << sblk_idx.div_ceil(2)); - for _ in 0..ndblks { - dblk_sizes.push(dblk_nelmts); + let mut direct: Vec = Vec::with_capacity(ndblk_addrs); + for &(ndblks, nelmts, first) in &levels[..direct_levels] { + for k in 0..ndblks { + direct.push(plan_dblk(&mut cursor, first + k * nelmts, nelmts)); } } - let n_direct_dblks = dblk_sizes.len(); - - // Super block addresses (for super blocks sblk_min..nsblks-1) - let n_sblk_addrs = nsblks.saturating_sub(sblk_min); - - // EAIB size - let aeib_size = 4 - + 1 - + 1 - + os - + idx_blk_elmts as usize * elem_size - + n_direct_dblks * os - + n_sblk_addrs * os - + 4; - - // Build AEHD - let mut aehd = Vec::with_capacity(aehd_size); - aehd.extend_from_slice(b"EAHD"); - aehd.push(0); // version - aehd.push(client_id); - aehd.push(elem_size as u8); - aehd.push(max_nelmts_bits); - aehd.push(idx_blk_elmts); - aehd.push(min_dblk_nelmts); - aehd.push(super_blk_min_nelmts); - aehd.push(max_dblk_nelmts_bits); - - // Count data blocks that will have chunks - let n_active_dblks: u64 = if remaining_after_inline > 0 { - let mut count = 0u64; - let mut ci = n_inline; - for &sz in &dblk_sizes { - if ci < num_elements { - count += 1; - ci += sz; - } + // (super block address, level, its data blocks) + let mut supers: Vec<(Option, usize, Vec)> = Vec::with_capacity(nsblk_addrs); + for (u, &(ndblks, nelmts, first)) in levels.iter().enumerate().skip(direct_levels) { + if !defined_in(first, ndblks.saturating_mul(nelmts)) { + supers.push((None, u, Vec::new())); + continue; } - count - } else { - 0 - }; - let blk_off_size = (max_nelmts_bits as usize).div_ceil(8); - let aedb_header_overhead = 4 + 1 + 1 + os + blk_off_size + 4; - let data_blk_total_size: u64 = if remaining_after_inline > 0 { - let mut total = 0u64; - let mut ci = n_inline; - for &sz in &dblk_sizes { - if ci < num_elements { - total += (aedb_header_overhead + sz * elem_size) as u64; - ci += sz; - } - } - total - } else { - 0 - }; - let max_idx_set: u64 = if remaining_after_inline > 0 { - let mut max_set = idx_blk_elmts as u64; - let mut ci = n_inline; - for &sz in &dblk_sizes { - if ci < num_elements { - max_set += sz as u64; - ci += sz; - } - } - max_set - } else { - idx_blk_elmts as u64 - }; + let sb_size = + 4 + 1 + 1 + os + arr_off_size + sblk_bitmap_len(ndblks, nelmts) + ndblks * os + 4; + let sb_addr = cursor; + cursor += sb_size as u64; + nsuper_blks += 1; + super_blk_size += sb_size as u64; + let dblks = (0..ndblks) + .map(|k| plan_dblk(&mut cursor, first + k * nelmts, nelmts)) + .collect(); + supers.push((Some(sb_addr), u, dblks)); + } + let slot = |i: usize| slots.get(i).and_then(Option::as_ref); let write_length = |buf: &mut Vec, val: u64| match length_size { 4 => buf.extend_from_slice(&(val as u32).to_le_bytes()), _ => buf.extend_from_slice(&val.to_le_bytes()), }; - let write_addr = |buf: &mut Vec, val: u64| match offset_size { - 4 => buf.extend_from_slice(&(val as u32).to_le_bytes()), - _ => buf.extend_from_slice(&val.to_le_bytes()), + let write_addr_opt = |buf: &mut Vec, addr: Option| match addr { + Some(a) => push_addr(buf, a, offset_size), + None => buf.extend(core::iter::repeat_n(0xFF, os)), + }; + let block_prefix = |buf: &mut Vec, sig: &[u8; 4], block_off: usize| { + buf.extend_from_slice(sig); + buf.push(0); // version + buf.push(client_id); + push_addr(buf, ea_base_address, offset_size); + buf.extend_from_slice(&(block_off as u64).to_le_bytes()[..arr_off_size]); + }; + // Serialise one data block (paged or not) onto `out`. + let write_dblk = |out: &mut Vec, db: &DataBlock| { + let at = out.len(); + block_prefix(out, b"EADB", db.start); + let first = idx_blk + db.start; + if db.nelmts > page_nelmts { + // Paged: the prefix carries only its own checksum; each page + // follows with one of its own. + let sum = jenkins_lookup3(&out[at..]); + out.extend_from_slice(&sum.to_le_bytes()); + for p in 0..db.nelmts / page_nelmts { + let page_at = out.len(); + for e in 0..page_nelmts { + let i = first + p * page_nelmts + e; + push_index_element(out, slot(i), offset_size, chunk_size_bytes); + } + let sum = jenkins_lookup3(&out[page_at..]); + out.extend_from_slice(&sum.to_le_bytes()); + } + } else { + for i in first..first + db.nelmts { + push_index_element(out, slot(i), offset_size, chunk_size_bytes); + } + let sum = jenkins_lookup3(&out[at..]); + out.extend_from_slice(&sum.to_le_bytes()); + } + debug_assert_eq!(out.len() - at, dblk_size(db.nelmts)); }; - write_length(&mut aehd, 0); - write_length(&mut aehd, 0); - write_length(&mut aehd, n_active_dblks); - write_length(&mut aehd, data_blk_total_size); - write_length(&mut aehd, num_elements as u64); - write_length(&mut aehd, max_idx_set); + // Header (EAHD). The six statistics are, in order: super blocks, their + // bytes, data blocks, their bytes, max index set, elements realised. + let mut out = Vec::with_capacity((cursor - ea_base_address) as usize); + out.extend_from_slice(b"EAHD"); + out.push(0); // version + out.push(client_id); + out.push(elem_size as u8); + out.push(MAX_NELMTS_BITS); + out.push(IDX_BLK_ELMTS); + out.push(DATA_BLK_MIN_ELMTS); + out.push(SUP_BLK_MIN_DATA_PTRS); + out.push(MAX_DBLK_PAGE_NELMTS_BITS); + write_length(&mut out, nsuper_blks); + write_length(&mut out, super_blk_size); + write_length(&mut out, ndata_blks); + write_length(&mut out, data_blk_size); + write_length(&mut out, max_idx_set as u64); + write_length(&mut out, realized); + push_addr(&mut out, aeib_address, offset_size); + let sum = jenkins_lookup3(&out); + out.extend_from_slice(&sum.to_le_bytes()); + debug_assert_eq!(out.len(), aehd_size); - write_addr(&mut aehd, aeib_address); - - let aehd_checksum = jenkins_lookup3(&aehd); - aehd.extend_from_slice(&aehd_checksum.to_le_bytes()); - debug_assert_eq!(aehd.len(), aehd_size); - - // Build AEIB - let mut aeib = Vec::with_capacity(aeib_size); - aeib.extend_from_slice(b"EAIB"); - aeib.push(0); - aeib.push(client_id); - - match offset_size { - 4 => aeib.extend_from_slice(&(ea_base_address as u32).to_le_bytes()), - 8 => aeib.extend_from_slice(&ea_base_address.to_le_bytes()), - _ => aeib.extend_from_slice(&ea_base_address.to_le_bytes()), + // Index block (EAIB): inline elements, data block and super block + // addresses. + let ib_start = out.len(); + out.extend_from_slice(b"EAIB"); + out.push(0); + out.push(client_id); + push_addr(&mut out, ea_base_address, offset_size); + for i in 0..idx_blk { + push_index_element(&mut out, slot(i), offset_size, chunk_size_bytes); } - - // Inline elements - #[allow(clippy::needless_range_loop)] - for i in 0..idx_blk_elmts as usize { - if i < n_inline { - write_chunk_element( - &mut aeib, - &chunks[i], - offset_size, - has_filters, - chunk_size_bytes, - ); - } else { - write_undefined_element(&mut aeib, offset_size, has_filters, chunk_size_bytes); - } + for db in &direct { + write_addr_opt(&mut out, db.addr); } + for (sb_addr, _, _) in &supers { + write_addr_opt(&mut out, *sb_addr); + } + let sum = jenkins_lookup3(&out[ib_start..]); + out.extend_from_slice(&sum.to_le_bytes()); + debug_assert_eq!(out.len() - ib_start, aeib_size); - // Data block addresses + build data blocks - let mut data_blocks_buf = Vec::new(); - let dblks_base = aeib_address + aeib_size as u64; - let mut dblk_cursor = dblks_base; - let mut chunk_idx = n_inline; - - for &nelmts in &dblk_sizes { - if chunk_idx >= num_elements { - match offset_size { - 4 => aeib.extend_from_slice(&u32::MAX.to_le_bytes()), - 8 => aeib.extend_from_slice(&u64::MAX.to_le_bytes()), - _ => aeib.extend_from_slice(&u64::MAX.to_le_bytes()), - } + for db in direct.iter().filter(|d| d.addr.is_some()) { + write_dblk(&mut out, db); + } + for (sb_addr, u, dblks) in &supers { + if sb_addr.is_none() { continue; } - - match offset_size { - 4 => aeib.extend_from_slice(&(dblk_cursor as u32).to_le_bytes()), - 8 => aeib.extend_from_slice(&dblk_cursor.to_le_bytes()), - _ => aeib.extend_from_slice(&dblk_cursor.to_le_bytes()), - } - - // Build EADB - let mut aedb = Vec::new(); - aedb.extend_from_slice(b"EADB"); - aedb.push(0); - aedb.push(client_id); - match offset_size { - 4 => aedb.extend_from_slice(&(ea_base_address as u32).to_le_bytes()), - 8 => aedb.extend_from_slice(&ea_base_address.to_le_bytes()), - _ => aedb.extend_from_slice(&ea_base_address.to_le_bytes()), - } - - let blk_off_size = (max_nelmts_bits as usize).div_ceil(8); - let blk_off_val = (chunk_idx - n_inline) as u64; - aedb.extend_from_slice(&blk_off_val.to_le_bytes()[..blk_off_size]); - - for slot in 0..nelmts { - if chunk_idx + slot < num_elements { - write_chunk_element( - &mut aedb, - &chunks[chunk_idx + slot], - offset_size, - has_filters, - chunk_size_bytes, - ); - } else { - write_undefined_element(&mut aedb, offset_size, has_filters, chunk_size_bytes); + let (ndblks, nelmts, first) = levels[*u]; + let sb_start = out.len(); + block_prefix(&mut out, b"EASB", first); + if nelmts > page_nelmts { + // Page-init bits, `npages` per data block, packed MSB-first + // (`H5VM_bit_set`): every page of an allocated data block is + // written. + let npages = nelmts / page_nelmts; + let mut bitmap = vec![0u8; sblk_bitmap_len(ndblks, nelmts)]; + for (k, db) in dblks.iter().enumerate() { + if db.addr.is_some() { + for p in 0..npages { + let bit = k * npages + p; + bitmap[bit / 8] |= 0x80 >> (bit % 8); + } + } } + out.extend_from_slice(&bitmap); } - - let aedb_checksum = jenkins_lookup3(&aedb); - aedb.extend_from_slice(&aedb_checksum.to_le_bytes()); - - dblk_cursor += aedb.len() as u64; - data_blocks_buf.extend_from_slice(&aedb); - chunk_idx += nelmts; - } - - // Super block addresses (all undefined) - for _ in 0..n_sblk_addrs { - match offset_size { - 4 => aeib.extend_from_slice(&u32::MAX.to_le_bytes()), - 8 => aeib.extend_from_slice(&u64::MAX.to_le_bytes()), - _ => aeib.extend_from_slice(&u64::MAX.to_le_bytes()), + for db in dblks { + write_addr_opt(&mut out, db.addr); + } + let sum = jenkins_lookup3(&out[sb_start..]); + out.extend_from_slice(&sum.to_le_bytes()); + for db in dblks.iter().filter(|d| d.addr.is_some()) { + write_dblk(&mut out, db); } } - - let aeib_checksum = jenkins_lookup3(&aeib); - aeib.extend_from_slice(&aeib_checksum.to_le_bytes()); - debug_assert_eq!(aeib.len(), aeib_size); - - let mut combined = aehd; - combined.extend_from_slice(&aeib); - combined.extend_from_slice(&data_blocks_buf); - combined -} - -fn write_chunk_element( - buf: &mut Vec, - chunk: &WrittenChunk, - offset_size: u8, - has_filters: bool, - chunk_size_bytes: usize, -) { - match offset_size { - 4 => buf.extend_from_slice(&(chunk.address as u32).to_le_bytes()), - 8 => buf.extend_from_slice(&chunk.address.to_le_bytes()), - _ => buf.extend_from_slice(&chunk.address.to_le_bytes()), - } - if has_filters { - let cs_bytes = chunk.compressed_size.to_le_bytes(); - buf.extend_from_slice(&cs_bytes[..chunk_size_bytes]); - buf.extend_from_slice(&chunk.filter_mask.to_le_bytes()); - } -} - -fn write_undefined_element( - buf: &mut Vec, - offset_size: u8, - has_filters: bool, - chunk_size_bytes: usize, -) { - let os = offset_size as usize; - // Use extend with repeat to avoid heap-allocating a temporary Vec on each call. - buf.extend(core::iter::repeat_n(0xFF, os)); - if has_filters { - buf.extend(core::iter::repeat_n(0x00, chunk_size_bytes)); - buf.extend_from_slice(&0u32.to_le_bytes()); - } + debug_assert_eq!(out.len() as u64, cursor - ea_base_address); + out } diff --git a/crates/clawhdf5/tests/chunk_index_interop.rs b/crates/clawhdf5/tests/chunk_index_interop.rs index cce2df6..0e7c845 100644 --- a/crates/clawhdf5/tests/chunk_index_interop.rs +++ b/crates/clawhdf5/tests/chunk_index_interop.rs @@ -390,3 +390,24 @@ fn we_write_paged_fixed_array() { cases.push(wcase("fa_2d_1100", &[110, 40], &[1, 4], None)); check_we_write(&cases); } + +/// An Extensible Array holds 4 elements in its index block and 240 in the +/// data blocks the index block addresses; everything after that lives under +/// super blocks, and from ~131K elements on in paged data blocks. Chunks past +/// index 243 used to be written but never indexed (read back as fill by us +/// and by libhdf5). +#[test] +fn we_write_extensible_array_past_index_block() { + let unl: &[u64] = &[u64::MAX]; + let mut cases: Vec = [1u64, 4, 5, 243, 244, 245, 300, 1000, 5000] + .iter() + .map(|&n| wcase(&format!("ea_{n}"), &[n * 4], &[4], Some(unl))) + .collect(); + let mut filtered = wcase("ea_300_deflate", &[300 * 4], &[4], Some(unl)); + filtered.deflate = true; + cases.push(filtered); + // Several super blocks and paged data blocks (level 13, the first with + // data blocks over 1024 elements, starts at element 4 + 131056). + cases.push(wcase("ea_140000", &[140_000], &[1], Some(unl))); + check_we_write(&cases); +} From 540fa08907913df9371e59d65f95dba0d357b2f0 Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:12:47 -0500 Subject: [PATCH 4/7] fix(format): write chunk indexes over the max extent, swizzled for EA The writer indexed chunks by their position in the current shape, the same mistake the reader had. With a finite maxshape larger than the shape the Fixed Array was sized for the shape, so libhdf5 looked up chunks past its end ("addr overflow"); with the unlimited dimension anywhere but first, e.g. maxshape (20, None), libhdf5 swizzles that dimension to the slowest position and read our Extensible Array scrambled. Two unlimited dimensions produced a file libhdf5 refused to open ("already found unlimited dimension"). Chunks are now placed with the shared chunk_grid linearisation: Fixed Array slots cover every chunk of the maximum extent (unwritten ones undefined), Extensible Array indexes are swizzled, Single Chunk is only used when the maximum extent is one chunk, and a maxshape that is smaller than the shape, has more than one unlimited dimension, or would need an absurd Fixed Array is an error instead of a bad file. build_chunked_data_from_precompressed now returns a Result. Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/chunk_grid.rs | 1 - crates/clawhdf5-format/src/chunked_write.rs | 137 +++++++++++++++++-- crates/clawhdf5-format/src/file_writer.rs | 4 +- crates/clawhdf5/tests/chunk_index_interop.rs | 117 ++++++++++++++++ 4 files changed, 243 insertions(+), 16 deletions(-) diff --git a/crates/clawhdf5-format/src/chunk_grid.rs b/crates/clawhdf5-format/src/chunk_grid.rs index c5a06fd..9e03b9d 100644 --- a/crates/clawhdf5-format/src/chunk_grid.rs +++ b/crates/clawhdf5-format/src/chunk_grid.rs @@ -147,7 +147,6 @@ impl ChunkGrid { /// Linear index of the chunk with scaled coordinates `scaled` /// (`offset / chunk_dim` per dimension, in dataset order). - #[allow(dead_code)] // used by the writer pub(crate) fn linear_index(&self, scaled: &[u64]) -> u64 { self.order .iter() diff --git a/crates/clawhdf5-format/src/chunked_write.rs b/crates/clawhdf5-format/src/chunked_write.rs index 66bc82d..3136f58 100644 --- a/crates/clawhdf5-format/src/chunked_write.rs +++ b/crates/clawhdf5-format/src/chunked_write.rs @@ -8,6 +8,7 @@ use alloc::{vec, vec::Vec}; use crate::checksum::jenkins_lookup3; use crate::chunk_cache::{CACHE_LINE_SIZE, align_to_cache_line}; +use crate::chunk_grid::ChunkGrid; use crate::ea_writer; use crate::error::FormatError; use crate::filter_pipeline::{ @@ -699,7 +700,8 @@ pub fn build_chunked_data_from_precompressed( pre: &PrecompressedChunks, base_address: u64, maxshape: Option<&[u64]>, -) -> ChunkedDataResult { +) -> Result { + let index = ChunkIndexPlan::new(&pre.shape, maxshape, &pre.chunk_dims)?; let offset_size: u8 = 8; let length_size: u8 = 8; let num_chunks = pre.chunks.len(); @@ -725,16 +727,15 @@ pub fn build_chunked_data_from_precompressed( } let chunk_dims_u32: Vec = pre.chunk_dims.iter().map(|&d| d as u32).collect(); - let use_extensible = maxshape.is_some_and(|ms| ms.contains(&u64::MAX)); let aligned_idx = align_to_cache_line(data_buf.len()); if aligned_idx > data_buf.len() { data_buf.resize(aligned_idx, 0u8); } - let layout_message = if use_extensible { + let layout_message = if let ChunkIndexPlan::ExtensibleArray(grid) = &index { let ea_address = base_address + data_buf.len() as u64; - let slots: Vec> = written_chunks.iter().cloned().map(Some).collect(); + let slots = index_slots(grid, &pre.shape, &pre.chunk_dims, &written_chunks, None)?; let ea_bytes = ea_writer::build_extensible_array_at( &slots, offset_size, @@ -749,7 +750,7 @@ pub fn build_chunked_data_from_precompressed( offset_size, element_size as u32, ) - } else if num_chunks == 1 { + } else if matches!(index, ChunkIndexPlan::SingleChunk) { let chunk_addr = written_chunks[0].address; let filtered_size = if pre.has_filters { Some(written_chunks[0].compressed_size) @@ -765,9 +766,15 @@ pub fn build_chunked_data_from_precompressed( offset_size, element_size as u32, ) - } else { + } else if let ChunkIndexPlan::FixedArray(grid, nslots) = &index { let fa_address = base_address + data_buf.len() as u64; - let slots: Vec> = written_chunks.iter().cloned().map(Some).collect(); + let slots = index_slots( + grid, + &pre.shape, + &pre.chunk_dims, + &written_chunks, + Some(*nslots), + )?; let fa_bytes = build_fixed_array_at( &slots, offset_size, @@ -783,15 +790,123 @@ pub fn build_chunked_data_from_precompressed( element_size as u32, FA_PAGE_BITS, ) + } else { + unreachable!("every chunk index plan is handled above") }; - ChunkedDataResult { + Ok(ChunkedDataResult { data_bytes: data_buf, layout_message, pipeline_message: pre.pipeline_message.clone(), + }) +} + +/// Most slots a Fixed Array index may have before we refuse to build it: its +/// data block holds one element per chunk of the *maximum* extent, so a huge +/// finite maxshape with small chunks would otherwise exhaust memory. +const MAX_FIXED_ARRAY_SLOTS: u64 = 1 << 26; + +/// Which chunk index a dataset gets, following the library's choice in +/// `H5D__layout_set_latest_indexing`: Extensible Array for exactly one +/// unlimited dimension, Fixed Array for a finite maxshape, Single Chunk when +/// the whole maximum extent is one chunk. +enum ChunkIndexPlan { + SingleChunk, + /// The grid and the number of array elements (chunks of the max extent). + FixedArray(ChunkGrid, usize), + ExtensibleArray(ChunkGrid), +} + +impl ChunkIndexPlan { + fn new( + shape: &[u64], + maxshape: Option<&[u64]>, + chunk_dims: &[u64], + ) -> Result { + let bad = |what: &str| FormatError::ChunkedReadError(format!("maxshape: {what}")); + if let Some(ms) = maxshape { + if ms.len() != shape.len() { + return Err(bad("rank differs from the shape")); + } + if ms.iter().zip(shape).any(|(&m, &s)| m < s) { + return Err(bad("smaller than the shape")); + } + } + let max = maxshape.unwrap_or(shape); + let nunlim = max.iter().filter(|&&d| d == u64::MAX).count(); + match nunlim { + 0 => { + let nslots = max + .iter() + .zip(chunk_dims) + .try_fold(1u64, |acc, (&m, &c)| acc.checked_mul(m.div_ceil(c.max(1)))) + .filter(|&n| n <= MAX_FIXED_ARRAY_SLOTS) + .ok_or_else(|| { + bad("too many chunks for a Fixed Array index; \ + use larger chunks or an unlimited dimension") + })?; + // A Single Chunk index needs that one chunk to exist; an + // empty dataset gets an all-unallocated Fixed Array instead. + let empty = shape.contains(&0); + if nslots == 1 && !empty { + Ok(Self::SingleChunk) + } else { + let grid = ChunkGrid::fixed_array(shape, Some(max), chunk_dims)?; + Ok(Self::FixedArray(grid, nslots as usize)) + } + } + 1 => Ok(Self::ExtensibleArray(ChunkGrid::extensible_array( + shape, + Some(max), + chunk_dims, + )?)), + _ => Err(bad( + "more than one unlimited dimension needs a B-tree v2 chunk index, \ + which the writer does not support", + )), + } } } +/// Place each written chunk at its linear index in `grid`. `chunks` are in +/// row-major order over the chunks of the current extent (`split_into_chunks`). +/// `len` fixes the slot count (Fixed Array); otherwise it is one past the +/// highest index used. +fn index_slots( + grid: &ChunkGrid, + shape: &[u64], + chunk_dims: &[u64], + chunks: &[WrittenChunk], + len: Option, +) -> Result>, FormatError> { + let rank = shape.len(); + let cur: Vec = shape + .iter() + .zip(chunk_dims) + .map(|(&s, &c)| s.div_ceil(c)) + .collect(); + let mut placed: Vec<(usize, &WrittenChunk)> = Vec::with_capacity(chunks.len()); + let mut scaled = vec![0u64; rank]; + for (i, chunk) in chunks.iter().enumerate() { + let mut rem = i as u64; + for d in (0..rank).rev() { + scaled[d] = rem % cur[d]; + rem /= cur[d]; + } + let idx = usize::try_from(grid.linear_index(&scaled)) + .map_err(|_| FormatError::Overflow("chunk index slot".into()))?; + placed.push((idx, chunk)); + } + let n = len.unwrap_or_else(|| placed.iter().map(|&(i, _)| i + 1).max().unwrap_or(0)); + let mut slots = vec![None; n]; + for (idx, chunk) in placed { + *slots + .get_mut(idx) + .ok_or_else(|| FormatError::Overflow("chunk index slot".into()))? = Some(chunk.clone()); + } + Ok(slots) +} + /// Build chunked data with absolute addresses. /// If `maxshape` has unlimited dims, uses Extensible Array index. pub fn build_chunked_data_at( @@ -824,11 +939,7 @@ pub fn build_chunked_data_at_ext( maxshape: Option<&[u64]>, ) -> Result { let pre = precompress_chunks(raw_data, shape, chunk_dims, element_size, options)?; - Ok(build_chunked_data_from_precompressed( - &pre, - base_address, - maxshape, - )) + build_chunked_data_from_precompressed(&pre, base_address, maxshape) } /// Write selected elements into an existing in-memory dataset buffer. diff --git a/crates/clawhdf5-format/src/file_writer.rs b/crates/clawhdf5-format/src/file_writer.rs index c262cd4..5350c88 100644 --- a/crates/clawhdf5-format/src/file_writer.rs +++ b/crates/clawhdf5-format/src/file_writer.rs @@ -1244,7 +1244,7 @@ impl FileWriter { &pre, dummy_cursor, d.maxshape.as_deref(), - ); + )?; dummy_cursor += result.data_bytes.len() as u64; let dense_blob = if ds_dense[i] { Some(build_dense_attrs(&d.attrs, 0)) @@ -1424,7 +1424,7 @@ impl FileWriter { .expect("chunked dataset missing precompressed cache"), base_address, d.maxshape.as_deref(), - ); + )?; cursor2 += result.data_bytes.len(); let oh = build_chunked_dataset_oh( &d.dt, diff --git a/crates/clawhdf5/tests/chunk_index_interop.rs b/crates/clawhdf5/tests/chunk_index_interop.rs index 0e7c845..41d7a5a 100644 --- a/crates/clawhdf5/tests/chunk_index_interop.rs +++ b/crates/clawhdf5/tests/chunk_index_interop.rs @@ -373,6 +373,73 @@ fn check_we_write(cases: &[WriteCase]) { "h5dump failed: {stderr}" ); } + + // Let libhdf5 grow every resizable dataset by two chunks per dimension + // (capped at the maxshape) and rewrite it, which updates our index in + // place and inserts new chunks into it. Then both readers must agree. + let script = format!( + r#" +import h5py, numpy as np +grown = {{}} +with h5py.File(r'{path_str}', 'r+') as f: + for name in f: + d = f[name] + if d.chunks is None: + continue + new = tuple(s + 2 * c if m is None else min(m, s + 2 * c) + for s, m, c in zip(d.shape, d.maxshape, d.chunks)) + if new == d.shape: + continue + old = d[()] + full = np.full(new, -7, 'i4') + full[tuple(slice(0, s) for s in old.shape)] = old + d.resize(new) + d[...] = full + grown[name] = (list(old.shape), list(new)) +with h5py.File(r'{path_str}', 'r') as f: + for name, (old, new) in grown.items(): + want = np.full(new, -7, 'i4') + want[tuple(slice(0, s) for s in old)] = np.arange(int(np.prod(old)), dtype='i4').reshape(old) + assert np.array_equal(f[name][()], want), name +for name, (old, new) in grown.items(): + print(name, ','.join(map(str, old)), ','.join(map(str, new))) +"# + ); + let out = run_python(&script); + let growable = cases + .iter() + .filter(|c| c.maxshape.as_ref().is_some_and(|m| *m != c.shape)) + .count(); + assert_eq!(out.lines().count(), growable, "libhdf5 grew: {out}"); + let dims = |s: &str| -> Vec { s.split(',').map(|x| x.parse().unwrap()).collect() }; + let file = File::open(&path).unwrap(); + for line in out.lines() { + let mut parts = line.split(' '); + let (name, old, new) = ( + parts.next().unwrap(), + dims(parts.next().unwrap()), + dims(parts.next().unwrap()), + ); + let got = file.dataset(name).unwrap().read_i32().unwrap(); + let n: usize = new.iter().product(); + let mut want = vec![-7i32; n]; + for (flat, w) in want.iter_mut().enumerate() { + let mut rem = flat; + let mut coords = vec![0usize; new.len()]; + for d in (0..new.len()).rev() { + coords[d] = rem % new[d]; + rem /= new[d]; + } + if coords.iter().zip(&old).all(|(c, o)| c < o) { + *w = coords.iter().zip(&old).fold(0, |acc, (c, o)| acc * o + c) as i32; + } + } + let bad = got.iter().zip(&want).filter(|(a, b)| a != b).count(); + assert!( + got.len() == n && bad == 0, + "{name}: after libhdf5 grew it, our reader got {bad} of {n} values wrong" + ); + } } /// A Fixed Array with more than 1024 elements must be paged, or libhdf5 @@ -411,3 +478,53 @@ fn we_write_extensible_array_past_index_block() { cases.push(wcase("ea_140000", &[140_000], &[1], Some(unl))); check_we_write(&cases); } + +/// A maxshape larger than the shape: the index must be laid out over the +/// chunks of the maximum extent (libhdf5 read our Fixed Array past its end: +/// "addr overflow"), and an Extensible Array whose unlimited dimension is not +/// the first must swizzle it to the slowest position (libhdf5 read our +/// `(20, None)` dataset scrambled). +#[test] +fn we_write_maxshape_larger_than_shape() { + const U: u64 = u64::MAX; + let mut cases = vec![ + // Fixed Array over the maximum extent. + wcase("fa2d_finite_max", &[20, 30], &[5, 5], Some(&[40, 60])), + wcase("fa1d_finite_max", &[40], &[4], Some(&[100])), + wcase("fa3d_edges", &[6, 7, 8], &[4, 3, 5], Some(&[10, 9, 20])), + wcase("fa_paged_max", &[30, 50], &[1, 1], Some(&[40, 60])), + wcase("fa_one_chunk_now", &[5], &[5], Some(&[50])), + // Extensible Array, unlimited dimension first (no swizzle) ... + wcase("ea2d_unl_fin", &[20, 30], &[5, 5], Some(&[U, 30])), + wcase("ea2d_unl_fin_max", &[20, 30], &[5, 5], Some(&[U, 60])), + // ... and not first (swizzled). + wcase("ea2d_fin_unl", &[20, 30], &[5, 5], Some(&[20, U])), + wcase("ea2d_fin_max_unl", &[20, 30], &[5, 5], Some(&[40, U])), + wcase("ea3d_mid", &[6, 7, 8], &[4, 3, 5], Some(&[10, U, 20])), + // Past the index block and into super blocks, swizzled. + wcase("ea2d_many", &[3, 2000], &[1, 1], Some(&[4, U])), + ]; + let mut filtered = wcase( + "ea3d_last_deflate", + &[6, 7, 8], + &[4, 3, 5], + Some(&[6, 8, U]), + ); + filtered.deflate = true; + cases.push(filtered); + check_we_write(&cases); +} + +/// More than one unlimited dimension needs a B-tree v2 chunk index; the +/// writer must not produce a file libhdf5 cannot open. +#[test] +fn two_unlimited_dims_are_refused() { + let mut b = FileBuilder::new(); + b.create_dataset("d") + .with_i32_data(&(0..600).collect::>()) + .with_shape(&[20, 30]) + .with_chunks(&[5, 5]) + .with_maxshape(&[u64::MAX, u64::MAX]); + let dir = tempfile::tempdir().unwrap(); + assert!(b.write(dir.path().join("unl_unl.h5")).is_err()); +} From 1dba7b465a95a8f2738b1a17241f99e6f17fbe6f Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:14:47 -0500 Subject: [PATCH 5/7] fix(format): index datasets with several unlimited dims by B-tree v2 A dataset with more than one unlimited dimension got an Extensible Array index, which libhdf5 refuses ("already found unlimited dimension"), so the whole file failed to open in h5py and h5dump. The previous commit turned that into a write error; this one writes what the library itself uses there: a version-2 B-tree chunk index (record type 10/11), as a single leaf of the library's 2048-byte node size, or a larger leaf when the records do not fit. The root's record count is 16-bit, so more than 65535 chunks is still refused rather than written wrong. Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/chunked_write.rs | 334 ++++++++++++++----- crates/clawhdf5/tests/chunk_index_interop.rs | 33 +- 2 files changed, 268 insertions(+), 99 deletions(-) diff --git a/crates/clawhdf5-format/src/chunked_write.rs b/crates/clawhdf5-format/src/chunked_write.rs index 3136f58..6ff0c07 100644 --- a/crates/clawhdf5-format/src/chunked_write.rs +++ b/crates/clawhdf5-format/src/chunked_write.rs @@ -444,6 +444,27 @@ fn serialize_v4_fixed_array( element_size: u32, max_bits: u8, ) -> Vec { + let mut buf = layout_v4_chunked_prefix(chunk_dims, element_size); + + // chunk index type = 3 (Fixed Array) + buf.push(3); + + // max_dblk_page_nelmts_bits — must match FAHD max_nelmts_bits + buf.push(max_bits); + + // Fixed Array header address + match offset_size { + 4 => buf.extend_from_slice(&(fixed_array_address as u32).to_le_bytes()), + 8 => buf.extend_from_slice(&fixed_array_address.to_le_bytes()), + _ => {} + } + + buf +} + +/// The part of a v4 chunked layout message before the chunk index type: +/// version, class, flags and the chunk dimensions (plus the element size). +fn layout_v4_chunked_prefix(chunk_dims: &[u32], element_size: u32) -> Vec { let mut buf = Vec::new(); buf.push(4); // version buf.push(2); // class = chunked @@ -483,20 +504,6 @@ fn serialize_v4_fixed_array( 4 => buf.extend_from_slice(&element_size.to_le_bytes()), _ => {} } - - // chunk index type = 3 (Fixed Array) - buf.push(3); - - // max_dblk_page_nelmts_bits — must match FAHD max_nelmts_bits - buf.push(max_bits); - - // Fixed Array header address - match offset_size { - 4 => buf.extend_from_slice(&(fixed_array_address as u32).to_le_bytes()), - 8 => buf.extend_from_slice(&fixed_array_address.to_le_bytes()), - _ => {} - } - buf } @@ -733,65 +740,91 @@ pub fn build_chunked_data_from_precompressed( data_buf.resize(aligned_idx, 0u8); } - let layout_message = if let ChunkIndexPlan::ExtensibleArray(grid) = &index { - let ea_address = base_address + data_buf.len() as u64; - let slots = index_slots(grid, &pre.shape, &pre.chunk_dims, &written_chunks, None)?; - let ea_bytes = ea_writer::build_extensible_array_at( - &slots, - offset_size, - length_size, - pre.has_filters, - ea_address, - ); - data_buf.extend_from_slice(&ea_bytes); - ea_writer::serialize_v4_extensible_array( - &chunk_dims_u32, - ea_address, - offset_size, - element_size as u32, - ) - } else if matches!(index, ChunkIndexPlan::SingleChunk) { - let chunk_addr = written_chunks[0].address; - let filtered_size = if pre.has_filters { - Some(written_chunks[0].compressed_size) - } else { - None - }; - let filter_mask = if pre.has_filters { Some(0u32) } else { None }; - serialize_v4_single_chunk( - &chunk_dims_u32, - chunk_addr, - filtered_size, - filter_mask, - offset_size, - element_size as u32, - ) - } else if let ChunkIndexPlan::FixedArray(grid, nslots) = &index { - let fa_address = base_address + data_buf.len() as u64; - let slots = index_slots( - grid, - &pre.shape, - &pre.chunk_dims, - &written_chunks, - Some(*nslots), - )?; - let fa_bytes = build_fixed_array_at( - &slots, - offset_size, - length_size, - pre.has_filters, - fa_address, - ); - data_buf.extend_from_slice(&fa_bytes); - serialize_v4_fixed_array( - &chunk_dims_u32, - fa_address, - offset_size, - element_size as u32, - FA_PAGE_BITS, - ) - } else { - unreachable!("every chunk index plan is handled above") + let layout_message = match &index { + ChunkIndexPlan::ExtensibleArray(grid) => { + let ea_address = base_address + data_buf.len() as u64; + let slots = index_slots(grid, &pre.shape, &pre.chunk_dims, &written_chunks, None)?; + let ea_bytes = ea_writer::build_extensible_array_at( + &slots, + offset_size, + length_size, + pre.has_filters, + ea_address, + ); + data_buf.extend_from_slice(&ea_bytes); + ea_writer::serialize_v4_extensible_array( + &chunk_dims_u32, + ea_address, + offset_size, + element_size as u32, + ) + } + ChunkIndexPlan::SingleChunk => { + let chunk_addr = written_chunks[0].address; + let filtered_size = if pre.has_filters { + Some(written_chunks[0].compressed_size) + } else { + None + }; + let filter_mask = if pre.has_filters { Some(0u32) } else { None }; + serialize_v4_single_chunk( + &chunk_dims_u32, + chunk_addr, + filtered_size, + filter_mask, + offset_size, + element_size as u32, + ) + } + ChunkIndexPlan::FixedArray(grid, nslots) => { + let fa_address = base_address + data_buf.len() as u64; + let slots = index_slots( + grid, + &pre.shape, + &pre.chunk_dims, + &written_chunks, + Some(*nslots), + )?; + let fa_bytes = build_fixed_array_at( + &slots, + offset_size, + length_size, + pre.has_filters, + fa_address, + ); + data_buf.extend_from_slice(&fa_bytes); + serialize_v4_fixed_array( + &chunk_dims_u32, + fa_address, + offset_size, + element_size as u32, + FA_PAGE_BITS, + ) + } + ChunkIndexPlan::BTreeV2 => { + let bt_address = base_address + data_buf.len() as u64; + let records: Vec<(Vec, &WrittenChunk)> = written_chunks + .iter() + .enumerate() + .map(|(i, c)| (scaled_coords(&pre.shape, &pre.chunk_dims, i), c)) + .collect(); + let (bt_bytes, node_size) = build_btree_v2_chunk_index_at( + pre.shape.len(), + &records, + offset_size, + length_size, + pre.has_filters, + bt_address, + )?; + data_buf.extend_from_slice(&bt_bytes); + serialize_v4_btree_v2( + &chunk_dims_u32, + bt_address, + offset_size, + element_size as u32, + node_size, + ) + } }; Ok(ChunkedDataResult { @@ -807,14 +840,15 @@ pub fn build_chunked_data_from_precompressed( const MAX_FIXED_ARRAY_SLOTS: u64 = 1 << 26; /// Which chunk index a dataset gets, following the library's choice in -/// `H5D__layout_set_latest_indexing`: Extensible Array for exactly one -/// unlimited dimension, Fixed Array for a finite maxshape, Single Chunk when -/// the whole maximum extent is one chunk. +/// `H5D__layout_set_latest_indexing`: version-2 B-tree for more than one +/// unlimited dimension, Extensible Array for exactly one, Fixed Array for a +/// finite maxshape, Single Chunk when the whole maximum extent is one chunk. enum ChunkIndexPlan { SingleChunk, /// The grid and the number of array elements (chunks of the max extent). FixedArray(ChunkGrid, usize), ExtensibleArray(ChunkGrid), + BTreeV2, } impl ChunkIndexPlan { @@ -860,10 +894,7 @@ impl ChunkIndexPlan { Some(max), chunk_dims, )?)), - _ => Err(bad( - "more than one unlimited dimension needs a B-tree v2 chunk index, \ - which the writer does not support", - )), + _ => Ok(Self::BTreeV2), } } } @@ -879,20 +910,9 @@ fn index_slots( chunks: &[WrittenChunk], len: Option, ) -> Result>, FormatError> { - let rank = shape.len(); - let cur: Vec = shape - .iter() - .zip(chunk_dims) - .map(|(&s, &c)| s.div_ceil(c)) - .collect(); let mut placed: Vec<(usize, &WrittenChunk)> = Vec::with_capacity(chunks.len()); - let mut scaled = vec![0u64; rank]; for (i, chunk) in chunks.iter().enumerate() { - let mut rem = i as u64; - for d in (0..rank).rev() { - scaled[d] = rem % cur[d]; - rem /= cur[d]; - } + let scaled = scaled_coords(shape, chunk_dims, i); let idx = usize::try_from(grid.linear_index(&scaled)) .map_err(|_| FormatError::Overflow("chunk index slot".into()))?; placed.push((idx, chunk)); @@ -907,6 +927,136 @@ fn index_slots( Ok(slots) } +/// Scaled coordinates (`offset / chunk_dim`) of the `i`-th chunk in the +/// row-major order `split_into_chunks` produces over the current extent. +fn scaled_coords(shape: &[u64], chunk_dims: &[u64], i: usize) -> Vec { + let rank = shape.len(); + let mut scaled = vec![0u64; rank]; + let mut rem = i as u64; + for d in (0..rank).rev() { + let n = shape[d].div_ceil(chunk_dims[d]); + scaled[d] = rem % n; + rem /= n; + } + scaled +} + +/// Node size the library gives a chunk index B-tree (`H5D_BT2_NODE_SIZE`), +/// with its split and merge percentages. +const BT2_NODE_SIZE: u32 = 2048; +const BT2_SPLIT_PERCENT: u8 = 100; +const BT2_MERGE_PERCENT: u8 = 40; +/// B-tree v2 record types for chunk indexes (`H5B2_CDSET_ID`, +/// `H5B2_CDSET_FILT_ID`). +const BT2_CHUNK_UNFILTERED: u8 = 10; +const BT2_CHUNK_FILTERED: u8 = 11; + +/// Build a version-2 B-tree chunk index (the library's index for datasets +/// with more than one unlimited dimension) at a known absolute address. +/// +/// `records` are `(scaled coordinates, chunk)` in lexicographic order of the +/// coordinates, which is the order the library's comparator +/// (`H5VM_vector_cmp_u`) keeps them in. The tree is a single leaf: the +/// library's 2048-byte node when the records fit, otherwise a leaf node +/// sized to hold them all (the root's record count is 16-bit, so at most +/// 65535 chunks). Returns the bytes and the node size the layout message +/// must record. +fn build_btree_v2_chunk_index_at( + rank: usize, + records: &[(Vec, &WrittenChunk)], + offset_size: u8, + length_size: u8, + has_filters: bool, + base_address: u64, +) -> Result<(Vec, u32), FormatError> { + let os = offset_size as usize; + let nrec = u16::try_from(records.len()).map_err(|_| { + FormatError::ChunkedReadError( + "more than 65535 chunks with more than one unlimited dimension: \ + use larger chunks" + .into(), + ) + })?; + let chunk_size_bytes = has_filters.then(|| { + let slots: Vec> = + records.iter().map(|(_, c)| Some((*c).clone())).collect(); + filtered_chunk_size_len(&slots) + }); + let record_size = os + chunk_size_bytes.map_or(0, |n| n + 4) + 8 * rank; + // Leaf: signature, version, type, records, checksum. + let leaf_len = 4 + 1 + 1 + records.len() * record_size + 4; + let node_size = u32::try_from(leaf_len) + .map_err(|_| FormatError::Overflow("B-tree v2 leaf size".into()))? + .max(BT2_NODE_SIZE); + let tree_type = if has_filters { + BT2_CHUNK_FILTERED + } else { + BT2_CHUNK_UNFILTERED + }; + + let hdr_len = 4 + 1 + 1 + 4 + 2 + 2 + 1 + 1 + os + 2 + length_size as usize + 4; + let leaf_address = base_address + hdr_len as u64; + + let mut out = Vec::with_capacity(hdr_len + node_size as usize); + out.extend_from_slice(b"BTHD"); + out.push(0); // version + out.push(tree_type); + out.extend_from_slice(&node_size.to_le_bytes()); + out.extend_from_slice(&(record_size as u16).to_le_bytes()); + out.extend_from_slice(&0u16.to_le_bytes()); // depth + out.push(BT2_SPLIT_PERCENT); + out.push(BT2_MERGE_PERCENT); + if records.is_empty() { + out.extend(core::iter::repeat_n(0xFF, os)); + } else { + push_addr(&mut out, leaf_address, offset_size); + } + out.extend_from_slice(&nrec.to_le_bytes()); + match length_size { + 4 => out.extend_from_slice(&(records.len() as u32).to_le_bytes()), + _ => out.extend_from_slice(&(records.len() as u64).to_le_bytes()), + } + let sum = jenkins_lookup3(&out); + out.extend_from_slice(&sum.to_le_bytes()); + debug_assert_eq!(out.len(), hdr_len); + if records.is_empty() { + return Ok((out, node_size)); + } + + let leaf_start = out.len(); + out.extend_from_slice(b"BTLF"); + out.push(0); // version + out.push(tree_type); + for (scaled, chunk) in records { + push_index_element(&mut out, Some(chunk), offset_size, chunk_size_bytes); + for &c in scaled { + out.extend_from_slice(&c.to_le_bytes()); + } + } + let sum = jenkins_lookup3(&out[leaf_start..]); + out.extend_from_slice(&sum.to_le_bytes()); + // The library reads whole nodes; pad the leaf out to the node size. + out.resize(leaf_start + node_size as usize, 0); + Ok((out, node_size)) +} + +/// Serialize a v4 layout message for a version-2 B-tree chunk index. +fn serialize_v4_btree_v2( + chunk_dims: &[u32], + btree_address: u64, + offset_size: u8, + element_size: u32, + node_size: u32, +) -> Vec { + let mut buf = layout_v4_chunked_prefix(chunk_dims, element_size); + buf.push(5); // chunk index type = 5 (version-2 B-tree) + buf.extend_from_slice(&node_size.to_le_bytes()); + buf.push(BT2_SPLIT_PERCENT); + buf.push(BT2_MERGE_PERCENT); + push_addr(&mut buf, btree_address, offset_size); + buf +} + /// Build chunked data with absolute addresses. /// If `maxshape` has unlimited dims, uses Extensible Array index. pub fn build_chunked_data_at( diff --git a/crates/clawhdf5/tests/chunk_index_interop.rs b/crates/clawhdf5/tests/chunk_index_interop.rs index 41d7a5a..2fb5aeb 100644 --- a/crates/clawhdf5/tests/chunk_index_interop.rs +++ b/crates/clawhdf5/tests/chunk_index_interop.rs @@ -515,16 +515,35 @@ fn we_write_maxshape_larger_than_shape() { check_we_write(&cases); } -/// More than one unlimited dimension needs a B-tree v2 chunk index; the -/// writer must not produce a file libhdf5 cannot open. +/// More than one unlimited dimension needs a version-2 B-tree chunk index, +/// as the library uses; an Extensible Array for `(None, None)` made libhdf5 +/// refuse the whole file ("already found unlimited dimension"). #[test] -fn two_unlimited_dims_are_refused() { +fn we_write_btree_v2_for_several_unlimited_dims() { + const U: u64 = u64::MAX; + let mut cases = vec![ + wcase("unl_unl", &[20, 30], &[5, 5], Some(&[U, U])), + wcase("unl_fin_unl", &[6, 7, 8], &[4, 3, 5], Some(&[U, 9, U])), + // More records than the library's 2048-byte node holds (84 here). + wcase("unl_unl_2400", &[40, 60], &[1, 1], Some(&[U, U])), + wcase("unl_unl_empty", &[0, 0], &[4, 4], Some(&[U, U])), + ]; + let mut filtered = wcase("unl_unl_deflate", &[6, 7, 8], &[4, 3, 5], Some(&[U, U, U])); + filtered.deflate = true; + cases.push(filtered); + check_we_write(&cases); +} + +/// A single-leaf B-tree has a 16-bit record count; beyond it the writer +/// refuses rather than writing a tree libhdf5 would misread. +#[test] +fn btree_v2_index_past_one_leaf_is_refused() { let mut b = FileBuilder::new(); b.create_dataset("d") - .with_i32_data(&(0..600).collect::>()) - .with_shape(&[20, 30]) - .with_chunks(&[5, 5]) + .with_i32_data(&vec![0i32; 70_000]) + .with_shape(&[70_000, 1]) + .with_chunks(&[1, 1]) .with_maxshape(&[u64::MAX, u64::MAX]); let dir = tempfile::tempdir().unwrap(); - assert!(b.write(dir.path().join("unl_unl.h5")).is_err()); + assert!(b.write(dir.path().join("too_many.h5")).is_err()); } From f5505fb03dce4944c2b9aa32a2618e4d6f485ae8 Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:15:55 -0500 Subject: [PATCH 6/7] fix(format): keep maxshape == shape datasets contiguous Any maxshape forced chunked storage, even one equal to the shape, which cannot grow. h5py and the library store such a dataset contiguously; we now do too unless chunks (or a filter) are requested. Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/file_writer.rs | 7 +++- crates/clawhdf5/tests/chunk_index_interop.rs | 44 ++++++++++++++++++++ 2 files changed, 50 insertions(+), 1 deletion(-) diff --git a/crates/clawhdf5-format/src/file_writer.rs b/crates/clawhdf5-format/src/file_writer.rs index 5350c88..644b64e 100644 --- a/crates/clawhdf5-format/src/file_writer.rs +++ b/crates/clawhdf5-format/src/file_writer.rs @@ -1124,7 +1124,12 @@ impl FileWriter { let is_chunked: Vec = all_ds .iter() .enumerate() - .map(|(i, d)| !is_vds[i] && (d.chunk_options.is_chunked() || d.maxshape.is_some())) + .map(|(i, d)| { + // Only a dataset that can grow needs chunks; a maxshape equal + // to the shape is as fixed as no maxshape at all. + let resizable = d.maxshape.as_ref().is_some_and(|m| *m != d.ds.dimensions); + !is_vds[i] && (d.chunk_options.is_chunked() || resizable) + }) .collect(); // Determine which datasets use compact storage let is_compact: Vec = all_ds diff --git a/crates/clawhdf5/tests/chunk_index_interop.rs b/crates/clawhdf5/tests/chunk_index_interop.rs index 2fb5aeb..b72b7e2 100644 --- a/crates/clawhdf5/tests/chunk_index_interop.rs +++ b/crates/clawhdf5/tests/chunk_index_interop.rs @@ -547,3 +547,47 @@ fn btree_v2_index_past_one_leaf_is_refused() { let dir = tempfile::tempdir().unwrap(); assert!(b.write(dir.path().join("too_many.h5")).is_err()); } + +/// A maxshape equal to the shape cannot grow, so it needs no chunks: the +/// dataset stays contiguous (as h5py makes it) unless chunks are requested. +#[test] +fn maxshape_equal_to_shape_stays_contiguous() { + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("ms_eq.h5"); + let data: Vec = (0..40).collect(); + let mut b = FileBuilder::new(); + b.create_dataset("plain") + .with_i32_data(&data) + .with_shape(&[40]) + .with_maxshape(&[40]); + b.create_dataset("chunked") + .with_i32_data(&data) + .with_shape(&[40]) + .with_maxshape(&[40]) + .with_chunks(&[8]); + b.write(&path).unwrap(); + + let file = File::open(&path).unwrap(); + let plain = file.dataset("plain").unwrap(); + assert_eq!(plain.read_i32().unwrap(), data); + assert_eq!(plain.max_dimensions().unwrap(), Some(vec![40])); + assert!( + plain.read_raw_ref().unwrap().is_some(), + "maxshape == shape should be contiguous" + ); + let chunked = file.dataset("chunked").unwrap(); + assert_eq!(chunked.read_i32().unwrap(), data); + assert!(chunked.read_raw_ref().unwrap().is_none()); + + skip_if_no_python!(); + let out = run_python(&format!( + "import h5py, numpy as np\n\ + f = h5py.File(r'{}', 'r')\n\ + for n in ('plain', 'chunked'):\n\ + \x20 d = f[n]\n\ + \x20 assert np.array_equal(d[()], np.arange(40, dtype='i4')), n\n\ + \x20 print(n, d.chunks, d.maxshape)\n", + path.display() + )); + assert_eq!(out, "plain None (40,)\nchunked (8,) (40,)"); +} From e7f2d8575da92c275dbdd30a0d0c12905e4c34a6 Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:18:13 -0500 Subject: [PATCH 7/7] fix(format): import format! for the no_std chunk index planner The maxshape checks added to chunked_write use format!, which a no_std build has to import from alloc (scripts/check-nostd.sh). Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/chunked_write.rs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/crates/clawhdf5-format/src/chunked_write.rs b/crates/clawhdf5-format/src/chunked_write.rs index 6ff0c07..682a1c6 100644 --- a/crates/clawhdf5-format/src/chunked_write.rs +++ b/crates/clawhdf5-format/src/chunked_write.rs @@ -4,7 +4,7 @@ extern crate alloc; #[cfg(not(feature = "std"))] -use alloc::{vec, vec::Vec}; +use alloc::{format, vec, vec::Vec}; use crate::checksum::jenkins_lookup3; use crate::chunk_cache::{CACHE_LINE_SIZE, align_to_cache_line};