crates/clawhdf5-format/src/data_read.rs

@@ -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

crates/clawhdf5/src/lazy.rs

@@ -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()?;

crates/clawhdf5/src/mmap_file.rs

@@ -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()?;

crates/clawhdf5/src/reader.rs

@@ -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()?;

crates/clawhdf5/tests/chunked_read_paths_interop.rs

@@ -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);
+}