diff --git a/CHANGELOG.md b/CHANGELOG.md index d679fbf..f19acc4 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -60,6 +60,15 @@ is now guarded and a panic becomes `clawhdf5.InternalError` (a `RuntimeError`) naming the object; with the implicit-index panic above restored, `ds[0:30]` raises it. +- **Wrong data: uninitialised padding in compound results of index lists.** + `ds[[0, 3, 6]]` joined one read per run with `np.concatenate`, which + copies structured dtypes field by field into an `np.empty` result, so the + padding bytes held whatever was in memory (pointers were seen) and leaked + through `tobytes()`, hashes and write-backs. The runs' bytes are now joined + in Rust, whole elements at a time, so the result carries the bytes read + from the file (h5py's, zero for files it wrote) and stays zero-copy. + The h5py comparisons now also compare every byte of structured values + (`test_compound_padding_bytes_match_h5py` and `assert_same`). - **CI builds and tests the Python package.** It was excluded from CI. `scripts/ci-test.sh` now lints `clawhdf5-py`, builds the wheel with maturin, unpacks it under `target/` and runs the pytest suite; skipped diff --git a/crates/clawhdf5-py/src/dataset.rs b/crates/clawhdf5-py/src/dataset.rs index 9109a32..36467d8 100644 --- a/crates/clawhdf5-py/src/dataset.rs +++ b/crates/clawhdf5-py/src/dataset.rs @@ -82,70 +82,60 @@ impl PyDataset { conv.empty(py, &out_shape)? } else { let (reads, list_axis) = plan.reads(dims); + let read_shape = plan.read_shape(); let file = &*self.file; let path = self.path.as_str(); let (vl, elem_size, unit) = (conv.is_vl(), conv.elem_size, conv.vl_unit); // Everything below touches only Rust data: release the GIL. - let read = || -> Result<_, ReadError> { + let read = || -> Result { let ds = file.dataset(path)?; - let sb = file.superblock(); let mut blocks = Vec::with_capacity(reads.len()); for (sel, shape) in reads { let raw = ds.read_selection(&sel)?; - let n: usize = shape.iter().product(); - let data = if vl { - if raw.len() != n * elem_size { - return Err(ReadError::Other(format!( - "read {} bytes of variable-length references, expected {}", - raw.len(), - n * elem_size - ))); - } - Elements::Vl( - resolve_vl( - file.as_bytes(), - &raw, - n, - sb.offset_size, - sb.length_size, - unit, - ) - .map_err(ReadError::Other)?, - ) - } else { - Elements::Bytes(raw) - }; - blocks.push((data, shape)); + let want = shape.iter().product::() * elem_size; + if raw.len() != want { + return Err(ReadError::Other(format!( + "read {} bytes, expected {want}", + raw.len() + ))); + } + blocks.push((raw, shape)); } - Ok(blocks) + // Several blocks only for a list index: join their bytes + // (every byte of every element, padding included) along + // that axis before anything becomes numpy. + let raw = match (blocks.len(), list_axis) { + (1, _) => blocks.pop().expect("one block").0, + (_, Some(axis)) => select::join_along(&blocks, axis, elem_size), + _ => { + return Err(ReadError::Other( + "several reads without an index list".into(), + )); + } + }; + if !vl { + return Ok(Elements::Bytes(raw)); + } + let sb = file.superblock(); + let n = read_shape.iter().product(); + resolve_vl( + file.as_bytes(), + &raw, + n, + sb.offset_size, + sb.length_size, + unit, + ) + .map(Elements::Vl) + .map_err(ReadError::Other) }; - let blocks: Vec<(Elements, Vec)> = py + let data = py .detach(|| { std::panic::catch_unwind(std::panic::AssertUnwindSafe(read)) .unwrap_or_else(|p| Err(ReadError::Panic(crate::panic_text(&*p)))) }) .map_err(|e| e.into_py(&self.path))?; - - let mut arrays = Vec::with_capacity(blocks.len()); - for (data, shape) in blocks { - arrays.push(conv.to_array(py, data, &shape, false)?); - } - let joined = if arrays.len() == 1 { - arrays.pop().expect("one block") - } else { - let axis = list_axis.expect("several reads only for a list index"); - // Name the dtype: left to itself numpy canonicalises a - // structured dtype here (drops padding, native byte order). - let kwargs = pyo3::types::PyDict::new(py); - kwargs.set_item("axis", axis)?; - kwargs.set_item("dtype", arrays[0].getattr("dtype")?)?; - kwargs.set_item("casting", "no")?; - py.import("numpy")?.call_method( - "concatenate", - (PyList::new(py, arrays)?,), - Some(&kwargs), - )? - }; + let joined = conv.to_array(py, data, &read_shape, false)?; // Drop the axes indexed by an integer (length 1 in the blocks). let mut shape = out_shape.clone(); if let crate::convert::Layout::Subarray(sub) = &conv.layout { diff --git a/crates/clawhdf5-py/src/select.rs b/crates/clawhdf5-py/src/select.rs index 1b1c796..dac62ab 100644 --- a/crates/clawhdf5-py/src/select.rs +++ b/crates/clawhdf5-py/src/select.rs @@ -56,6 +56,12 @@ impl Plan { .collect() } + /// The shape of the result before the integer-indexed axes are dropped + /// (they have length 1 here): the shape of the joined reads. + pub fn read_shape(&self) -> Vec { + self.axes.iter().map(|a| a.len() as usize).collect() + } + /// Whether the selection is empty. pub fn is_empty(&self) -> bool { self.axes.iter().any(|a| a.len() == 0) @@ -64,7 +70,7 @@ impl Plan { /// The hyperslab reads that make up this selection, each with the shape /// of its block (index axes kept at length 1). More than one only when an /// axis is indexed by a list: one read per run of consecutive indices, - /// concatenated along `list_axis` afterwards. + /// joined along `list_axis` afterwards (`join_along`). pub fn reads(&self, dims: &[u64]) -> (Vec<(Selection, Vec)>, Option) { let list_axis = self.axes.iter().position(|a| matches!(a, Axis::List(_))); let runs: Vec<(u64, u64)> = match list_axis.map(|i| &self.axes[i]) { @@ -119,6 +125,33 @@ fn consecutive_runs(idx: &[u64]) -> Vec<(u64, u64)> { runs } +/// Join row-major blocks of `elem_size`-byte elements whose shapes differ +/// only along `axis` into one buffer, in order along that axis. Whole +/// elements are copied, so compound padding keeps the bytes that were read. +pub(crate) fn join_along( + blocks: &[(Vec, Vec)], + axis: usize, + elem_size: usize, +) -> Vec { + let Some((_, first)) = blocks.first() else { + return Vec::new(); + }; + let outer: usize = first[..axis].iter().product(); + let inner: usize = first[axis + 1..].iter().product::() * elem_size; + let total: usize = blocks.iter().map(|(_, s)| s[axis]).sum(); + let mut out = vec![0u8; outer * total * inner]; + let mut at = 0; + for (bytes, shape) in blocks { + let len = shape[axis] * inner; + for o in 0..outer { + let dst = (o * total) * inner + at; + out[dst..dst + len].copy_from_slice(&bytes[o * len..(o + 1) * len]); + } + at += len; + } + out +} + /// Parse `key` for a dataset of shape `dims`. pub(crate) fn parse(key: &Bound<'_, PyAny>, dims: &[u64]) -> PyResult { let items: Vec> = match key.cast::() { @@ -311,6 +344,18 @@ mod tests { assert_eq!(consecutive_runs(&[]), vec![]); } + #[test] + fn blocks_join_along_the_list_axis() { + // Two 2x1 and 2x2 blocks of 1-byte elements, joined along axis 1. + let a = (vec![1, 2], vec![2, 1]); + let b = (vec![3, 4, 5, 6], vec![2, 2]); + assert_eq!(join_along(&[a, b], 1, 1), vec![1, 3, 4, 2, 5, 6]); + // Along axis 0 it is concatenation; 2-byte elements stay whole. + let a = (vec![1, 2, 3, 4], vec![1, 2]); + let b = (vec![5, 6, 7, 8], vec![1, 2]); + assert_eq!(join_along(&[a, b], 0, 2), vec![1, 2, 3, 4, 5, 6, 7, 8]); + } + #[test] fn full_selection_reads_everything() { let plan = Plan { diff --git a/crates/clawhdf5-py/tests/test_read_vs_h5py.py b/crates/clawhdf5-py/tests/test_read_vs_h5py.py index 66fc643..c5876d2 100644 --- a/crates/clawhdf5-py/tests/test_read_vs_h5py.py +++ b/crates/clawhdf5-py/tests/test_read_vs_h5py.py @@ -212,14 +212,19 @@ def assert_same(ours, theirs, what=""): if theirs.dtype == object: for a, b in zip(ours.ravel(), theirs.ravel()): assert_same(a, b, what) - elif theirs.dtype.names is None and theirs.dtype.kind == "V": - assert ours.tobytes() == theirs.tobytes(), what + elif theirs.dtype.kind == "V": + # Structured and opaque: every byte, padding included (h5py's + # padding is zero; uninitialised memory there would leak). + if theirs.dtype.names is not None: + np.testing.assert_array_equal(ours, theirs, err_msg=what) + assert ours.tobytes() == theirs.tobytes(), f"{what}: bytes differ" else: np.testing.assert_array_equal(ours, theirs, err_msg=what) elif isinstance(theirs, np.generic): assert ours.dtype == theirs.dtype, what if theirs.dtype.names is not None: np.testing.assert_array_equal(np.asarray(ours), np.asarray(theirs), err_msg=what) + assert ours.tobytes() == theirs.tobytes(), f"{what}: bytes differ" else: assert ours == theirs or (ours != ours and theirs != theirs), what else: @@ -511,3 +516,18 @@ def test_a_library_panic_is_an_ordinary_exception(): clawhdf5._panic_for_test() except Exception: # noqa: BLE001 - the point of the test pass + + +def test_compound_padding_bytes_match_h5py(pair): + """Every byte of a padded compound, padding included, is h5py's, for + index lists with many runs as well as slices. Joining the runs with + np.concatenate left the padding uninitialised: process memory ended up + in tobytes().""" + ours, theirs, _ = pair + keys = [[0, 3, 6], [1, 2, 5, 9], [0, 2, 4, 6, 8], slice(None), slice(1, 9, 3), 4, [9]] + for name in ["cmp/padded", "cmp/padded_chunked"]: + for _ in range(20): # garbage varies between runs; zeros do not + for key in keys: + assert ours[name][key].tobytes() == theirs[name][key].tobytes(), f"{name}[{key!r}]" + for key in [(slice(None), [0, 2]), ([0, 2, 3], slice(None)), ([1, 3], 1)]: + assert ours["cmp/nested_2d"][key].tobytes() == theirs["cmp/nested_2d"][key].tobytes(), key