Merge branch 'perf/p2b-chunked-full-reads' into feat/p2b-scale

This commit is contained in:
osobh
2026-09-26 11:57:11 -05:00
12 changed files with 2095 additions and 387 deletions
+46
View File
@@ -2,6 +2,52 @@
## Unreleased
### Chunked full reads (2026-09-26)
- **Chunks are decoded straight into the output, into reused buffers.** A
full read of a chunked dataset faulted in about three times its size in
fresh pages: every chunk was decoded into a new buffer per filter stage
(the cached reader behind `read_*` decoded 128 chunks at a time before
placing any; the uncached one behind `MmapFile`, `LazyFile` and
`verify_provenance` decoded the whole dataset first), then assembled into
a byte buffer, which the typed readers copied once more. Now each chunk
is decoded into buffers the thread keeps between chunks and reads
(`clawhdf5_format::filters::DecodeScratch`,
`decompress_chunk_exact_with`: deflate inflates into a kept buffer with a
reset inflater, shuffle into the other one, Fletcher32 is checked in
place; other filters go through the registry as before) and copied
directly to its place in the output. Chunks still go into the file's
chunk cache when the whole dataset fits. Selection reads decode the
chunks they touch the same way.
- **Typed full reads of chunked data skip the byte buffer.** `read_f32`,
`read_f64`, `read_i32`, `read_i64` and `read_u64` (on `File`, `MmapFile`
and `LazyFile`) of a chunked dataset stored as that type in native byte
order decode every chunk into the returned `Vec` (huge-page backed when
large, like the byte readers' output); other types and byte orders
convert as before. New public
`clawhdf5_format::data_read::read_chunked_native`.
- **Reading threads no longer wait for a busy rayon pool.** A full read
handed its chunks to rayon and the calling thread slept until the pool had
decoded them, so with a small pool (2-4 threads) readers outside it queued
behind its workers. The calling thread now decodes too, and pool workers
join in only when free; a helper the pool starts after the read has
finished returns at once. A single read still spreads over the default
pool. This replaces the one-thread-pool special case below for full
reads. Chunks are placed from several threads only when the chunk index
puts them on the chunk grid at distinct places; a corrupt index is read
one chunk at a time, and the error reported is still the first failing
chunk's. New test `crates/clawhdf5/tests/busy_decode_pool.rs`.
- **Fixed:** in a filtered dataset, a chunk stored with every filter skipped
(filter mask) and shorter than a chunk read with zeros in place of its
missing part through `File`'s `read_*`; it is now an error naming the
chunk, as `MmapFile`/`LazyFile` already made it.
- New h5py comparison `crates/clawhdf5/tests/chunked_read_paths_interop.rs`:
every chunked read path (cached and uncached full reads, `MmapFile`,
`LazyFile`, small, strided and point selections, with and without the
`parallel` feature) for 1-8-byte integers and 2-8-byte floats in both
byte orders, through deflate, shuffle, Fletcher32, LZF, SZIP and Blosc,
with partial edge chunks, sparse datasets and fill values, and datasets
larger than the chunk cache.
### Concurrent reads (2026-09-26)
- **Full reads of chunked datasets scale with threads again when rayon's
pool has one thread.** Each full read handed its chunks to rayon to
File diff suppressed because it is too large Load Diff
+114
View File
@@ -756,6 +756,120 @@ pub fn read_selection_native<T: NativeElement>(
crate::gather::gather::<T>(raw, dims, elem_size, selection).map(Some)
}
/// The bytes of a slice of [`NativeElement`]s.
#[cfg(feature = "std")]
fn bytes_of_mut<T: NativeElement>(values: &mut [T]) -> &mut [u8] {
// SAFETY: `T: NativeElement` has no padding and every bit pattern is a
// valid value, so its storage may be viewed, and written, as bytes; the
// byte slice covers exactly the values' storage and borrows it
// exclusively for its lifetime.
unsafe {
core::slice::from_raw_parts_mut(
values.as_mut_ptr().cast::<u8>(),
core::mem::size_of_val(values),
)
}
}
/// `count` zeroed values of `T`, from zeroed pages where the allocator can
/// (see [`crate::chunked_read::alloc_output`]) and backed by huge pages when
/// large. A size taken from the file surfaces as an error, not an abort.
#[cfg(feature = "std")]
fn alloc_zeroed_values<T: NativeElement>(count: usize) -> Result<Vec<T>, FormatError> {
if count == 0 || core::mem::size_of::<T>() == 0 {
return Ok(Vec::new());
}
let failed = || {
FormatError::Overflow(format!(
"cannot allocate {count} values of {} bytes for dataset output",
core::mem::size_of::<T>()
))
};
let layout = core::alloc::Layout::array::<T>(count).map_err(|_| failed())?;
// SAFETY: `layout` has non-zero size (count > 0, T not zero-sized).
let ptr = unsafe { std::alloc::alloc_zeroed(layout) };
if ptr.is_null() {
return Err(failed());
}
crate::bulk_alloc::advise_huge_pages(ptr, layout.size());
// SAFETY: allocated by the global allocator with the layout of
// `[T; count]`, which is what `Vec<T>` with capacity `count` frees; all
// bytes are zero, a valid `T` (`NativeElement`: any bit pattern is).
Ok(unsafe { Vec::from_raw_parts(ptr.cast::<T>(), count, count) })
}
/// Read a whole chunked dataset that stores `T` natively
/// ([`NativeElement::is_native`]) straight into a `Vec<T>`: each chunk is
/// decoded and copied to its place in the typed output, with no byte buffer
/// to convert from afterwards. Unallocated chunks read as the dataset's fill
/// value, as [`crate::fill_value::read_full_with_fill`] makes them.
///
/// `Ok(None)` when this does not apply — the datatype is not `T`'s native
/// representation (another type, another byte order: the caller converts
/// through the byte readers and the `read_as_*` functions), the layout is
/// not chunked, no storage is allocated, or the data lives in external
/// files. `cache` is the file's chunk cache, used as
/// [`crate::chunked_read::read_chunked_data_cached`] uses it.
#[cfg(feature = "std")]
#[allow(clippy::too_many_arguments)]
pub fn read_chunked_native<T: NativeElement>(
messages: &[crate::object_header::HeaderMessage],
file_data: &[u8],
layout: &DataLayout,
dataspace: &Dataspace,
datatype: &Datatype,
pipeline: Option<&FilterPipeline>,
offset_size: u8,
length_size: u8,
cache: Option<&ChunkCache>,
) -> Result<Option<Vec<T>>, FormatError> {
use crate::fill_value;
use crate::message_type::MessageType;
if !T::is_native(datatype)
|| !matches!(layout, DataLayout::Chunked { .. })
|| !fill_value::has_storage(layout)
|| messages
.iter()
.any(|m| m.msg_type == MessageType::ExternalDataFiles)
{
return Ok(None);
}
let size = core::mem::size_of::<T>();
let mut values = crate::chunked_read::read_chunked_full(
file_data,
layout,
dataspace,
datatype,
pipeline,
offset_size,
length_size,
cache,
|total_bytes| {
if !total_bytes.is_multiple_of(size) {
return Err(FormatError::DataSizeMismatch {
expected: total_bytes.next_multiple_of(size),
actual: total_bytes,
});
}
alloc_zeroed_values::<T>(total_bytes / size)
},
|values| bytes_of_mut(values),
)?;
let fill = fill_value::dataset_fill_value_in(file_data, messages, offset_size, length_size)?;
fill_value::apply_to_unallocated_chunks(
bytes_of_mut(&mut values),
file_data,
layout,
dataspace,
size,
fill.as_deref(),
offset_size,
length_size,
)?;
Ok(Some(values))
}
/// Convert raw bytes to `f64` values.
pub fn read_as_f64(raw: &[u8], datatype: &Datatype) -> Result<Vec<f64>, FormatError> {
// Array datatypes read as a flat sequence of their base elements, and
+345 -21
View File
@@ -148,6 +148,183 @@ pub fn decompress_chunk_exact(
Ok(data)
}
/// Buffers a chunk decoder keeps between chunks, so decoding a dataset's
/// chunks one after another reuses the same memory instead of allocating
/// (and faulting in) fresh buffers for every chunk and every filter stage.
///
/// Use one per thread with [`decompress_chunk_exact_with`]. Buffers larger
/// than [`DecodeScratch::RETAIN_BYTES`] are released by
/// [`DecodeScratch::trim`], so a scratch kept for a long time (a
/// thread-local, say) does not hold on to a huge chunk's memory.
#[derive(Default)]
pub struct DecodeScratch {
a: Vec<u8>,
b: Vec<u8>,
#[cfg(feature = "deflate")]
inflater: Option<flate2::Decompress>,
}
impl core::fmt::Debug for DecodeScratch {
fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
f.debug_struct("DecodeScratch")
.field("a_capacity", &self.a.capacity())
.field("b_capacity", &self.b.capacity())
.finish_non_exhaustive()
}
}
/// Which buffer holds the data between two filter stages.
#[derive(Clone, Copy)]
enum Stage {
/// `compressed[..len]`: still (a prefix of) the stored bytes.
Stored(usize),
A,
B,
}
impl DecodeScratch {
/// Largest buffer [`trim`](Self::trim) keeps (1 MiB): enough for common
/// chunk sizes (a 256 x 256 `f32` chunk is 256 KiB) while bounding what
/// every long-lived thread (rayon's workers never exit) holds on to, at
/// two buffers each.
pub const RETAIN_BYTES: usize = 1 << 20;
/// An empty scratch; buffers are allocated on first use.
pub fn new() -> Self {
Self::default()
}
/// Release any buffer larger than [`Self::RETAIN_BYTES`].
pub fn trim(&mut self) {
for buf in [&mut self.a, &mut self.b] {
if buf.capacity() > Self::RETAIN_BYTES {
*buf = Vec::new();
}
}
}
}
/// [`decompress_chunk_exact`] into reusable buffers: the decoded chunk is
/// returned as a slice of `scratch` (or of `compressed`, when every filter
/// that was applied only appended a checksum), valid until `scratch` is used
/// again.
///
/// Deflate, shuffle and Fletcher32 — h5py's and libhdf5's usual pipeline —
/// decode without allocating once `scratch` has grown to the chunk size: the
/// inflater writes into a kept buffer (and its state is reset, not
/// rebuilt), shuffle interleaves into the other buffer, and Fletcher32 checks
/// the checksum and drops it in place. Every other filter goes through the
/// filter registry as [`decompress_chunk_masked`] does, and its output
/// replaces a scratch buffer. The result is byte for byte what
/// [`decompress_chunk_exact`] returns, with the same errors.
pub fn decompress_chunk_exact_with<'s>(
compressed: &'s [u8],
pipeline: &FilterPipeline,
chunk_size: usize,
element_size: u32,
filter_mask: u32,
coords: &[u64],
scratch: &'s mut DecodeScratch,
) -> Result<&'s [u8], FormatError> {
// Same per-stage bounds as `decompress_chunk_masked`.
let mut bounds = [0usize; 32];
let mut bounds_vec = Vec::new();
let bounds: &mut [usize] = if pipeline.filters.len() <= bounds.len() {
&mut bounds[..pipeline.filters.len()]
} else {
bounds_vec.resize(pipeline.filters.len(), 0);
&mut bounds_vec
};
let mut size = chunk_size;
for (i, filter) in pipeline.filters.iter().enumerate() {
bounds[i] = size;
if !filter_skipped(filter_mask, i) {
size = filter_output_bound(filter.filter_id, size);
}
}
let mut stage = Stage::Stored(compressed.len());
for (i, filter) in pipeline.filters.iter().enumerate().rev() {
if filter_skipped(filter_mask, i) {
continue;
}
let ctx = FilterContext {
filter,
element_size: element_size as usize,
max_output: bounds[i],
};
// The stage's input, and the buffer its output goes to (the one
// not holding the input).
let (input, out): (&[u8], &mut Vec<u8>) = match stage {
Stage::Stored(len) => (&compressed[..len], &mut scratch.a),
Stage::A => (&scratch.a, &mut scratch.b),
Stage::B => (&scratch.b, &mut scratch.a),
};
let next = match stage {
Stage::Stored(_) | Stage::B => Stage::A,
Stage::A => Stage::B,
};
match filter.filter_id {
// Built in and never overridable (`register_filter` refuses
// built-in IDs), so the registry would pick exactly these.
FILTER_FLETCHER32 => {
// Check and drop the checksum where the data is.
let payload = fletcher32_payload(input)?;
stage = match stage {
Stage::Stored(_) => Stage::Stored(payload),
Stage::A => {
scratch.a.truncate(payload);
Stage::A
}
Stage::B => {
scratch.b.truncate(payload);
Stage::B
}
};
continue;
}
FILTER_SHUFFLE => shuffle_decompress_into(input, ctx.element_size, out),
#[cfg(all(
feature = "deflate",
not(all(target_os = "macos", feature = "system-zlib-decompress"))
))]
FILTER_DEFLATE => {
let limit = if ctx.max_output != 0 {
ctx.max_output
} else {
MAX_DECOMPRESS_SIZE
};
let size_hint = if ctx.max_output != 0 {
ctx.max_output
} else {
input.len().saturating_mul(4).min(1 << 20)
};
let inflater = scratch
.inflater
.get_or_insert_with(|| flate2::Decompress::new(true));
inflater.reset(true);
inflate_bounded_into(inflater, input, size_hint, limit, out)
.map_err(FormatError::DecompressionError)?;
}
_ => *out = filter_registry::decode(input, &ctx)?,
}
stage = next;
}
let data: &[u8] = match stage {
Stage::Stored(len) => &compressed[..len],
Stage::A => &scratch.a,
Stage::B => &scratch.b,
};
if chunk_size != 0 && data.len() != chunk_size {
return Err(FormatError::ChunkedReadError(format!(
"chunk at {coords:?} decoded to {} bytes, expected {chunk_size}",
data.len()
)));
}
Ok(data)
}
/// Apply a filter pipeline to compress a chunk.
/// Filters are applied in FORWARD order for compression.
pub fn compress_chunk(
@@ -892,33 +1069,55 @@ pub(crate) fn inflate_bounded(
size_hint: usize,
limit: usize,
) -> Result<Vec<u8>, String> {
use flate2::{Decompress, FlushDecompress, Status};
let mut out = Vec::new();
inflate_bounded_into(
&mut flate2::Decompress::new(true),
data,
size_hint,
limit,
&mut out,
)?;
Ok(out)
}
/// [`inflate_bounded`] with a fresh or reset `inflater`, into `out`: its
/// contents are replaced and its allocation reused.
#[cfg(feature = "deflate")]
fn inflate_bounded_into(
inflater: &mut flate2::Decompress,
data: &[u8],
size_hint: usize,
limit: usize,
out: &mut Vec<u8>,
) -> Result<(), String> {
use flate2::{FlushDecompress, Status};
// One byte of headroom past the limit distinguishes an over-size stream
// from one that legitimately ends exactly at the limit.
let max_capacity = limit.saturating_add(1);
let mut out = Vec::new();
out.try_reserve_exact(size_hint.clamp(1, max_capacity))
// A kept buffer may already be larger than `max_capacity`; the decoder
// can then write past the limit, which the check below still refuses.
out.clear();
let want = size_hint.clamp(1, max_capacity);
out.try_reserve_exact(want)
.map_err(|e| format!("deflate: cannot allocate output: {e}"))?;
let mut inflater = Decompress::new(true);
loop {
let (in_before, out_before) = (inflater.total_in(), inflater.total_out());
let status = inflater
.decompress_vec(
&data[in_before as usize..],
&mut out,
FlushDecompress::Finish,
)
.decompress_vec(&data[in_before as usize..], out, FlushDecompress::Finish)
.map_err(|e| format!("deflate: {e}"))?;
if out.len() > limit {
return Err("deflate: output exceeds size limit".into());
}
match status {
Status::StreamEnd => return Ok(out),
Status::StreamEnd => return Ok(()),
Status::Ok | Status::BufError if out.len() == out.capacity() => {
// Out of room: double, up to the limit.
let grow = out.capacity().min(max_capacity - out.capacity()).max(1);
let grow = out
.capacity()
.min(max_capacity.saturating_sub(out.capacity()))
.max(1);
out.try_reserve_exact(grow)
.map_err(|e| format!("deflate: cannot allocate output: {e}"))?;
}
@@ -1217,16 +1416,30 @@ fn zstd_compress(data: &[u8], level: u32) -> Result<Vec<u8>, FormatError> {
/// On disk: all byte-0s of each element together, then all byte-1s, etc.
/// Output: elements in natural order.
fn shuffle_decompress(data: &[u8], element_size: usize) -> Result<Vec<u8>, FormatError> {
let mut result = Vec::new();
shuffle_decompress_into(data, element_size, &mut result);
Ok(result)
}
/// [`shuffle_decompress`] into `result`, replacing its contents and reusing
/// its allocation.
fn shuffle_decompress_into(data: &[u8], element_size: usize, result: &mut Vec<u8>) {
if element_size <= 1 {
return Ok(data.to_vec());
result.clear();
result.extend_from_slice(data);
return;
}
// Like libhdf5, only whole elements are shuffled; trailing bytes (e.g. a
// Fletcher32 checksum appended before the shuffle) are stored as-is.
let whole = data.len() - data.len() % element_size;
let (data, tail) = data.split_at(whole);
let num_elements = data.len() / element_size;
let mut result = vec![0u8; whole];
result.reserve_exact(tail.len());
// Every byte of `result[..whole]` is overwritten below, so a reused
// buffer keeps its old bytes instead of being zeroed first; only growth
// is zero-filled.
result.truncate(whole);
result.reserve_exact(whole + tail.len() - result.len());
result.resize(whole, 0);
// The shuffled stream is `element_size` byte planes of `num_elements`
// bytes each; un-shuffling interleaves them. This is on the read path of
@@ -1245,10 +1458,10 @@ fn shuffle_decompress(data: &[u8], element_size: usize) -> Result<Vec<u8>, Forma
}
}
match element_size {
2 => interleave::<2>(data, num_elements, &mut result),
4 => interleave::<4>(data, num_elements, &mut result),
8 => interleave::<8>(data, num_elements, &mut result),
16 => interleave::<16>(data, num_elements, &mut result),
2 => interleave::<2>(data, num_elements, result),
4 => interleave::<4>(data, num_elements, result),
8 => interleave::<8>(data, num_elements, result),
16 => interleave::<16>(data, num_elements, result),
_ => {
for (i, element) in result.chunks_exact_mut(element_size).enumerate() {
for (j, byte) in element.iter_mut().enumerate() {
@@ -1258,8 +1471,6 @@ fn shuffle_decompress(data: &[u8], element_size: usize) -> Result<Vec<u8>, Forma
}
}
result.extend_from_slice(tail);
Ok(result)
}
/// Shuffle (compress direction): group bytes by position within each element.
@@ -1411,6 +1622,12 @@ fn fletcher32_compute(data: &[u8]) -> u32 {
/// Verify Fletcher32 checksum and strip it from the data.
/// The last 4 bytes are the stored checksum.
fn fletcher32_verify(data: &[u8]) -> Result<Vec<u8>, FormatError> {
fletcher32_payload(data).map(|len| data[..len].to_vec())
}
/// Verify the Fletcher32 checksum that ends `data`; the length of the data
/// before it.
fn fletcher32_payload(data: &[u8]) -> Result<usize, FormatError> {
if data.len() < 4 {
return Err(FormatError::FilterError(
"fletcher32: data too short for checksum".into(),
@@ -1430,7 +1647,7 @@ fn fletcher32_verify(data: &[u8]) -> Result<Vec<u8>, FormatError> {
computed,
});
}
Ok(payload.to_vec())
Ok(payload.len())
}
/// Append Fletcher32 checksum to data.
@@ -1545,6 +1762,113 @@ fn pcodec_decompress(
#[cfg(test)]
mod tests {
/// `decompress_chunk_exact_with` returns exactly what
/// `decompress_chunk_exact` returns — data or error — for every pipeline
/// shape, filter mask and chunk size, with one scratch reused across all
/// of them in an order that grows, shrinks and swaps its buffers.
#[test]
fn decode_with_scratch_matches_the_allocating_decoder() {
let f = |filter_id: u16, client_data: Vec<u32>| FilterDescription {
filter_id,
name: None,
flags: 0,
client_data,
};
let mut pipelines = vec![
vec![f(FILTER_SHUFFLE, vec![4])],
vec![f(FILTER_FLETCHER32, vec![])],
// NetCDF-4's order: the checksum is taken before shuffle.
vec![f(FILTER_FLETCHER32, vec![]), f(FILTER_SHUFFLE, vec![4])],
];
#[cfg(feature = "deflate")]
pipelines.extend([
vec![f(FILTER_DEFLATE, vec![4])],
vec![f(FILTER_SHUFFLE, vec![4]), f(FILTER_DEFLATE, vec![4])],
// h5py's order with `fletcher32=True`: checksum last.
vec![
f(FILTER_SHUFFLE, vec![4]),
f(FILTER_DEFLATE, vec![4]),
f(FILTER_FLETCHER32, vec![]),
],
vec![
f(FILTER_FLETCHER32, vec![]),
f(FILTER_SHUFFLE, vec![4]),
f(FILTER_DEFLATE, vec![1]),
],
]);
#[cfg(feature = "lzf")]
pipelines.push(vec![
f(FILTER_SHUFFLE, vec![4]),
f(crate::filter_pipeline::FILTER_LZF, vec![]),
]);
let mut scratch = DecodeScratch::new();
for elements in [1usize, 7, 4096, 3, 65536, 100] {
let data: Vec<u8> = (0..elements as u32)
.flat_map(|i| (i.wrapping_mul(2654435761) >> (i % 13)).to_le_bytes())
.collect();
for filters in &pipelines {
let pipeline = FilterPipeline {
version: 2,
filters: filters.clone(),
};
let n = filters.len() as u32;
for mask in 0..(1u32 << n) {
// Encode only the filters the mask says were applied.
let mut stored = data.clone();
for (i, filter) in filters.iter().enumerate() {
if mask & (1 << i) == 0 {
let ctx = FilterContext {
filter,
element_size: 4,
max_output: 0,
};
stored = filter_registry::encode(&stored, &ctx).unwrap();
}
}
let mut cases = vec![(stored.clone(), data.len())];
// Corrupt: last byte flipped, truncated, wrong size.
let mut flipped = stored.clone();
*flipped.last_mut().unwrap() ^= 0x5a;
cases.push((flipped, data.len()));
cases.push((stored[..stored.len() / 2].to_vec(), data.len()));
cases.push((stored.clone(), data.len() + 4));
cases.push((stored.clone(), 0));
for (bytes, size) in cases {
let want = decompress_chunk_exact(&bytes, &pipeline, size, 4, mask, &[3]);
let got = decompress_chunk_exact_with(
&bytes,
&pipeline,
size,
4,
mask,
&[3],
&mut scratch,
)
.map(<[u8]>::to_vec);
match (&want, &got) {
(Ok(w), Ok(g)) => assert_eq!(w, g, "{filters:?} mask {mask}"),
(Err(w), Err(g)) => {
assert_eq!(w.to_string(), g.to_string(), "{filters:?}")
}
_ => panic!("{filters:?} mask {mask} size {size}: {want:?} vs {got:?}"),
}
}
}
}
}
// Long-lived scratch gives back a huge chunk's buffers.
let big = vec![0u8; DecodeScratch::RETAIN_BYTES + 8];
let shuffle = FilterPipeline {
version: 2,
filters: vec![f(FILTER_SHUFFLE, vec![4])],
};
decompress_chunk_exact_with(&big, &shuffle, big.len(), 4, 0, &[0], &mut scratch).unwrap();
scratch.trim();
assert!(scratch.a.capacity() <= DecodeScratch::RETAIN_BYTES);
assert!(scratch.b.capacity() <= DecodeScratch::RETAIN_BYTES);
}
/// A chunk whose pipeline decodes to fewer bytes than the chunk holds is
/// an error naming the chunk, never a short buffer the reader pads.
#[test]
+180
View File
@@ -41,6 +41,132 @@ pub fn pool_can_parallelise() -> bool {
rayon::current_num_threads() > 1
}
/// How many rayon workers [`run_with_helpers`] should ask to help with
/// `items` work items, given that the calling thread works too: the pool's
/// other threads (all of them when the caller is not one), at most one per
/// item beyond the caller's first.
pub(crate) fn helper_count(items: usize) -> usize {
let pool = rayon::current_num_threads();
// A one-thread pool means "decode on the calling thread" (the setting
// benchmarks use to compare with h5py, where each call decodes on its
// caller): no helper, so one read never uses two cores.
if pool <= 1 {
return 0;
}
let others = if rayon::current_thread_index().is_some() {
pool.saturating_sub(1)
} else {
pool
};
others.min(items.saturating_sub(1))
}
/// Run `body` on the calling thread and on up to `helpers` rayon workers at
/// once, returning when the caller's call has finished and every worker that
/// started one has too. `body` shares its work out itself (typically by
/// claiming items from an atomic counter until none are left).
///
/// The caller never waits for a worker to *become* free: helpers are queued
/// on the pool, and one that only gets to run after the caller has finished
/// returns without calling `body`. So a busy or small pool can only fail to
/// speed a read up, never hold it back — with `par_iter`, the calling thread
/// (not a pool worker) handed all the work to the pool and slept, and N
/// threads reading through a 2-worker pool decoded on 2 cores.
///
/// A panic in `body`, on any thread, is resumed on the caller once every
/// helper that started has stopped.
pub(crate) fn run_with_helpers(helpers: usize, body: &(dyn Fn() + Sync)) {
use std::panic::{AssertUnwindSafe, catch_unwind, resume_unwind};
use std::sync::{Arc, Condvar, Mutex, PoisonError};
if helpers == 0 {
body();
return;
}
type Body = dyn Fn() + Sync + 'static;
struct Shared {
/// `body`, its lifetime erased. Only dereferenced by a helper that
/// registered in `state` while it was open (see below).
body: *const Body,
/// (closed, helpers inside `body`).
state: Mutex<(bool, usize)>,
idle: Condvar,
panic: Mutex<Option<Box<dyn core::any::Any + Send>>>,
}
// SAFETY: `body` points to a `Sync` closure, so calling it from other
// threads is allowed; the pointer is only used under the protocol below,
// which keeps it from outliving the closure.
unsafe impl Send for Shared {}
unsafe impl Sync for Shared {}
fn help(shared: &Shared) {
{
let mut state = shared.state.lock().unwrap_or_else(PoisonError::into_inner);
if state.0 {
return;
}
state.1 += 1;
}
// SAFETY: registered while open, so the caller of `run_with_helpers`
// is still inside it (it closes, then waits until no helper is
// registered, before returning), and `body` is alive.
let body = unsafe { &*shared.body };
if let Err(payload) = catch_unwind(AssertUnwindSafe(body)) {
shared
.panic
.lock()
.unwrap_or_else(PoisonError::into_inner)
.get_or_insert(payload);
}
let mut state = shared.state.lock().unwrap_or_else(PoisonError::into_inner);
state.1 -= 1;
if state.1 == 0 {
shared.idle.notify_all();
}
}
let body_ptr: *const (dyn Fn() + Sync + '_) = body;
// SAFETY: only the lifetime changes (same fat-pointer layout). The
// pointer is dereferenced only while this function is running: see
// `help` and the wait below.
let body_ptr: *const Body = unsafe { core::mem::transmute(body_ptr) };
let shared = Arc::new(Shared {
body: body_ptr,
state: Mutex::new((false, 0)),
idle: Condvar::new(),
panic: Mutex::new(None),
});
for _ in 0..helpers {
let shared = Arc::clone(&shared);
rayon::spawn(move || help(&shared));
}
let caller = catch_unwind(AssertUnwindSafe(body));
{
// Close, then wait for the helpers inside `body`; later ones return
// at once. This must happen even if `body` panicked on this thread.
let mut state = shared.state.lock().unwrap_or_else(PoisonError::into_inner);
state.0 = true;
while state.1 > 0 {
state = shared
.idle
.wait(state)
.unwrap_or_else(PoisonError::into_inner);
}
}
if let Err(payload) = caller {
resume_unwind(payload);
}
let helper_panic = shared
.panic
.lock()
.unwrap_or_else(PoisonError::into_inner)
.take();
if let Some(payload) = helper_panic {
resume_unwind(payload);
}
}
/// Decompress chunks in parallel using lane-partitioned assignment.
///
/// Instead of naive `par_iter`, chunks are deterministically assigned to lanes
@@ -258,6 +384,60 @@ mod tests {
(file, infos)
}
/// Every item is processed exactly once, whatever mix of caller and
/// helpers ends up doing it.
#[test]
fn run_with_helpers_shares_all_work() {
use core::sync::atomic::{AtomicUsize, Ordering};
for helpers in [0, 1, 3, 16] {
let n = 1000;
let next = AtomicUsize::new(0);
let done: Vec<AtomicUsize> = (0..n).map(|_| AtomicUsize::new(0)).collect();
run_with_helpers(helpers, &|| {
loop {
let i = next.fetch_add(1, Ordering::Relaxed);
if i >= n {
break;
}
done[i].fetch_add(1, Ordering::Relaxed);
}
});
assert!(done.iter().all(|d| d.load(Ordering::Relaxed) == 1));
}
}
/// A panic in the shared body reaches the caller whichever thread it
/// happened on, and only after the helpers inside the body have left it
/// (they borrow the caller's stack).
#[test]
fn run_with_helpers_propagates_panics() {
use core::sync::atomic::{AtomicUsize, Ordering};
use std::panic::{AssertUnwindSafe, catch_unwind};
let caller = std::thread::current().id();
for panic_on_caller in [true, false] {
let inside = AtomicUsize::new(0);
let calls = AtomicUsize::new(0);
let result = catch_unwind(AssertUnwindSafe(|| {
run_with_helpers(4, &|| {
inside.fetch_add(1, Ordering::SeqCst);
calls.fetch_add(1, Ordering::SeqCst);
let on_caller = std::thread::current().id() == caller;
std::thread::sleep(std::time::Duration::from_millis(20));
inside.fetch_sub(1, Ordering::SeqCst);
if on_caller == panic_on_caller {
panic!("boom");
}
});
}));
// A helper may never have run (the pool was slow to start it),
// in which case nothing panicked when `panic_on_caller` is false.
if panic_on_caller || calls.load(Ordering::SeqCst) > 1 {
assert!(result.is_err());
}
assert_eq!(inside.load(Ordering::SeqCst), 0);
}
}
/// Every parallel decoder refuses a chunk that decodes short, naming it.
#[test]
fn short_decoded_chunk_is_an_error() {
+8 -5
View File
@@ -24,7 +24,7 @@ use crate::data_read::extract_selection_from_buffer;
use crate::dataspace::Dataspace;
use crate::error::FormatError;
use crate::filter_pipeline::FilterPipeline;
use crate::filters::{all_filters_skipped, decompress_chunk_exact};
use crate::filters::{all_filters_skipped, decompress_chunk_exact_with};
use crate::selection::Selection;
/// The smallest axis-aligned box containing every selected element, as
@@ -305,6 +305,8 @@ pub fn read_selection(
let rank = dims.len();
let chunk_shape: Vec<u64> = chunk_dims.iter().map(|&d| d as u64).collect();
let chunk_bytes = crate::chunked_read::checked_chunk_byte_len(&chunk_dims, elem_size)?;
// Chunks are decoded into this thread's reusable buffers.
crate::chunked_read::with_scratch(|scratch| -> Result<(), FormatError> {
for chunk in &chunks {
if chunk.offsets.len() < rank || chunk.address == u64::MAX {
continue;
@@ -328,18 +330,17 @@ pub fn read_selection(
})?;
// Mirrors the full-read path: filter-mask bit i set means
// filter i was not applied to this chunk.
let decoded;
let data: &[u8] = match pipeline {
Some(pl) if !all_filters_skipped(pl, chunk.filter_mask) => {
decoded = decompress_chunk_exact(
decompress_chunk_exact_with(
raw,
pl,
chunk_bytes,
elem_size as u32,
chunk.filter_mask,
&chunk.offsets[..rank],
)?;
&decoded
scratch,
)?
}
_ => raw,
};
@@ -353,6 +354,8 @@ pub fn read_selection(
elem_size,
);
}
Ok(())
})?;
}
_ => return Ok(None),
}
+42
View File
@@ -349,6 +349,9 @@ impl<'f, R: HDF5Read> LazyDataset<'f, R> {
/// Read all data as `f64` values.
pub fn read_f64(&self) -> Result<Vec<f64>, Error> {
if let Some(values) = self.read_chunked_native::<f64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_f64(&raw, &dt)?)
@@ -396,6 +399,9 @@ impl<'f, R: HDF5Read> LazyDataset<'f, R> {
/// Read all data as `f32` values.
pub fn read_f32(&self) -> Result<Vec<f32>, Error> {
if let Some(values) = self.read_chunked_native::<f32>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_f32(&raw, &dt)?)
@@ -403,6 +409,9 @@ impl<'f, R: HDF5Read> LazyDataset<'f, R> {
/// Read all data as `i32` values.
pub fn read_i32(&self) -> Result<Vec<i32>, Error> {
if let Some(values) = self.read_chunked_native::<i32>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_i32(&raw, &dt)?)
@@ -410,6 +419,9 @@ impl<'f, R: HDF5Read> LazyDataset<'f, R> {
/// Read all data as `i64` values.
pub fn read_i64(&self) -> Result<Vec<i64>, Error> {
if let Some(values) = self.read_chunked_native::<i64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_i64(&raw, &dt)?)
@@ -417,6 +429,9 @@ impl<'f, R: HDF5Read> LazyDataset<'f, R> {
/// Read all data as `u64` values.
pub fn read_u64(&self) -> Result<Vec<u64>, Error> {
if let Some(values) = self.read_chunked_native::<u64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_u64(&raw, &dt)?)
@@ -550,6 +565,33 @@ impl<'f, R: HDF5Read> LazyDataset<'f, R> {
.transpose()
}
/// A chunked dataset that stores `T` natively, decoded straight into a
/// `Vec<T>` (no byte buffer to convert); `None` for any other dataset
/// (see [`data_read::read_chunked_native`]).
fn read_chunked_native<T: data_read::NativeElement>(&self) -> Result<Option<Vec<T>>, Error> {
let dl = self.data_layout()?;
if !matches!(dl, DataLayout::Chunked { .. }) {
return Ok(None);
}
let dt = self.datatype()?;
if !T::is_native(&dt) {
return Ok(None);
}
let ds = self.dataspace()?;
let pipeline = self.filter_pipeline()?;
Ok(data_read::read_chunked_native::<T>(
&self.header.messages,
self.file.hdf5_bytes(),
&dl,
&ds,
&dt,
pipeline.as_ref(),
self.file.offset_size(),
self.file.length_size(),
None,
)?)
}
fn read_raw(&self) -> Result<Vec<u8>, Error> {
let dt = self.datatype()?;
let ds = self.dataspace()?;
+42
View File
@@ -276,6 +276,9 @@ impl<'f> MmapDataset<'f> {
/// Read all data as `f64` values.
pub fn read_f64(&self) -> Result<Vec<f64>, Error> {
if let Some(values) = self.read_chunked_native::<f64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_f64(&raw, &dt)?)
@@ -310,6 +313,9 @@ impl<'f> MmapDataset<'f> {
/// Read all data as `f32` values.
pub fn read_f32(&self) -> Result<Vec<f32>, Error> {
if let Some(values) = self.read_chunked_native::<f32>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_f32(&raw, &dt)?)
@@ -317,6 +323,9 @@ impl<'f> MmapDataset<'f> {
/// Read all data as `i32` values.
pub fn read_i32(&self) -> Result<Vec<i32>, Error> {
if let Some(values) = self.read_chunked_native::<i32>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_i32(&raw, &dt)?)
@@ -324,6 +333,9 @@ impl<'f> MmapDataset<'f> {
/// Read all data as `i64` values.
pub fn read_i64(&self) -> Result<Vec<i64>, Error> {
if let Some(values) = self.read_chunked_native::<i64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_i64(&raw, &dt)?)
@@ -331,6 +343,9 @@ impl<'f> MmapDataset<'f> {
/// Read all data as `u64` values.
pub fn read_u64(&self) -> Result<Vec<u64>, Error> {
if let Some(values) = self.read_chunked_native::<u64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
let dt = self.datatype()?;
Ok(data_read::read_as_u64(&raw, &dt)?)
@@ -496,6 +511,33 @@ impl<'f> MmapDataset<'f> {
.transpose()
}
/// A chunked dataset that stores `T` natively, decoded straight into a
/// `Vec<T>` (no byte buffer to convert); `None` for any other dataset
/// (see [`data_read::read_chunked_native`]).
fn read_chunked_native<T: data_read::NativeElement>(&self) -> Result<Option<Vec<T>>, Error> {
let dl = self.data_layout()?;
if !matches!(dl, DataLayout::Chunked { .. }) {
return Ok(None);
}
let dt = self.datatype()?;
if !T::is_native(&dt) {
return Ok(None);
}
let ds = self.dataspace()?;
let pipeline = self.filter_pipeline()?;
Ok(data_read::read_chunked_native::<T>(
&self.header.messages,
self.file.hdf5_bytes(),
&dl,
&ds,
&dt,
pipeline.as_ref(),
self.file.offset_size(),
self.file.length_size(),
None,
)?)
}
fn read_raw(&self) -> Result<Vec<u8>, Error> {
let dt = self.datatype()?;
let ds = self.dataspace()?;
+43
View File
@@ -507,6 +507,9 @@ impl<'f> Dataset<'f> {
if let Ok(Some(bytes)) = self.read_raw_ref() {
return Ok(data_read::read_as_f64(bytes, &dt)?);
}
if let Some(values) = self.read_chunked_native::<f64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
Ok(data_read::read_as_f64(&raw, &dt)?)
}
@@ -524,6 +527,9 @@ impl<'f> Dataset<'f> {
if let Ok(Some(bytes)) = self.read_raw_ref() {
return Ok(data_read::read_as_f32(bytes, &dt)?);
}
if let Some(values) = self.read_chunked_native::<f32>()? {
return Ok(values);
}
let raw = self.read_raw()?;
Ok(data_read::read_as_f32(&raw, &dt)?)
}
@@ -536,6 +542,9 @@ impl<'f> Dataset<'f> {
if let Ok(Some(bytes)) = self.read_raw_ref() {
return Ok(data_read::read_as_i32(bytes, &dt)?);
}
if let Some(values) = self.read_chunked_native::<i32>()? {
return Ok(values);
}
let raw = self.read_raw()?;
Ok(data_read::read_as_i32(&raw, &dt)?)
}
@@ -548,6 +557,9 @@ impl<'f> Dataset<'f> {
if let Ok(Some(bytes)) = self.read_raw_ref() {
return Ok(data_read::read_as_i64(bytes, &dt)?);
}
if let Some(values) = self.read_chunked_native::<i64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
Ok(data_read::read_as_i64(&raw, &dt)?)
}
@@ -560,6 +572,9 @@ impl<'f> Dataset<'f> {
if let Ok(Some(bytes)) = self.read_raw_ref() {
return Ok(data_read::read_as_u64(bytes, &dt)?);
}
if let Some(values) = self.read_chunked_native::<u64>()? {
return Ok(values);
}
let raw = self.read_raw()?;
Ok(data_read::read_as_u64(&raw, &dt)?)
}
@@ -1050,6 +1065,34 @@ impl<'f> Dataset<'f> {
.transpose()
}
/// A chunked dataset that stores `T` natively, decoded straight into a
/// `Vec<T>` through the file's chunk cache (no byte buffer to convert);
/// `None` for any other dataset (see
/// [`data_read::read_chunked_native`]).
fn read_chunked_native<T: data_read::NativeElement>(&self) -> Result<Option<Vec<T>>, Error> {
let dl = self.data_layout()?;
if !matches!(dl, DataLayout::Chunked { .. }) {
return Ok(None);
}
let dt = self.datatype()?;
if !T::is_native(&dt) {
return Ok(None);
}
let ds = self.dataspace()?;
let pipeline = self.filter_pipeline()?;
Ok(data_read::read_chunked_native::<T>(
&self.header.messages,
self.file.data.as_bytes(),
&dl,
&ds,
&dt,
pipeline.as_ref(),
self.file.offset_size(),
self.file.length_size(),
Some(&self.file.chunk_cache),
)?)
}
fn read_raw(&self) -> Result<Vec<u8>, Error> {
let dt = self.datatype()?;
let ds = self.dataspace()?;
+112
View File
@@ -0,0 +1,112 @@
//! Full reads of chunked datasets must not wait for a busy rayon pool of
//! any size.
//!
//! A full read handed its chunks to the rayon pool (`par_iter`) and the
//! calling thread — not a pool worker — slept until the pool had decoded
//! them. With a small pool (2-4 threads) and more reading threads than
//! workers, every reader queued behind the same few workers
//! (`docs/known-issues.md`, "Concurrent and contiguous read performance").
//! Now the calling thread decodes too, and pool workers only help when they
//! are free. The test keeps both workers of a two-thread pool busy and
//! requires reads to finish anyway, with the right values.
//!
//! One test in its own binary: it configures the process-wide rayon pool.
#![cfg(feature = "parallel")]
use std::sync::mpsc;
use std::time::Duration;
use clawhdf5::{File, FileBuilder};
const N: usize = 4096; // 64 chunks of 64 elements
fn values() -> Vec<f64> {
(0..N).map(|i| i as f64 * 0.25 - 7.0).collect()
}
fn build() -> File {
let mut b = FileBuilder::new();
b.create_dataset("data")
.with_f64_data(&values())
.with_shape(&[N as u64])
.with_chunks(&[64])
.with_deflate(1)
.with_provenance("test-suite", "2026-09-26T00:00:00Z", None);
File::from_bytes(b.finish().unwrap()).unwrap()
}
/// Run `f` on a fresh thread; `None` if it has not finished within `limit`.
fn finishes_within<T: Send + 'static>(
limit: Duration,
f: impl FnOnce() -> T + Send + 'static,
) -> Option<T> {
let (tx, rx) = mpsc::channel();
std::thread::spawn(move || {
let _ = tx.send(f());
});
rx.recv_timeout(limit).ok()
}
#[test]
fn full_reads_do_not_wait_for_a_busy_small_pool() {
const WORKERS: usize = 2;
rayon::ThreadPoolBuilder::new()
.num_threads(WORKERS)
.build_global()
.expect("this test binary configures the global pool first");
// Built first: the writer compresses on the pool too.
let file = std::sync::Arc::new(build());
// Occupy every worker of the pool until the reads are done.
let (started_tx, started_rx) = mpsc::channel();
let (release_tx, release_rx) = mpsc::channel::<()>();
let release_rx = std::sync::Arc::new(std::sync::Mutex::new(release_rx));
for _ in 0..WORKERS {
let (started_tx, release_rx) = (started_tx.clone(), release_rx.clone());
rayon::spawn(move || {
started_tx.send(()).unwrap();
let _ = release_rx.lock().unwrap().recv();
});
}
for _ in 0..WORKERS {
started_rx.recv().unwrap();
}
let limit = Duration::from_secs(20);
// Several readers at once, as in the `concurrent_read` benchmark: the
// cached full read (`read_*`), the typed one and the uncached reader
// behind `verify_provenance`.
let readers: Vec<_> = (0..4)
.map(|_| {
let file = std::sync::Arc::clone(&file);
std::thread::spawn(move || {
finishes_within(limit, move || {
let ds = file.dataset("data").unwrap();
(
ds.read_f64().unwrap(),
ds.read_f32().unwrap(),
ds.verify_provenance().unwrap(),
)
})
})
})
.collect();
let results: Vec<_> = readers.into_iter().map(|h| h.join().unwrap()).collect();
// Free the workers before asserting, so a failure does not hang the
// blocked reader threads forever.
for _ in 0..WORKERS {
release_tx.send(()).unwrap();
}
let want = values();
let want_f32: Vec<f32> = want.iter().map(|&v| v as f32).collect();
for result in results {
let (f64s, f32s, verified) =
result.expect("a full read waited for the busy two-thread rayon pool");
assert_eq!(f64s, want);
assert_eq!(f32s, want_f32);
assert_eq!(verified, clawhdf5::provenance::VerifyResult::Ok);
}
}
@@ -0,0 +1,361 @@
//! Every full and selection read path of chunked datasets against h5py.
//!
//! h5py (libhdf5) writes chunked datasets of every numeric type the typed
//! readers cover, in both byte orders, through deflate, shuffle,
//! Fletcher32, LZF, SZIP and Blosc, in 1-3 dimensional shapes whose chunks
//! do not divide them (partial edge chunks), plus sparse datasets whose
//! unwritten chunks read as a fill value. Next to each it stores the values
//! as contiguous `f64`, and checks that h5py reads the chunked dataset back
//! as those values.
//!
//! clawhdf5 must read every chunked dataset as those values through every
//! reader: `File` (the chunk-cached reader, twice so the second read can hit
//! the cache, and the typed readers that decode straight into their output),
//! `File::from_bytes`, `MmapFile`, `LazyFile`, and the selection readers
//! (a small hyperslab, a strided one covering most of the dataset, points).
//! With `--features parallel` the same reads decode chunks on several
//! threads. A filter this build does not include must be an error, never
//! data.
//!
//! Skipped when python3 with h5py is unavailable, unless
//! `CLAWHDF5_REQUIRE_INTEROP=1`.
use std::process::Command;
use clawhdf5::{File, LazyFile, MmapFile, Selection};
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")
}
/// Whether the interop test can run; panics instead of skipping when
/// `CLAWHDF5_REQUIRE_INTEROP=1`.
fn have_python() -> bool {
let ok = Command::new(python())
.args(["-c", "import h5py, numpy"])
.output()
.map(|o| o.status.success())
.unwrap_or(false);
if ok {
return true;
}
assert!(
!interop_required(),
"CLAWHDF5_REQUIRE_INTEROP=1 but python3 with h5py is not available"
);
eprintln!("SKIP: python3 with h5py not available");
false
}
/// Writes `c/<name>` (chunked) and `e/<name>` (the expected values, `f64`,
/// contiguous) for every case; prints one line per case:
/// `name filter_ids(comma-separated or -)`.
const GENERATE: &str = r#"
import sys
import numpy as np, h5py
try:
import hdf5plugin
except ImportError:
hdf5plugin = None
path = sys.argv[1]
dtypes = []
for code in ['i1', 'u1', 'i2', 'u2', 'i4', 'u4', 'i8', 'u8', 'f2', 'f4', 'f8']:
orders = ['|'] if code[1] == '1' else ['<', '>']
dtypes += [np.dtype(o + code) for o in orders]
filters = {
'none': {},
'gzip': dict(compression='gzip', compression_opts=4),
'shuffle_gzip': dict(shuffle=True, compression='gzip', compression_opts=1),
# h5py puts the checksum last (applied after deflate).
'shuffle_gzip_fletcher': dict(shuffle=True, compression='gzip', fletcher32=True),
'fletcher': dict(fletcher32=True),
'shuffle_lzf': dict(shuffle=True, compression='lzf'),
}
if h5py.h5z.filter_avail(h5py.h5z.FILTER_SZIP):
filters['szip'] = dict(compression='szip', compression_opts=('nn', 8))
if hdf5plugin is not None:
filters['blosc'] = dict(**hdf5plugin.Blosc(cname='lz4', clevel=5,
shuffle=hdf5plugin.Blosc.SHUFFLE))
shapes = [((37, 23), (8, 5)), ((101,), (16,)), ((9, 10, 11), (4, 3, 5))]
def values(n, dt):
i = np.arange(n, dtype=np.int64)
if dt.kind == 'f':
v = ((i * 7 + 3) % 1000) / 8.0 - 60.0
elif dt.kind == 'i':
v = (i * 7 + 3) % 200 - 100
else:
v = (i * 7 + 3) % 250
return v.astype(dt)
def filter_ids(dset):
plist = dset.id.get_create_plist()
ids = [str(plist.get_filter(k)[0]) for k in range(plist.get_nfilters())]
return ','.join(ids) or '-'
with h5py.File(path, 'w') as f:
n = 0
for dt in dtypes:
for fname, fopts in filters.items():
for s, (shape, chunks) in enumerate(shapes):
name = f"{dt.str.replace('|', 'x').replace('<', 'le').replace('>', 'be')}_{fname}_{s}"
data = values(int(np.prod(shape)), dt).reshape(shape)
try:
d = f.create_dataset('c/' + name, data=data, chunks=chunks, **fopts)
except (ValueError, TypeError) as e:
continue # a filter that refuses this type
f.create_dataset('e/' + name, data=data.astype('<f8'))
assert np.array_equal(d[()], data), name
print(name, filter_ids(d))
n += 1
# Larger than the file's 16 MiB chunk cache, so `File` decodes into
# its reusable buffers instead of caching chunks.
if dt.str in ('<f4', '>i4'):
name = f"{dt.str.replace('<', 'le').replace('>', 'be')}_large"
shape = (2600, 2048)
data = values(shape[0] * shape[1], dt).reshape(shape)
d = f.create_dataset('c/' + name, data=data, chunks=(256, 256),
shuffle=True, compression='gzip', compression_opts=1)
f.create_dataset('e/' + name, data=data.astype('<f8'))
assert np.array_equal(d[()], data), name
print(name, filter_ids(d))
# Sparse: only some chunks written; the rest read as the fill value
# (non-default, and default 0).
for fill in [None, 42]:
for fname in ['none', 'shuffle_gzip']:
name = f"{dt.str.replace('|', 'x').replace('<', 'le').replace('>', 'be')}_sparse_{fname}_{fill}"
shape, chunks = (37, 23), (8, 5)
kw = dict(filters[fname])
if fill is not None:
kw['fillvalue'] = np.array(fill, dtype=dt)
d = f.create_dataset('c/' + name, shape=shape, dtype=dt, chunks=chunks, **kw)
full = values(37 * 23, dt).reshape(shape)
d[3:17, 4:12] = full[3:17, 4:12]
d[30:, 20:] = full[30:, 20:]
expect = d[()]
want = np.full(shape, 0 if fill is None else fill, dtype=dt)
want[3:17, 4:12] = full[3:17, 4:12]
want[30:, 20:] = full[30:, 20:]
assert np.array_equal(expect, want), name
f.create_dataset('e/' + name, data=want.astype('<f8'))
print(name, filter_ids(d))
"#;
/// A filter this build can decode.
fn filters_available(ids: &str) -> bool {
ids == "-"
|| ids.split(',').all(|id| {
clawhdf5_format::filter_registry::is_filter_available(id.parse().expect("filter id"))
})
}
/// The expected values converted as libhdf5 converts them for each typed
/// reader (every case's values are exact in `f32` and within `i32`).
struct Expected {
f64s: Vec<f64>,
}
impl Expected {
fn f32s(&self) -> Vec<f32> {
self.f64s.iter().map(|&v| v as f32).collect()
}
/// Truncation toward zero; negative values read as unsigned are 0.
fn i32s(&self) -> Vec<i32> {
self.f64s.iter().map(|&v| v as i32).collect()
}
fn i64s(&self) -> Vec<i64> {
self.f64s.iter().map(|&v| v as i64).collect()
}
fn u64s(&self) -> Vec<u64> {
self.f64s.iter().map(|&v| v as u64).collect()
}
fn select(&self, idx: &[usize]) -> Expected {
Expected {
f64s: idx.iter().map(|&i| self.f64s[i]).collect(),
}
}
}
/// The typed readers' results for one dataset, compared with `want`.
macro_rules! check_typed {
($ds:expr, $want:expr, $name:expr, $path:expr) => {{
let (ds, want, name, path) = (&$ds, &$want, $name, $path);
assert_eq!(ds.read_f64().unwrap(), want.f64s, "{name} {path} f64");
assert_eq!(ds.read_f32().unwrap(), want.f32s(), "{name} {path} f32");
assert_eq!(ds.read_i32().unwrap(), want.i32s(), "{name} {path} i32");
assert_eq!(ds.read_i64().unwrap(), want.i64s(), "{name} {path} i64");
assert_eq!(ds.read_u64().unwrap(), want.u64s(), "{name} {path} u64");
}};
}
/// Row-major indices of the elements `sel` picks from a dataset of `dims`.
fn selected(sel: &Selection, dims: &[u64]) -> Vec<usize> {
let strides: Vec<usize> = (0..dims.len())
.map(|d| dims[d + 1..].iter().product::<u64>() as usize)
.collect();
match sel {
Selection::Hyperslab {
start,
stride,
count,
block,
} => {
// Per dimension, the selected coordinates in order.
let axes: Vec<Vec<u64>> = (0..dims.len())
.map(|d| {
(0..count[d])
.flat_map(|c| (0..block[d]).map(move |b| start[d] + c * stride[d] + b))
.collect()
})
.collect();
let mut out = vec![0usize];
for (d, axis) in axes.iter().enumerate() {
let stride = strides[d];
out = out
.iter()
.flat_map(|&base| axis.iter().map(move |&x| base + x as usize * stride))
.collect();
}
out
}
Selection::Points(points) => points
.iter()
.map(|p| p.iter().zip(&strides).map(|(&x, &s)| x as usize * s).sum())
.collect(),
_ => unreachable!(),
}
}
/// A small box (read through the partial-read path), a strided selection
/// covering most of the dataset (read in full, then selected), and points.
fn selections(dims: &[u64]) -> Vec<Selection> {
let rank = dims.len();
let small = Selection::Hyperslab {
start: dims.iter().map(|&d| d / 3).collect(),
stride: vec![1; rank],
count: dims.iter().map(|&d| (d / 4).max(1)).collect(),
block: vec![1; rank],
};
let strided = Selection::Hyperslab {
start: vec![0; rank],
stride: vec![2; rank],
count: dims.iter().map(|&d| d.div_ceil(2)).collect(),
block: vec![1; rank],
};
let points = Selection::Points(vec![
vec![0; rank],
dims.iter().map(|&d| d - 1).collect(),
dims.iter().map(|&d| d / 2).collect(),
dims.iter().map(|&d| (d * 2) / 3).collect(),
]);
vec![small, strided, points]
}
#[test]
fn chunked_reads_match_h5py_on_every_path() {
if !have_python() {
return;
}
let dir = tempfile::tempdir().unwrap();
let path = dir.path().join("chunked_paths.h5");
let out = Command::new(python())
.arg("-c")
.arg(GENERATE)
.arg(&path)
.output()
.expect("run python");
assert!(
out.status.success(),
"generator failed: {}",
String::from_utf8_lossy(&out.stderr)
);
let cases: Vec<(String, String)> = String::from_utf8(out.stdout)
.unwrap()
.lines()
.map(|l| {
let (name, ids) = l.split_once(' ').unwrap();
(name.to_string(), ids.to_string())
})
.collect();
assert!(cases.len() > 300, "only {} cases", cases.len());
if interop_required() {
// The generator must have covered the plugin filters too.
for f in ["_szip_", "_blosc_", "_shuffle_lzf_", "_sparse_"] {
assert!(cases.iter().any(|(n, _)| n.contains(f)), "no {f} case");
}
}
let file = File::open(&path).unwrap();
let owned = File::from_bytes(std::fs::read(&path).unwrap()).unwrap();
let mmap = MmapFile::open(&path).unwrap();
let lazy = LazyFile::open_mmap(&path).unwrap();
let (mut checked, mut refused) = (0, 0);
for (name, ids) in &cases {
let chunked = format!("c/{name}");
let want = Expected {
f64s: file
.dataset(&format!("e/{name}"))
.unwrap()
.read_f64()
.unwrap(),
};
let ds = file.dataset(&chunked).unwrap();
if !filters_available(ids) {
// Unsupported filter: an error from every reader, never wrong
// data. (An optional filter, as Blosc is, may have declined every
// chunk — they are then stored as-is and read fine.)
let results = [
ds.read_f64(),
owned.dataset(&chunked).unwrap().read_f64(),
mmap.dataset(&chunked).unwrap().read_f64(),
lazy.dataset(&chunked).unwrap().read_f64(),
];
let errors = results.iter().filter(|r| r.is_err()).count();
assert!(errors == 0 || errors == results.len(), "{name}");
for values in results.into_iter().flatten() {
assert_eq!(values, want.f64s, "{name}");
}
if errors > 0 {
assert!(ds.read_f32().is_err(), "{name}");
refused += 1;
continue;
}
}
// Twice: the second read of a small dataset comes from the cache.
check_typed!(ds, want, name, "File");
check_typed!(ds, want, name, "File (cached)");
check_typed!(owned.dataset(&chunked).unwrap(), want, name, "from_bytes");
check_typed!(mmap.dataset(&chunked).unwrap(), want, name, "MmapFile");
check_typed!(lazy.dataset(&chunked).unwrap(), want, name, "LazyFile");
let dims = ds.shape().unwrap();
for sel in selections(&dims) {
let want = want.select(&selected(&sel, &dims));
assert_eq!(
ds.read_f64_selection(&sel).unwrap(),
want.f64s,
"{name} {sel:?}"
);
assert_eq!(
ds.read_f32_selection(&sel).unwrap(),
want.f32s(),
"{name} {sel:?}"
);
assert_eq!(
ds.read_i64_selection(&sel).unwrap(),
want.i64s(),
"{name} {sel:?}"
);
}
checked += 1;
}
eprintln!("{checked} datasets read on every path, {refused} refused (filter not built in)");
assert!(checked > 300);
}
+16 -5
View File
@@ -26,8 +26,9 @@ performance" below.
## Concurrent and contiguous read performance (measured 2026-09-26)
**Status:** open for chunked full reads (one cause fixed 2026-09-26); the
contiguous item is fixed (2026-09-26). Measured on
**Status:** open for chunked full reads at 16 threads until re-measured
(the causes identified below are fixed as of 2026-09-26); the contiguous
item is fixed (2026-09-26). Measured on
tank with `concurrent_read` against h5py 3.16 / HDF5 2.0 (`BENCHMARKS.md`,
"Concurrent reads"):
- **Partly fixed 2026-09-26.** Full reads of chunked datasets from several threads
@@ -48,9 +49,19 @@ tank with `concurrent_read` against h5py 3.16 / HDF5 2.0 (`BENCHMARKS.md`,
readers outside it still wait on its workers. Datasets larger than the
cache's budget were already read without inserting into it, and skipping
its lookups entirely gained only a few percent at 16 threads. Remaining
per-read overhead, not yet addressed: each full `read_f32` of a chunked
dataset faults in about three times its size in fresh pages (the output,
the `f32` copy of it, and a new buffer per decoded chunk).
per-read overhead: each full `read_f32` of a chunked dataset faults in
about three times its size in fresh pages (the output, the `f32` copy of
it, and a new buffer per decoded chunk).
**Fixed 2026-09-26** (both causes; see `CHANGELOG.md`, "Chunked full
reads"): chunks are decoded into buffers each thread reuses and copied
straight into the output, and the typed readers decode into the `Vec<T>`
they return, so a full read no longer faults in a buffer per chunk or a
second copy of its output; and the reading
thread decodes its own chunks with pool workers helping when free, so no
reader waits on a small or busy pool
(`crates/clawhdf5/tests/busy_decode_pool.rs`). The 16-thread comparison
with h5py processes has not been re-measured yet (tank was busy with
other work); this item stays open until it is.
- Contiguous datasets read 4x slower than h5py on one thread (2.5 vs
9.8 GB/s full, 0.12x for 256 x 256 hyperslabs).
**Fixed 2026-09-26** (re-measured on tank at `408f69e`: 13665 MB/s