From 190918a478bfc8bc412f9c58e5a70b087e560d0b Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:57:32 -0500 Subject: [PATCH 1/4] feat(format): decode hyperslab selection versions 1 and 2 in VDS mappings libhdf5 serializes a VDS hyperslab as version 1 (irregular, 4-byte block corners) for the default format bounds, and as version 2 (regular, 8-byte) for unlimited selections in the 1.10 format. Only version 3 was accepted, so every h5py VDS written with default libver failed with "only version-3 hyperslab selections are supported" (5 libhdf5 test files in the sweep). Decode all three versions following H5S__hyper_deserialize, including irregular hyperslabs (a union of blocks, enumerated in row-major order as libhdf5 iterates them) and the all-ones "unlimited" count/block marker. SerializedSelection exposes the raw form for unlimited-mapping support. Test: vds_interop::vds_version1_irregular_hyperslab_selections compares default-libver h5py VDS reads (contiguous, strided and 2-D block mappings) with libhdf5's values. Co-Authored-By: Claude Opus 5.5 (1M context) --- CHANGELOG.md | 9 + crates/clawhdf5-format/src/selection.rs | 498 +++++++++++++++++++----- crates/clawhdf5/tests/vds_interop.rs | 156 ++++++++ docs/known-issues.md | 3 +- 4 files changed, 562 insertions(+), 104 deletions(-) create mode 100644 crates/clawhdf5/tests/vds_interop.rs diff --git a/CHANGELOG.md b/CHANGELOG.md index 8b87415..a349fba 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -273,6 +273,15 @@ - CI keeps zlib-ng building and tested; the arm64 job no longer needs cmake. ### Correctness +- `clawhdf5-format` virtual datasets (VDS), checked against HDF5 2.0 through + h5py (`crates/clawhdf5/tests/vds_interop.rs`): + - Hyperslab selection versions 1 and 2 were refused ("only version-3 + hyperslab selections are supported"). Version 1 is what libhdf5 writes for + every VDS created with the default format bounds (h5py's default), so + those could not be read at all; version 2 is its encoding of an unlimited + selection. Both are decoded now, as are irregular hyperslabs (a union of + blocks, read in row-major order as libhdf5 iterates them). + `SerializedSelection` exposes the raw form, including unlimited counts. - `clawhdf5-format` reader — **values returned wrong with no error:** - Fixed Array and Extensible Array chunk indexes were laid out by the dataset's current shape instead of its max shape (23 libhdf5 test files, diff --git a/crates/clawhdf5-format/src/selection.rs b/crates/clawhdf5-format/src/selection.rs index e1f98aa..de31e46 100644 --- a/crates/clawhdf5-format/src/selection.rs +++ b/crates/clawhdf5-format/src/selection.rs @@ -229,44 +229,47 @@ impl Selection { /// self-describing in length, so the count lets a caller walk a packed list /// of selections — as the Virtual Dataset global-heap block does). /// - /// Only the forms needed for VDS assembly are decoded: `ALL`, `NONE`, and - /// **regular** hyperslabs serialized at **version 3** (the encoding HDF5 - /// 1.10+/2.0 emit). Point selections, irregular hyperslabs, and older - /// hyperslab versions return an error rather than mis-decoding. + /// Decodes `ALL`, `NONE`, and hyperslabs at every version libhdf5 writes + /// (1: irregular, 4-byte coordinates — the default-format encoding; 2: + /// regular, 8-byte; 3: either, variable width). A regular hyperslab maps + /// to [`Selection::Hyperslab`]; an *irregular* one (a union of blocks) + /// maps to a single-block hyperslab when it has one block, and otherwise to + /// [`Selection::Points`] listing the union in row-major order (the order + /// libhdf5 iterates it in). Unlimited counts/blocks decode as `u64::MAX` + /// (see [`SerializedSelection::decode`] for the raw form). Point + /// selections are refused: libhdf5 does not allow them in virtual datasets + /// either. pub fn decode_serialized(data: &[u8]) -> Result<(Selection, usize), FormatError> { - if data.len() < 8 { - return Err(FormatError::UnexpectedEof { - expected: 8, - available: data.len(), - }); - } - let sel_type = u32::from_le_bytes([data[0], data[1], data[2], data[3]]); - let version = u32::from_le_bytes([data[4], data[5], data[6], data[7]]); - - match sel_type { - // ALL / NONE: type(4) + version(4) + reserved(4) + length(4) = 16 bytes. - 3 | 0 => { - if data.len() < 16 { - return Err(FormatError::UnexpectedEof { - expected: 16, - available: data.len(), - }); - } - let sel = if sel_type == 3 { - Selection::All + let (raw, len) = SerializedSelection::decode(data)?; + let sel = match raw { + SerializedSelection::All => Selection::All, + SerializedSelection::None => Selection::None, + SerializedSelection::Regular { + start, + stride, + count, + block, + } => Selection::Hyperslab { + start, + stride, + count, + block, + }, + SerializedSelection::Blocks { rank, starts, ends } => { + if starts.len() == rank { + let block = starts.iter().zip(&ends).map(|(&s, &e)| e - s + 1).collect(); + Selection::Hyperslab { + start: starts, + stride: vec![1; rank], + count: vec![1; rank], + block, + } } else { - Selection::None - }; - Ok((sel, 16)) + Selection::Points(blocks_union_coords(rank, &starts, &ends)?) + } } - 2 => decode_hyperslab_serialized(data, version), - 1 => Err(FormatError::ChunkedReadError( - "VDS point selections are not supported".into(), - )), - _ => Err(FormatError::ChunkedReadError( - "unknown dataspace selection type".into(), - )), - } + }; + Ok((sel, len)) } /// Enumerate the selected element indices of a **1-D** dataspace of the @@ -314,6 +317,11 @@ impl Selection { "VDS selection rank does not match dataspace rank".into(), )); } + if count.iter().chain(block.iter()).any(|&v| v == UNLIMITED) { + return Err(FormatError::ChunkedReadError( + "unlimited selection must be clipped before it is enumerated".into(), + )); + } // Selected coordinates along each dimension, in order. let mut per_dim: Vec> = Vec::with_capacity(rank); for d in 0..rank { @@ -400,84 +408,279 @@ impl Selection { } } -/// Decode an `H5S_SEL_HYPER` selection in its serialized form. Only version-3 -/// **regular** hyperslabs are supported. -fn decode_hyperslab_serialized( - data: &[u8], - version: u32, -) -> Result<(Selection, usize), FormatError> { - if version != 3 { - return Err(FormatError::ChunkedReadError( - "only version-3 hyperslab selections are supported".into(), - )); +/// Hyperslab count/block value meaning "unlimited" (`H5S_UNLIMITED`). +pub const UNLIMITED: u64 = u64::MAX; + +/// Largest number of elements an irregular selection is expanded to when it +/// is converted to a point list by [`Selection::decode_serialized`]. +const MAX_EXPANDED_POINTS: u64 = 1 << 26; + +/// A selection exactly as `H5S_select_serialize` stores it, before it is +/// applied to any dataspace. +/// +/// Unlike [`Selection`] this keeps an irregular hyperslab as its list of +/// blocks, and a regular hyperslab's count/block may be [`UNLIMITED`] (the +/// unlimited selections used by unlimited and "printf" virtual dataset +/// mappings). +#[derive(Debug, Clone, PartialEq)] +pub enum SerializedSelection { + /// `H5S_SEL_ALL`. + All, + /// `H5S_SEL_NONE`. + None, + /// A regular hyperslab. `count[d]` or `block[d]` may be [`UNLIMITED`]. + Regular { + start: Vec, + stride: Vec, + count: Vec, + block: Vec, + }, + /// An irregular hyperslab: the union of `starts.len() / rank` blocks, each + /// given by its first (`starts`) and last (`ends`, inclusive) coordinate, + /// flattened block-major. + Blocks { + rank: usize, + starts: Vec, + ends: Vec, + }, +} + +fn sel_err(msg: &str) -> FormatError { + FormatError::ChunkedReadError(msg.into()) +} + +/// Bounds-checked little-endian reader over a serialized selection. +struct SelReader<'a> { + data: &'a [u8], + pos: usize, +} + +impl SelReader<'_> { + fn take(&mut self, n: usize) -> Result<&[u8], FormatError> { + let end = self.pos.checked_add(n).filter(|&e| e <= self.data.len()); + let end = end.ok_or(FormatError::UnexpectedEof { + expected: self.pos.saturating_add(n), + available: self.data.len(), + })?; + let s = &self.data[self.pos..end]; + self.pos = end; + Ok(s) } - // type(4) ver(4) flags(1) enc_size(1) rank(4) [start,stride,count,block]*rank - if data.len() < 14 { - return Err(FormatError::UnexpectedEof { - expected: 14, - available: data.len(), - }); + + fn uint(&mut self, size: usize) -> Result { + let bytes = self.take(size)?; + Ok(bytes + .iter() + .enumerate() + .fold(0u64, |v, (i, &b)| v | (b as u64) << (i * 8))) } - let flags = data[8]; - let enc_size = data[9] as usize; - // Bit 0 set => regular hyperslab. Irregular hyperslabs list explicit blocks. - if flags & 0x01 == 0 { - return Err(FormatError::ChunkedReadError( - "irregular VDS hyperslab selections are not supported".into(), - )); + + fn remaining(&self) -> usize { + self.data.len() - self.pos } - if enc_size != 2 && enc_size != 4 && enc_size != 8 { - return Err(FormatError::ChunkedReadError( - "unsupported hyperslab coordinate encoding size".into(), - )); - } - let rank = u32::from_le_bytes([data[10], data[11], data[12], data[13]]) as usize; - // HDF5 caps dataspace rank at 32 (H5S_MAX_RANK). Reject anything larger so a - // corrupt rank can't drive a huge allocation or read loop. - if rank > 32 { - return Err(FormatError::ChunkedReadError( - "hyperslab selection rank exceeds maximum (32)".into(), - )); - } - let mut pos = 14; - let read_coord = |data: &[u8], pos: usize| -> Result { - if pos + enc_size > data.len() { - return Err(FormatError::UnexpectedEof { - expected: pos + enc_size, - available: data.len(), - }); +} + +impl SerializedSelection { + /// Decode a serialized selection, returning it and the number of bytes it + /// occupies. Mirrors libhdf5's `H5S_select_deserialize`: `ALL`/`NONE` and + /// hyperslab versions 1-3 are decoded; point selections (which libhdf5 + /// refuses in virtual datasets) and malformed input are errors. + pub fn decode(data: &[u8]) -> Result<(SerializedSelection, usize), FormatError> { + let mut r = SelReader { data, pos: 0 }; + let sel_type = r.uint(4)?; + let version = r.uint(4)?; + match sel_type { + // ALL / NONE: type(4) + version(4) + reserved(4) + length(4). + 0 | 3 => { + r.take(8)?; + let sel = if sel_type == 3 { + SerializedSelection::All + } else { + SerializedSelection::None + }; + Ok((sel, r.pos)) + } + 2 => { + let sel = decode_hyperslab(&mut r, version)?; + Ok((sel, r.pos)) + } + 1 => Err(sel_err( + "VDS point selections are not supported (libhdf5 rejects them too)", + )), + _ => Err(sel_err("unknown dataspace selection type")), } - let mut v = 0u64; - for (i, &b) in data[pos..pos + enc_size].iter().enumerate() { - v |= (b as u64) << (i * 8); + } + + /// The single dimension in which this selection is unlimited, if any. + pub fn unlimited_dim(&self) -> Option { + match self { + SerializedSelection::Regular { count, block, .. } => count + .iter() + .zip(block) + .position(|(&c, &b)| c == UNLIMITED || b == UNLIMITED), + _ => None, } - Ok(v) + } + + /// The rank the selection was serialized with (`None` for ALL/NONE, which + /// carry no rank). + pub fn rank(&self) -> Option { + match self { + SerializedSelection::Regular { start, .. } => Some(start.len()), + SerializedSelection::Blocks { rank, .. } => Some(*rank), + _ => None, + } + } +} + +/// `H5S__hyper_deserialize`: after the type and version words. +fn decode_hyperslab(r: &mut SelReader, version: u64) -> Result { + const REGULAR: u8 = 0x01; + let (flags, enc_size) = match version { + // v1: reserved(4) + length(4), always irregular, 4-byte coordinates. + 1 => { + r.take(8)?; + (0u8, 4usize) + } + // v2: flags(1) + length(4), 8-byte coordinates. + 2 => { + let flags = r.take(1)?[0]; + r.take(4)?; + (flags, 8) + } + // v3: flags(1) + encoding size(1). + 3 => { + let flags = r.take(1)?[0]; + let enc = r.take(1)?[0] as usize; + (flags, enc) + } + _ => return Err(sel_err("unsupported hyperslab selection version")), }; - let (mut start, mut stride, mut count, mut block) = ( - Vec::with_capacity(rank), - Vec::with_capacity(rank), - Vec::with_capacity(rank), - Vec::with_capacity(rank), - ); - for _ in 0..rank { - start.push(read_coord(data, pos)?); - pos += enc_size; - stride.push(read_coord(data, pos)?); - pos += enc_size; - count.push(read_coord(data, pos)?); - pos += enc_size; - block.push(read_coord(data, pos)?); - pos += enc_size; + if flags & !REGULAR != 0 { + return Err(sel_err("unknown hyperslab selection flags")); } - Ok(( - Selection::Hyperslab { + if !matches!(enc_size, 2 | 4 | 8) { + return Err(sel_err("unsupported hyperslab coordinate encoding size")); + } + let rank = r.uint(4)? as usize; + // HDF5 caps dataspace rank at 32 (H5S_MAX_RANK). Reject anything else so a + // corrupt rank can't drive a huge allocation or read loop. + if rank == 0 || rank > 32 { + return Err(sel_err("hyperslab selection rank must be 1..=32")); + } + // The all-ones value of the encoding width means "unlimited". + let unlim_raw = if enc_size == 8 { + u64::MAX + } else { + (1u64 << (enc_size * 8)) - 1 + }; + + if flags & REGULAR != 0 { + let (mut start, mut stride, mut count, mut block) = ( + Vec::with_capacity(rank), + Vec::with_capacity(rank), + Vec::with_capacity(rank), + Vec::with_capacity(rank), + ); + for _ in 0..rank { + start.push(r.uint(enc_size)?); + stride.push(r.uint(enc_size)?); + let c = r.uint(enc_size)?; + count.push(if c == unlim_raw { UNLIMITED } else { c }); + let b = r.uint(enc_size)?; + block.push(if b == unlim_raw { UNLIMITED } else { b }); + } + let unlimited = count + .iter() + .zip(&block) + .filter(|&(&c, &b)| c == UNLIMITED || b == UNLIMITED) + .count(); + if unlimited > 1 { + return Err(sel_err( + "hyperslab selection is unlimited in more than one dimension", + )); + } + for d in 0..rank { + // Overlapping blocks are not a valid regular hyperslab. + if count[d] > 1 && block[d] != UNLIMITED && block[d] > stride[d] { + return Err(sel_err("regular hyperslab blocks overlap")); + } + } + return Ok(SerializedSelection::Regular { start, stride, count, block, - }, - pos, - )) + }); + } + + // Irregular: number of blocks, then each block's start and end corners. + let nblocks = r.uint(enc_size)?; + let per_block = (rank * 2 * enc_size) as u64; + // Untrusted count: it must fit in what is left of the buffer. + if nblocks + .checked_mul(per_block) + .is_none_or(|need| need > r.remaining() as u64) + { + return Err(FormatError::UnexpectedEof { + expected: r + .pos + .saturating_add(nblocks.saturating_mul(per_block) as usize), + available: r.data.len(), + }); + } + let n = nblocks as usize * rank; + let (mut starts, mut ends) = (Vec::with_capacity(n), Vec::with_capacity(n)); + for _ in 0..nblocks { + for _ in 0..rank { + starts.push(r.uint(enc_size)?); + } + for _ in 0..rank { + ends.push(r.uint(enc_size)?); + } + } + if starts.iter().zip(&ends).any(|(s, e)| e < s) { + return Err(sel_err("hyperslab block ends before it starts")); + } + Ok(SerializedSelection::Blocks { rank, starts, ends }) +} + +/// The coordinates of the union of the given blocks, in row-major order. +fn blocks_union_coords( + rank: usize, + starts: &[u64], + ends: &[u64], +) -> Result>, FormatError> { + let mut total = 0u64; + for (s, e) in starts.chunks_exact(rank).zip(ends.chunks_exact(rank)) { + let vol = s + .iter() + .zip(e) + .try_fold(1u64, |acc, (&s, &e)| acc.checked_mul(e - s + 1)); + total = vol + .and_then(|v| total.checked_add(v)) + .filter(|&t| t <= MAX_EXPANDED_POINTS) + .ok_or_else(|| sel_err("irregular hyperslab selection is too large to expand"))?; + } + let mut out = Vec::with_capacity(total as usize); + for (s, e) in starts.chunks_exact(rank).zip(ends.chunks_exact(rank)) { + let mut cur = s.to_vec(); + 'block: loop { + out.push(cur.clone()); + for d in (0..rank).rev() { + if cur[d] < e[d] { + cur[d] += 1; + continue 'block; + } + cur[d] = s[d]; + } + break; + } + } + // Lexicographic order of coordinates is row-major order. + out.sort_unstable(); + out.dedup(); + Ok(out) } // --------------------------------------------------------------------------- @@ -642,11 +845,100 @@ mod tests { } #[test] - fn decode_irregular_hyperslab_rejected() { + fn decode_truncated_irregular_hyperslab_is_error() { + // Irregular, rank 1, but the block count is missing. let bytes = [0x02u8, 0, 0, 0, 0x03, 0, 0, 0, 0x00, 0x02, 0x01, 0, 0, 0]; assert!(Selection::decode_serialized(&bytes).is_err()); } + /// Version 1 as libhdf5 writes it for the default (earliest) format bounds: + /// type, version, reserved(4), length(4), rank(4), nblocks(4), then each + /// block's start and inclusive end corner as 4-byte values. + fn v1_blocks(rank: u32, blocks: &[(&[u32], &[u32])]) -> Vec { + let mut b = Vec::new(); + for w in [2u32, 1, 0, 0, rank, blocks.len() as u32] { + b.extend_from_slice(&w.to_le_bytes()); + } + for (s, e) in blocks { + for v in s.iter().chain(e.iter()) { + b.extend_from_slice(&v.to_le_bytes()); + } + } + b + } + + #[test] + fn decode_v1_irregular_single_block() { + // Exactly what h5py/HDF5 2.0 writes for `[0:4]` with default libver. + let bytes = v1_blocks(1, &[(&[0], &[3])]); + let (sel, used) = Selection::decode_serialized(&bytes).unwrap(); + assert_eq!(used, bytes.len()); + assert_eq!(sel.iter_linear_1d(8).unwrap(), vec![0, 1, 2, 3]); + } + + #[test] + fn decode_v1_irregular_union_is_row_major() { + // Blocks given out of order and overlapping still enumerate once each, + // in row-major order (libhdf5 iterates the union, not the list). + let bytes = v1_blocks(2, &[(&[1, 0], &[1, 1]), (&[0, 2], &[1, 2])]); + let (sel, used) = Selection::decode_serialized(&bytes).unwrap(); + assert_eq!(used, bytes.len()); + // (0,2) (1,0) (1,1) (1,2) in a 2x3 space. + assert_eq!(sel.iter_linear(&[2, 3]).unwrap(), vec![2, 3, 4, 5]); + } + + #[test] + fn decode_v2_regular_with_unlimited_count() { + // v2: flags(1) + length(4), then 8-byte start/stride/count/block. + let mut b = Vec::new(); + b.extend_from_slice(&2u32.to_le_bytes()); + b.extend_from_slice(&2u32.to_le_bytes()); + b.push(0x01); + b.extend_from_slice(&36u32.to_le_bytes()); + b.extend_from_slice(&1u32.to_le_bytes()); + for v in [0u64, 10, u64::MAX, 10] { + b.extend_from_slice(&v.to_le_bytes()); + } + let (raw, used) = SerializedSelection::decode(&b).unwrap(); + assert_eq!(used, b.len()); + assert_eq!(raw.unlimited_dim(), Some(0)); + assert_eq!( + raw, + SerializedSelection::Regular { + start: vec![0], + stride: vec![10], + count: vec![UNLIMITED], + block: vec![10], + } + ); + // An unclipped unlimited selection cannot be enumerated. + let (sel, _) = Selection::decode_serialized(&b).unwrap(); + assert!(sel.iter_linear_1d(100).is_err()); + } + + #[test] + fn decode_v3_two_byte_all_ones_is_unlimited() { + let bytes = [ + 0x02, 0, 0, 0, 0x03, 0, 0, 0, 0x01, 0x02, 0x01, 0, 0, 0, // + 0x00, 0x00, 0x01, 0x00, 0x01, 0x00, 0xFF, 0xFF, + ]; + let (raw, _) = SerializedSelection::decode(&bytes).unwrap(); + assert_eq!(raw.unlimited_dim(), Some(0)); + } + + #[test] + fn decode_irregular_block_count_beyond_buffer_is_error() { + let mut b = v1_blocks(1, &[(&[0], &[3])]); + b[20..24].copy_from_slice(&u32::MAX.to_le_bytes()); + assert!(Selection::decode_serialized(&b).is_err()); + } + + #[test] + fn decode_point_selection_is_refused() { + let bytes = [1u8, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]; + assert!(Selection::decode_serialized(&bytes).is_err()); + } + #[test] fn iter_linear_2d_block_row_major() { // A 2x2 block at the top-left of a 4x4 space => linear 0,1,4,5. diff --git a/crates/clawhdf5/tests/vds_interop.rs b/crates/clawhdf5/tests/vds_interop.rs new file mode 100644 index 0000000..43eaee4 --- /dev/null +++ b/crates/clawhdf5/tests/vds_interop.rs @@ -0,0 +1,156 @@ +//! Virtual Dataset (VDS) reads checked against libhdf5 (through h5py). +//! +//! Each test has h5py build virtual datasets and their source files in a temp +//! directory, record what libhdf5 reads back (shape and values) next to them, +//! and then compares that with what clawhdf5 reads from the same files. +//! +//! Skipped when python3 or h5py is unavailable, unless +//! `CLAWHDF5_REQUIRE_INTEROP=1`. + +use std::path::Path; +use std::process::Command; + +use clawhdf5::File; + +/// The Python interpreter to drive interop checks with (`CLAWHDF5_PYTHON` +/// lets these run against a virtualenv holding h5py). +fn python() -> String { + std::env::var("CLAWHDF5_PYTHON").unwrap_or_else(|_| "python3".to_string()) +} + +/// When `CLAWHDF5_REQUIRE_INTEROP=1` (set in CI), a missing Python dependency +/// is a test failure instead of a silent skip. +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, numpy"]) + .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; + } + }; +} + +/// Prelude for every generator script: `expect(file, dset, tag)` records what +/// libhdf5 reads for `file:dset` as `.expect` (shape line, values line). +const PRELUDE: &str = r#" +import h5py, numpy as np +def expect(fn, dset, tag): + with h5py.File(fn, "r") as f: + d = f[dset] + a = d[...] + with open(tag + ".expect", "w") as out: + out.write(" ".join(str(n) for n in d.shape) + "\n") + out.write(" ".join(repr(float(v)) for v in a.ravel()) + "\n") +"#; + +/// Run `body` (after [`PRELUDE`]) with `dir` as the working directory, so +/// relative source file names land next to the virtual file. +fn generate(dir: &Path, body: &str) { + let script = format!("{PRELUDE}\n{body}"); + let out = Command::new(python()) + .args(["-c", &script]) + .current_dir(dir) + .output() + .expect("failed to run python"); + assert!( + out.status.success(), + "generator failed:\nSTDOUT: {}\nSTDERR: {}", + String::from_utf8_lossy(&out.stdout), + String::from_utf8_lossy(&out.stderr) + ); +} + +/// What libhdf5 read for `tag`: (shape, values as f64). +fn expected(dir: &Path, tag: &str) -> (Vec, Vec) { + let text = std::fs::read_to_string(dir.join(format!("{tag}.expect"))).unwrap(); + let mut lines = text.lines(); + let parse_line = |l: Option<&str>| -> Vec { + l.unwrap_or("") + .split_whitespace() + .map(str::to_string) + .collect() + }; + let shape = parse_line(lines.next()) + .iter() + .map(|s| s.parse().unwrap()) + .collect(); + let values = parse_line(lines.next()) + .iter() + .map(|s| s.parse().unwrap()) + .collect(); + (shape, values) +} + +/// Assert clawhdf5 reads `file:dset` exactly as libhdf5 did for `tag`. +fn assert_matches_libhdf5(dir: &Path, file: &str, dset: &str, tag: &str) { + let (shape, values) = expected(dir, tag); + let f = File::open(dir.join(file)).unwrap(); + let ds = f.dataset(dset).unwrap(); + assert_eq!( + ds.shape().unwrap(), + shape, + "{tag}: shape differs from libhdf5" + ); + let got = ds + .read_f64() + .unwrap_or_else(|e| panic!("{tag}: read failed: {e}")); + assert_eq!(got.len(), values.len(), "{tag}: element count differs"); + for (i, (g, e)) in got.iter().zip(&values).enumerate() { + assert!( + g == e || (g.is_nan() && e.is_nan()), + "{tag}: element {i} is {g}, libhdf5 reads {e}\n ours: {got:?}\n libhdf5: {values:?}" + ); + } +} + +// --------------------------------------------------------------------------- +// Selection encodings +// --------------------------------------------------------------------------- + +/// Files written with the default (earliest) format bounds serialize every +/// VDS hyperslab as a version-1 *irregular* selection (4-byte block corners), +/// and a strided selection as many blocks. These were refused outright. +#[test] +fn vds_version1_irregular_hyperslab_selections() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + generate( + dir.path(), + r#" +with h5py.File("src.h5", "w") as s: + s.create_dataset("a", data=np.arange(12.0)) + s.create_dataset("m", data=np.arange(20.0).reshape(4, 5)) +with h5py.File("v1.h5", "w") as f: # default libver: hyperslab version 1 + f.create_dataset("local", data=np.arange(10.0) * -1) + lay = h5py.VirtualLayout(shape=(12,), dtype="f8") + lay[0:4] = h5py.VirtualSource(".", "local", shape=(10,))[2:6] + lay[4:10] = h5py.VirtualSource("src.h5", "a", shape=(12,))[::2] + lay[10:12] = h5py.VirtualSource("src.h5", "a", shape=(12,))[10:12] + f.create_virtual_dataset("strided", lay) + lay = h5py.VirtualLayout(shape=(4, 6), dtype="f8") + lay[:, 0:2] = h5py.VirtualSource("src.h5", "m", shape=(4, 5))[:, 3:5] + lay[:, 2:6:2] = h5py.VirtualSource("src.h5", "m", shape=(4, 5))[:, 0:2] + lay[:, 3:6:2] = h5py.VirtualSource("src.h5", "m", shape=(4, 5))[:, 2:4] + f.create_virtual_dataset("grid", lay) +expect("v1.h5", "strided", "strided") +expect("v1.h5", "grid", "grid") +"#, + ); + assert_matches_libhdf5(dir.path(), "v1.h5", "strided", "strided"); + assert_matches_libhdf5(dir.path(), "v1.h5", "grid", "grid"); +} diff --git a/docs/known-issues.md b/docs/known-issues.md index 6dfa369..c15e824 100644 --- a/docs/known-issues.md +++ b/docs/known-issues.md @@ -71,7 +71,8 @@ the VDS item, which is marked. - **Virtual datasets:** - **Wrong data:** unmapped regions read as 0 instead of the fill value. - `%b` printf-style source names are not expanded. - - Hyperslab selection versions 1 and 2 are refused. + - ~~Hyperslab selection versions 1 and 2 are refused.~~ Fixed 2026-09-25: + versions 1-3 and irregular hyperslabs are decoded. - **Files with a user block:** the base address is not applied. - **Old-style shared messages (version 1)** read the wrong address. - **Groups and links:** From 2c6c6c176e0e504ed7a91a11a98ee704ad4bb89f Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:59:15 -0500 Subject: [PATCH 2/4] fix(format): decode the version-1 VDS mapping list HDF5 2.0 writes With a 2.0 low version bound, libhdf5 stores the VDS mapping list as heap block version 1: every entry starts with a flags byte (0x04 same file, no file name; 0x01/0x02 file/dataset name shared with an earlier entry, whose index is stored in place of the name). The parser treated only a leading 0x04 byte as special, so a 0x00 flags byte read as an empty (same-file) name and shared names were read as garbage. Decode it as H5D__virtual_load_layout does, refusing unknown flags, forward references and block versions above 1. Test: vds_interop::vds_mapping_block_version1_shared_names (h5py libver=("v200","v200") with repeated long names; failed before with "unknown dataspace selection type") plus the exact heap block as a unit test. Co-Authored-By: Claude Opus 5.5 (1M context) --- CHANGELOG.md | 4 + crates/clawhdf5-format/src/data_layout.rs | 132 ++++++++++++++++++++-- crates/clawhdf5/tests/vds_interop.rs | 36 ++++++ docs/known-issues.md | 2 + 4 files changed, 162 insertions(+), 12 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index a349fba..fd8624f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -282,6 +282,10 @@ selection. Both are decoded now, as are irregular hyperslabs (a union of blocks, read in row-major order as libhdf5 iterates them). `SerializedSelection` exposes the raw form, including unlimited counts. + - The version-1 mapping list HDF5 2.0 writes (low version bound 2.0) was + misparsed: each entry's flags byte was read as the start of the source + file name, and names shared with an earlier entry (stored as that entry's + index) were not followed. Now decoded as `H5D__virtual_load_layout` does. - `clawhdf5-format` reader — **values returned wrong with no error:** - Fixed Array and Extensible Array chunk indexes were laid out by the dataset's current shape instead of its max shape (23 libhdf5 test files, diff --git a/crates/clawhdf5-format/src/data_layout.rs b/crates/clawhdf5-format/src/data_layout.rs index 8429c5d..b58843c 100644 --- a/crates/clawhdf5-format/src/data_layout.rs +++ b/crates/clawhdf5-format/src/data_layout.rs @@ -72,21 +72,33 @@ pub enum DataLayout { }, } +/// Version-1 VDS mapping flag: the source file name is stored by an earlier +/// entry, whose index follows in place of the name. +const VDS_SOURCE_FILE_SHARED: u8 = 0x01; +/// Version-1 VDS mapping flag: likewise for the source dataset name. +const VDS_SOURCE_DSET_SHARED: u8 = 0x02; +/// Version-1 VDS mapping flag: the source is in the virtual file itself +/// (`"."`); no file name is stored. +const VDS_SOURCE_SAME_FILE: u8 = 0x04; +const VDS_ALL_FLAGS: u8 = VDS_SOURCE_FILE_SHARED | VDS_SOURCE_DSET_SHARED | VDS_SOURCE_SAME_FILE; + /// Parse VDS mappings from global-heap object data. /// /// The global-heap block holding a VDS mapping list is laid out as -/// (reverse-engineered and validated against HDF5 2.0): +/// (`H5D__virtual_store_layout` / `H5D__virtual_load_layout` in libhdf5): /// /// ```text /// version(1) · nused(length_size, LE) · entry[nused] · checksum(4) /// ``` /// /// Each entry is: -/// - source file name — a null-terminated string in **block version 0**; in -/// **block version 1** a same-file reference is encoded as a single `0x04` -/// marker byte (the source file is the virtual file itself) in place of the -/// name; -/// - source dataset name (null-terminated string); +/// - **block version 1 only:** a flags byte. `0x04`: the source is in the +/// virtual file itself and no file name is stored; `0x01`/`0x02`: the +/// source file/dataset name is that of an earlier entry, whose index +/// (`length_size` bytes) is stored instead of the name. libhdf5 2.0 writes +/// version 1 when the file's low version bound is 2.0 and it saves space; +/// - source file name (null-terminated string, unless flagged above); +/// - source dataset name (null-terminated string, unless flagged above); /// - source selection (serialized `H5S` dataspace selection — self-describing /// in length); /// - virtual selection (serialized `H5S` dataspace selection). @@ -112,7 +124,7 @@ pub fn parse_vds_mappings( // `nused` is untrusted; don't pre-allocate from it. Each entry consumes at // least a few bytes, so the loop is naturally bounded by the heap data and // a bogus `nused` simply errors out on the first short read. - let mut mappings = Vec::new(); + let mut mappings: Vec = Vec::new(); // Reads one self-describing selection at `pos`, returning its raw bytes and // advancing past it — bounds-checked so a corrupt selection can't overrun. let read_selection = |heap_data: &[u8], pos: &mut usize| -> Result, FormatError> { @@ -132,17 +144,57 @@ pub fn parse_vds_mappings( Ok(bytes) }; - for _ in 0..nused { - // Source file name (with the version-1 same-file marker handled). - let source_file = if version >= 1 && heap_data.get(pos) == Some(&0x04) { + if version > 1 { + return Err(FormatError::ChunkedReadError( + "unsupported VDS mapping block version".into(), + )); + } + for i in 0..nused { + // Version 1 prefixes each entry with a flags byte; a name may then be + // omitted (same file) or replaced by the index of an earlier entry + // holding the same name (`H5D__virtual_load_layout`). + let flags = if version >= 1 { + let f = *heap_data.get(pos).ok_or(FormatError::UnexpectedEof { + expected: pos + 1, + available: heap_data.len(), + })?; pos += 1; + if f & !VDS_ALL_FLAGS != 0 { + return Err(FormatError::ChunkedReadError( + "unknown VDS mapping flags".into(), + )); + } + f + } else { + 0 + }; + // Index of an earlier entry, for a shared name. + let earlier = |pos: &mut usize| -> Result { + let idx = read_length(heap_data, *pos, length_size)?; + *pos += ls; + if idx >= i { + return Err(FormatError::ChunkedReadError( + "VDS mapping shares a name with a later entry".into(), + )); + } + Ok(idx as usize) + }; + + let source_file = if flags & VDS_SOURCE_SAME_FILE != 0 { String::from(".") + } else if flags & VDS_SOURCE_FILE_SHARED != 0 { + let idx = earlier(&mut pos)?; + mappings[idx].source_file.clone() } else { read_null_terminated_string(heap_data, &mut pos)? }; - // Source dataset name. - let source_dataset = read_null_terminated_string(heap_data, &mut pos)?; + let source_dataset = if flags & VDS_SOURCE_DSET_SHARED != 0 { + let idx = earlier(&mut pos)?; + mappings[idx].source_dataset.clone() + } else { + read_null_terminated_string(heap_data, &mut pos)? + }; // Source selection, then virtual selection (both self-describing length). let source_selection = read_selection(heap_data, &mut pos)?; @@ -849,6 +901,62 @@ mod tests { assert_eq!(v1.iter_linear_1d(8).unwrap(), vec![4, 5, 6, 7]); } + #[test] + fn parse_vds_mappings_v1_shared_names() { + // Written by HDF5 2.0 (h5py, libver=("v200", "v200")) for three + // mappings from `a_rather_long_source_file.h5:a_rather_long_dataset_name` + // and one from the same file: the entries carry flags 0x00, 0x03, 0x03 + // and 0x06, so names after the first are stored as entry indices. + let blob: &[u8] = &[ + 0x01, 0x04, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x61, 0x5f, 0x72, 0x61, + 0x74, 0x68, 0x65, 0x72, 0x5f, 0x6c, 0x6f, 0x6e, 0x67, 0x5f, 0x73, 0x6f, 0x75, 0x72, + 0x63, 0x65, 0x5f, 0x66, 0x69, 0x6c, 0x65, 0x2e, 0x68, 0x35, 0x00, 0x61, 0x5f, 0x72, + 0x61, 0x74, 0x68, 0x65, 0x72, 0x5f, 0x6c, 0x6f, 0x6e, 0x67, 0x5f, 0x64, 0x61, 0x74, + 0x61, 0x73, 0x65, 0x74, 0x5f, 0x6e, 0x61, 0x6d, 0x65, 0x00, 0x02, 0x00, 0x00, 0x00, + 0x03, 0x00, 0x00, 0x00, 0x01, 0x02, 0x01, 0x00, 0x00, 0x00, 0x00, 0x00, 0x01, 0x00, + 0x01, 0x00, 0x04, 0x00, 0x02, 0x00, 0x00, 0x00, 0x03, 0x00, 0x00, 0x00, 0x01, 0x02, + 0x02, 0x00, 0x00, 0x00, 0x00, 0x00, 0x01, 0x00, 0x01, 0x00, 0x01, 0x00, 0x00, 0x00, + 0x01, 0x00, 0x01, 0x00, 0x04, 0x00, 0x03, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, + 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x02, 0x00, 0x00, 0x00, 0x03, + 0x00, 0x00, 0x00, 0x01, 0x02, 0x01, 0x00, 0x00, 0x00, 0x04, 0x00, 0x01, 0x00, 0x01, + 0x00, 0x04, 0x00, 0x02, 0x00, 0x00, 0x00, 0x03, 0x00, 0x00, 0x00, 0x01, 0x02, 0x02, + 0x00, 0x00, 0x00, 0x01, 0x00, 0x01, 0x00, 0x01, 0x00, 0x01, 0x00, 0x00, 0x00, 0x01, + 0x00, 0x01, 0x00, 0x04, 0x00, 0x03, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, + 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x02, 0x00, 0x00, 0x00, 0x03, 0x00, + 0x00, 0x00, 0x01, 0x02, 0x01, 0x00, 0x00, 0x00, 0x08, 0x00, 0x01, 0x00, 0x01, 0x00, + 0x04, 0x00, 0x02, 0x00, 0x00, 0x00, 0x03, 0x00, 0x00, 0x00, 0x01, 0x02, 0x02, 0x00, + 0x00, 0x00, 0x02, 0x00, 0x01, 0x00, 0x01, 0x00, 0x01, 0x00, 0x00, 0x00, 0x01, 0x00, + 0x01, 0x00, 0x04, 0x00, 0x06, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x00, 0x02, + 0x00, 0x00, 0x00, 0x03, 0x00, 0x00, 0x00, 0x01, 0x02, 0x01, 0x00, 0x00, 0x00, 0x00, + 0x00, 0x01, 0x00, 0x01, 0x00, 0x04, 0x00, 0x02, 0x00, 0x00, 0x00, 0x03, 0x00, 0x00, + 0x00, 0x01, 0x02, 0x02, 0x00, 0x00, 0x00, 0x03, 0x00, 0x01, 0x00, 0x01, 0x00, 0x01, + 0x00, 0x00, 0x00, 0x01, 0x00, 0x01, 0x00, 0x04, 0x00, 0x8e, 0xa7, 0xea, 0x7a, + ]; + let mappings = parse_vds_mappings(blob, 8).unwrap(); + let names: Vec<(&str, &str)> = mappings + .iter() + .map(|m| (m.source_file.as_str(), m.source_dataset.as_str())) + .collect(); + let (file, dset) = ("a_rather_long_source_file.h5", "a_rather_long_dataset_name"); + assert_eq!( + names, + vec![(file, dset), (file, dset), (file, dset), (".", dset)] + ); + } + + #[test] + fn parse_vds_mappings_v1_forward_reference_is_error() { + // Entry 0 claiming to share entry 0's file name must not index past + // the entries decoded so far. + let mut blob = vec![0x01u8, 1, 0, 0, 0, 0, 0, 0, 0, 0x01]; + blob.extend_from_slice(&[0u8; 8]); + blob.extend_from_slice(b"d\0"); + assert!(parse_vds_mappings(&blob, 8).is_err()); + // Unknown flag bits are refused. + let blob = [0x01u8, 1, 0, 0, 0, 0, 0, 0, 0, 0x08, b'd', 0]; + assert!(parse_vds_mappings(&blob, 8).is_err()); + } + #[test] fn parse_vds_mappings_external_v0() { // Block version 0 with an explicit (external) source file name. diff --git a/crates/clawhdf5/tests/vds_interop.rs b/crates/clawhdf5/tests/vds_interop.rs index 43eaee4..f3ce7e3 100644 --- a/crates/clawhdf5/tests/vds_interop.rs +++ b/crates/clawhdf5/tests/vds_interop.rs @@ -154,3 +154,39 @@ expect("v1.h5", "grid", "grid") assert_matches_libhdf5(dir.path(), "v1.h5", "strided", "strided"); assert_matches_libhdf5(dir.path(), "v1.h5", "grid", "grid"); } + +// --------------------------------------------------------------------------- +// Mapping list encoding +// --------------------------------------------------------------------------- + +/// With a 2.0 low version bound libhdf5 writes the mapping list as block +/// version 1: a flags byte per entry, and repeated names stored as the index +/// of the entry that first spelled them out. The flags byte was mistaken for +/// an empty (same-file) name. +#[test] +fn vds_mapping_block_version1_shared_names() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + generate( + dir.path(), + r#" +name = "a_rather_long_dataset_name" +with h5py.File("a_rather_long_source_file.h5", "w") as s: + s.create_dataset(name, data=np.arange(12.0) + 100) +with h5py.File("shared.h5", "w", libver=("v200", "v200")) as f: + f.create_dataset(name, data=np.arange(12.0) * -1) + lay = h5py.VirtualLayout(shape=(4, 4), dtype="f8") + for i in range(3): + src = h5py.VirtualSource("a_rather_long_source_file.h5", name, shape=(12,)) + lay[i] = src[4 * i:4 * i + 4] + lay[3] = h5py.VirtualSource(".", name, shape=(12,))[0:4] + f.create_virtual_dataset("v", lay) +# the heap block must really be version 1 for this test to mean anything +raw = open("shared.h5", "rb").read() +gcol = raw.index(b"GCOL") +assert raw[gcol + 32] == 1, "expected a version-1 VDS mapping block" +expect("shared.h5", "v", "shared") +"#, + ); + assert_matches_libhdf5(dir.path(), "shared.h5", "v", "shared"); +} diff --git a/docs/known-issues.md b/docs/known-issues.md index c15e824..eb7d80f 100644 --- a/docs/known-issues.md +++ b/docs/known-issues.md @@ -73,6 +73,8 @@ the VDS item, which is marked. - `%b` printf-style source names are not expanded. - ~~Hyperslab selection versions 1 and 2 are refused.~~ Fixed 2026-09-25: versions 1-3 and irregular hyperslabs are decoded. + - ~~The version-1 mapping list written with a 2.0 low bound (flags byte, + shared names) was misparsed.~~ Found and fixed 2026-09-25. - **Files with a user block:** the base address is not applied. - **Old-style shared messages (version 1)** read the wrong address. - **Groups and links:** From e94a52a88b214f9fdaa93c2817b16fe61ebcbfe0 Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 22:03:47 -0500 Subject: [PATCH 3/4] fix(format): read unmapped VDS elements as the virtual dataset's fill value MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Elements of a virtual dataset that no mapping supplies (unmapped regions, a missing source file, a missing source dataset) read as 0 instead of the fill value libhdf5 returns — silent wrong data for any VDS created with a non-zero fillvalue (read-matrix cases 0471/0472: -1 and 7 read as 0). A missing source dataset was an error; libhdf5 reads it as fill. Move VDS assembly into a new vds module following H5Dvirtual.c: vds::read_virtual_dataset takes the dataset's fill value and a VdsFileResolver that can refuse a name, and reports how many elements were unmapped. Sources are read with their own fill value, and a source whose datatype differs from the virtual dataset's is an error (libhdf5 converts). File passes the dataset's fill value, resolves source names against the virtual file's directory, and refuses names that leave it with an error instead of reading them as fill. read_selection on a VDS goes through the same fill-aware path. The raw-read API (read_raw_data_full*) has no fill value, so it now errors for a VDS with unmapped elements instead of guessing zeros. Tests: vds_interop::vds_unmapped_regions_read_as_fill_value (external, same-file, missing file/dataset, sparse source with its own fill, int fill; earliest and latest format) and vds_source_outside_directory_is_an_error_not_fill, both against h5py; integration_test::v4_virtual_dataset_raw_api_refuses_to_guess_the_fill_value. Co-Authored-By: Claude Opus 5.5 (1M context) --- CHANGELOG.md | 15 + crates/clawhdf5-format/src/data_read.rs | 181 ++----- crates/clawhdf5-format/src/lib.rs | 1 + crates/clawhdf5-format/src/vds.rs | 497 ++++++++++++++++++ .../clawhdf5-format/tests/integration_test.rs | 75 ++- crates/clawhdf5/src/reader.rs | 76 ++- crates/clawhdf5/tests/vds_interop.rs | 99 ++++ docs/known-issues.md | 8 +- 8 files changed, 785 insertions(+), 167 deletions(-) create mode 100644 crates/clawhdf5-format/src/vds.rs diff --git a/CHANGELOG.md b/CHANGELOG.md index fd8624f..de54b91 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -275,6 +275,21 @@ ### Correctness - `clawhdf5-format` virtual datasets (VDS), checked against HDF5 2.0 through h5py (`crates/clawhdf5/tests/vds_interop.rs`): + - **Wrong data:** elements no mapping supplies — unmapped regions, and + mappings whose source file or dataset is missing — read as 0 instead of + the virtual dataset's fill value (e.g. h5py `fillvalue=-1`). Assembly moved + to the new `vds` module: `vds::read_virtual_dataset` takes the fill value + and a resolver that can refuse a name (`VdsFileResolver`), and `File` + passes the dataset's fill value. A missing source *dataset* read as an + error; it is fill now, as in libhdf5. Source datasets are read with their + own fill value for unallocated chunks, and a source whose datatype differs + from the virtual dataset's is an error (libhdf5 converts; we do not). + `File` now refuses a source name that leaves the virtual file's directory + (`../x.h5`, absolute paths), or any external source of a `File::from_bytes` + file, with an error — these used to read as fill. + **Behaviour change:** the raw-read API (`read_raw_data_full*`), which has + no fill value, now returns an error for a virtual dataset with unmapped + elements instead of zeros. - Hyperslab selection versions 1 and 2 were refused ("only version-3 hyperslab selections are supported"). Version 1 is what libhdf5 writes for every VDS created with the default format bounds (h5py's default), so diff --git a/crates/clawhdf5-format/src/data_read.rs b/crates/clawhdf5-format/src/data_read.rs index 82e6544..76c2b0e 100644 --- a/crates/clawhdf5-format/src/data_read.rs +++ b/crates/clawhdf5-format/src/data_read.rs @@ -191,14 +191,9 @@ fn read_raw_data_full_impl( offset_size, length_size, ), - DataLayout::Virtual { - global_heap_address, - global_heap_index, - .. - } => read_virtual_data( + DataLayout::Virtual { .. } => read_virtual_data( file_data, - *global_heap_address, - *global_heap_index, + layout, dataspace, datatype, offset_size, @@ -465,158 +460,54 @@ pub fn read_raw_data_selection( } } -/// Assemble a **Virtual Dataset (VDS)** from its source mappings. +/// Assemble a **Virtual Dataset (VDS)** through the raw-read API, which has no +/// access to the dataset's fill value message. /// -/// Supports virtual datasets of any rank. Same-file sources are read directly; -/// **external-file** sources are read through the caller-supplied `resolver`, -/// which maps a stored source file name to that file's bytes. Each mapping's -/// selected source elements are scattered into the virtual buffer at the -/// positions given by the virtual selection (both enumerated in row-major -/// order, as HDF5 pairs them). Unmapped regions are left at the zero fill value. -/// -/// A mapping whose external source file the resolver cannot supply (`None`) is -/// skipped, leaving its region at fill — matching HDF5's tolerance of missing -/// sources. An external source with no resolver at all is a hard error. +/// Delegates to [`crate::vds::read_virtual_dataset`]. Because the fill value +/// is unknown here, a virtual dataset with any element no mapping supplies +/// (an unmapped region, or a missing source file or dataset) is an error +/// rather than a guess at the fill value; so is one whose extent libhdf5 +/// would report differently from the stored dataspace (unlimited mappings). +/// Use [`crate::vds::read_virtual_dataset`] to read those. #[allow(clippy::too_many_arguments)] fn read_virtual_data( file_data: &[u8], - global_heap_address: Option, - global_heap_index: u32, + layout: &DataLayout, dataspace: &Dataspace, datatype: &Datatype, offset_size: u8, length_size: u8, resolver: Option<&VdsSourceResolver>, ) -> Result, FormatError> { - use crate::data_layout::parse_vds_mappings; - use crate::global_heap::GlobalHeapCollection; - use crate::selection::Selection; - - let elem_size = datatype.type_size() as usize; - let mut out = crate::chunked_read::alloc_output(crate::chunked_read::checked_byte_len( - dataspace.checked_num_elements()?, - elem_size, - )?)?; - - let virtual_dims = &dataspace.dimensions; - - let addr = global_heap_address.ok_or_else(|| { - FormatError::ChunkedReadError("virtual dataset has no mapping global heap".into()) - })?; - let coll = GlobalHeapCollection::parse(file_data, addr as usize, length_size)?; - let obj = - coll.get_object(global_heap_index as u16) - .ok_or(FormatError::GlobalHeapObjectNotFound { - collection_address: addr, - index: global_heap_index as u16, - })?; - let mappings = parse_vds_mappings(&obj.data, length_size)?; - - for m in &mappings { - let same_file = m.source_file.is_empty() || m.source_file == "."; - - // Resolve the bytes of the file holding this source dataset. - let external; - let src_file_data: &[u8] = if same_file { - file_data - } else { - let r = resolver.ok_or_else(|| { - FormatError::ChunkedReadError( - "external-file virtual dataset sources require a file resolver".into(), - ) - })?; - match r(&m.source_file) { - Some(bytes) => { - external = bytes; - &external - } - // Source file unavailable: leave this region at fill value. - None => continue, - } - }; - - let (vsel, _) = Selection::decode_serialized(&m.virtual_selection)?; - let (ssel, _) = Selection::decode_serialized(&m.source_selection)?; - - let (src_raw, src_dims) = - read_named_dataset_raw(src_file_data, &m.source_dataset, offset_size, length_size)?; - - let vidx = vsel.iter_linear(virtual_dims)?; - let sidx = ssel.iter_linear(&src_dims)?; - if vidx.len() != sidx.len() { - return Err(FormatError::ChunkedReadError( - "virtual/source selection element counts differ".into(), - )); - } - - for (&v, &s) in vidx.iter().zip(sidx.iter()) { - let (vo, so) = (v as usize * elem_size, s as usize * elem_size); - if vo + elem_size > out.len() || so + elem_size > src_raw.len() { - return Err(FormatError::ChunkedReadError( - "virtual dataset selection out of bounds".into(), - )); - } - out[vo..vo + elem_size].copy_from_slice(&src_raw[so..so + elem_size]); - } - } - - Ok(out) -} - -/// Read a named dataset's raw (decoded) bytes and its dimensions, navigating -/// from the superblock. Used to pull VDS source datasets out of the same file. -fn read_named_dataset_raw( - file_data: &[u8], - path: &str, - _offset_size: u8, - _length_size: u8, -) -> Result<(Vec, Vec), FormatError> { - use crate::filter_pipeline::FilterPipeline; - use crate::group_v2::resolve_path_any; - use crate::message_type::MessageType; - use crate::object_header::ObjectHeader; - use crate::signature::find_signature; - use crate::superblock::Superblock; - - let sig = find_signature(file_data)?; - let sb = Superblock::parse(file_data, sig)?; - let addr = resolve_path_any(file_data, &sb, path)?; - let hdr = ObjectHeader::parse(file_data, addr as usize, sb.offset_size, sb.length_size)?; - - let find = |t: MessageType| hdr.messages.iter().find(|m| m.msg_type == t); - let ds_msg = find(MessageType::Dataspace) - .ok_or_else(|| FormatError::ChunkedReadError("VDS source has no dataspace".into()))?; - let dataspace = Dataspace::parse(&ds_msg.data, sb.length_size)?; - let dt_msg = find(MessageType::Datatype) - .ok_or_else(|| FormatError::ChunkedReadError("VDS source has no datatype".into()))?; - let (datatype, _) = Datatype::parse(&dt_msg.data)?; - let dl_msg = find(MessageType::DataLayout) - .ok_or_else(|| FormatError::ChunkedReadError("VDS source has no data layout".into()))?; - let layout = DataLayout::parse(&dl_msg.data, sb.offset_size, sb.length_size)?; - // A virtual dataset whose source is itself another virtual dataset could - // form a cycle (A -> B -> A) and recurse into a stack overflow. Nested - // virtual sources are exotic and unsupported, so stop here cleanly. - if matches!(layout, DataLayout::Virtual { .. }) { + let wrapped = + resolver.map(|r| move |name: &str| -> Result>, FormatError> { Ok(r(name)) }); + let wrapped_ref = wrapped.as_ref().map(|w| w as &crate::vds::VdsFileResolver); + let v = crate::vds::read_virtual_dataset( + file_data, + layout, + dataspace, + datatype, + None, + offset_size, + length_size, + wrapped_ref, + )?; + if v.dims != dataspace.dimensions { return Err(FormatError::ChunkedReadError( - "virtual dataset source is itself virtual (unsupported)".into(), + "virtual dataset extent differs from its stored dataspace; \ + read it with vds::read_virtual_dataset" + .into(), )); } - let pipeline = find(MessageType::FilterPipeline) - .map(|m| FilterPipeline::parse(&m.data)) - .transpose()?; - - let raw = read_raw_data_full( - file_data, - &layout, - &dataspace, - &datatype, - pipeline.as_ref(), - sb.offset_size, - sb.length_size, - )?; - Ok((raw, dataspace.dimensions.clone())) + if v.unmapped > 0 { + return Err(FormatError::ChunkedReadError( + "virtual dataset has elements no source supplies, which read as its \ + fill value; read it with vds::read_virtual_dataset and the fill value" + .into(), + )); + } + Ok(v.data) } - /// Extract selected elements from a full dataset buffer. pub fn extract_selection_from_buffer( full_data: &[u8], diff --git a/crates/clawhdf5-format/src/lib.rs b/crates/clawhdf5-format/src/lib.rs index 4a6c2a1..879708d 100644 --- a/crates/clawhdf5-format/src/lib.rs +++ b/crates/clawhdf5-format/src/lib.rs @@ -100,6 +100,7 @@ pub mod signature; pub mod superblock; pub mod symbol_table; pub mod type_builders; +pub mod vds; pub mod vl_data; #[cfg(feature = "provenance")] diff --git a/crates/clawhdf5-format/src/vds.rs b/crates/clawhdf5-format/src/vds.rs new file mode 100644 index 0000000..3ac0232 --- /dev/null +++ b/crates/clawhdf5-format/src/vds.rs @@ -0,0 +1,497 @@ +//! Virtual Dataset (VDS) assembly, following libhdf5's `H5Dvirtual.c`. +//! +//! A virtual dataset stores no data of its own: a list of mappings (kept in +//! the global heap) pairs a selection of the virtual dataspace with a +//! selection of a *source* dataset, in the same file (`"."`) or another one. +//! Reading it means reading each source and scattering the selected source +//! elements into the virtual buffer, pairing the two selections element by +//! element in row-major order. Elements no mapping supplies — unmapped +//! regions, and mappings whose source file or dataset does not exist — read +//! as the virtual dataset's **fill value**, as in libhdf5. +//! +//! Source files other than the virtual file itself are obtained through a +//! caller-supplied [`VdsFileResolver`], since this crate has no filesystem. + +#[cfg(not(feature = "std"))] +use alloc::{format, string::String, vec, vec::Vec}; + +use crate::data_layout::{DataLayout, VdsMapping, parse_vds_mappings}; +use crate::dataspace::Dataspace; +use crate::datatype::Datatype; +use crate::error::FormatError; +use crate::selection::{SerializedSelection, UNLIMITED}; + +/// Resolves the name of an external VDS source file, as stored in the +/// mapping, to that file's bytes. +/// +/// `Ok(None)` means the file does not exist; its mappings then read as the +/// fill value, as libhdf5 does for a missing source. `Err` refuses the name +/// (e.g. a path the caller will not follow) and fails the read, so that a +/// refused source is never passed off as fill. +pub type VdsFileResolver<'a> = dyn Fn(&str) -> Result>, FormatError> + 'a; + +/// A fully assembled virtual dataset. +#[derive(Debug, Clone, PartialEq)] +pub struct VirtualData { + /// The virtual dataset's extent. + pub dims: Vec, + /// Raw element bytes, row-major, in the virtual dataset's datatype. + pub data: Vec, + /// Number of elements no mapping supplied; they hold the fill value. + pub unmapped: u64, +} + +fn vds_err(msg: impl Into) -> FormatError { + FormatError::ChunkedReadError(msg.into()) +} + +/// One mapping with its selections decoded. +struct Mapping { + file: String, + dataset: String, + vsel: SerializedSelection, + ssel: SerializedSelection, +} + +/// Load and decode the mapping list of a virtual layout. +fn load_mappings( + file_data: &[u8], + layout: &DataLayout, + length_size: u8, +) -> Result, FormatError> { + let DataLayout::Virtual { + global_heap_address, + global_heap_index, + .. + } = layout + else { + return Err(vds_err("not a virtual dataset layout")); + }; + let Some(addr) = *global_heap_address else { + return Ok(Vec::new()); + }; + let coll = + crate::global_heap::GlobalHeapCollection::parse(file_data, addr as usize, length_size)?; + let index = u16::try_from(*global_heap_index) + .map_err(|_| vds_err("VDS mapping heap index out of range"))?; + let obj = coll + .get_object(index) + .ok_or(FormatError::GlobalHeapObjectNotFound { + collection_address: addr, + index, + })?; + parse_vds_mappings(&obj.data, length_size)? + .into_iter() + .map(|m: VdsMapping| { + let (vsel, _) = SerializedSelection::decode(&m.virtual_selection)?; + let (ssel, _) = SerializedSelection::decode(&m.source_selection)?; + Ok(Mapping { + file: m.source_file, + dataset: m.source_dataset, + vsel, + ssel, + }) + }) + .collect() +} + +/// The virtual dataset's extent as libhdf5 reports it (`H5Dget_space`). +/// +/// For a virtual dataset whose mappings are all of fixed size this is the +/// stored dataspace. +pub fn virtual_dataset_extent( + file_data: &[u8], + layout: &DataLayout, + dataspace: &Dataspace, + _offset_size: u8, + length_size: u8, + _resolver: Option<&VdsFileResolver>, +) -> Result, FormatError> { + let mappings = load_mappings(file_data, layout, length_size)?; + check_fixed(&mappings)?; + Ok(dataspace.dimensions.clone()) +} + +fn check_fixed(mappings: &[Mapping]) -> Result<(), FormatError> { + if mappings + .iter() + .any(|m| m.vsel.unlimited_dim().is_some() || m.ssel.unlimited_dim().is_some()) + { + return Err(vds_err( + "unlimited virtual dataset mappings are not supported", + )); + } + Ok(()) +} + +/// Read a whole virtual dataset. +/// +/// `fill` is the virtual dataset's fill value (from its fill value message; +/// `None` for the default of zeros); every element no mapping supplies holds +/// it. External source files are read through `resolver`; without one, a +/// mapping to another file is an error. +#[allow(clippy::too_many_arguments)] +pub fn read_virtual_dataset( + file_data: &[u8], + layout: &DataLayout, + dataspace: &Dataspace, + datatype: &Datatype, + fill: Option<&[u8]>, + offset_size: u8, + length_size: u8, + resolver: Option<&VdsFileResolver>, +) -> Result { + let mappings = load_mappings(file_data, layout, length_size)?; + check_fixed(&mappings)?; + let dims = virtual_dataset_extent( + file_data, + layout, + dataspace, + offset_size, + length_size, + resolver, + )?; + + let elem_size = datatype.type_size() as usize; + let total = dims + .iter() + .try_fold(1u64, |acc, &d| acc.checked_mul(d)) + .ok_or_else(|| FormatError::Overflow("virtual dataset extent".into()))?; + let mut data = crate::chunked_read::alloc_output(crate::chunked_read::checked_byte_len( + total, elem_size, + )?)?; + if let Some(fill) = fill.filter(|f| f.len() == elem_size && f.iter().any(|&b| b != 0)) { + for element in data.chunks_exact_mut(elem_size) { + element.copy_from_slice(fill); + } + } + let mut mapped = vec![false; usize::try_from(total).map_err(|_| vds_err("VDS too large"))?]; + + let mut sources = Sources::new(file_data, resolver); + for m in &mappings { + let Some(src) = sources.dataset(&m.file, &m.dataset, datatype)? else { + continue; // missing source file or dataset: fill + }; + let vidx = selection_indices(&m.vsel, &dims, None)?; + let sidx = selection_indices(&m.ssel, &src.dims, None)?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + + let unmapped = mapped.iter().filter(|&&m| !m).count() as u64; + Ok(VirtualData { + dims, + data, + unmapped, + }) +} + +/// Copy source element `sidx[i]` to virtual element `vidx[i]` for every `i`. +fn scatter( + out: &mut [u8], + mapped: &mut [bool], + src: &[u8], + vidx: &[u64], + sidx: &[u64], + elem_size: usize, +) -> Result<(), FormatError> { + if vidx.len() != sidx.len() { + return Err(vds_err("virtual/source selection element counts differ")); + } + for (&v, &s) in vidx.iter().zip(sidx) { + let (vo, so) = (v as usize * elem_size, s as usize * elem_size); + if vo + elem_size > out.len() || so + elem_size > src.len() { + return Err(vds_err("virtual dataset selection out of bounds")); + } + out[vo..vo + elem_size].copy_from_slice(&src[so..so + elem_size]); + mapped[v as usize] = true; + } + Ok(()) +} + +/// Row-major linear indices of the elements `sel` selects in a dataspace of +/// shape `dims`, in the order libhdf5 iterates them (row-major). +/// +/// `clip` = `(dim, limit)` drops every coordinate `>= limit` in `dim` — the +/// clipping libhdf5 applies to an unlimited selection +/// (`H5S_hyper_clip_unlim`); an unlimited selection must be clipped. +fn selection_indices( + sel: &SerializedSelection, + dims: &[u64], + clip: Option<(usize, u64)>, +) -> Result, FormatError> { + let overflow = || FormatError::Overflow("VDS selection index overflow".into()); + let rank = dims.len(); + let total = dims + .iter() + .try_fold(1u64, |acc, &d| acc.checked_mul(d)) + .ok_or_else(overflow)?; + let mut row_stride = vec![1u64; rank]; + for d in (0..rank.saturating_sub(1)).rev() { + row_stride[d] = row_stride[d + 1] + .checked_mul(dims[d + 1]) + .ok_or_else(overflow)?; + } + if sel.rank().is_some_and(|r| r != rank) { + return Err(vds_err("VDS selection rank does not match dataspace rank")); + } + let limit = |d: usize| match clip { + Some((cd, l)) if cd == d => l, + _ => u64::MAX, + }; + + match sel { + SerializedSelection::All => Ok((0..total).collect()), + SerializedSelection::None => Ok(Vec::new()), + SerializedSelection::Regular { + start, + stride, + count, + block, + } => { + // Selected coordinates along each dimension, in order. + let mut per_dim: Vec> = Vec::with_capacity(rank); + for d in 0..rank { + let lim = limit(d); + if (count[d] == UNLIMITED || block[d] == UNLIMITED) && lim == u64::MAX { + return Err(vds_err("unlimited VDS selection was not clipped")); + } + let mut coords = Vec::new(); + let mut ci = 0u64; + 'blocks: while ci < count[d] { + let base = ci + .checked_mul(stride[d]) + .and_then(|o| start[d].checked_add(o)) + .ok_or_else(overflow)?; + if base >= lim { + break; + } + let mut bi = 0u64; + while bi < block[d] { + let coord = base.checked_add(bi).ok_or_else(overflow)?; + if coord >= lim { + break 'blocks; + } + // Past the extent is malformed; bail before the list + // can grow without bound. + if coord >= dims[d] { + return Err(vds_err("VDS selection exceeds the dataspace extent")); + } + coords.push(coord); + bi += 1; + } + ci += 1; + } + per_dim.push(coords); + } + if per_dim.iter().any(|c| c.is_empty()) { + return Ok(Vec::new()); + } + let n = per_dim + .iter() + .try_fold(1usize, |acc, c| acc.checked_mul(c.len())) + .ok_or_else(overflow)?; + let mut out = Vec::with_capacity(n); + let mut idx = vec![0usize; rank]; + loop { + let lin: u64 = (0..rank).map(|d| per_dim[d][idx[d]] * row_stride[d]).sum(); + out.push(lin); + // Mixed-radix increment, last dimension fastest. + let mut d = rank; + loop { + if d == 0 { + return Ok(out); + } + d -= 1; + idx[d] += 1; + if idx[d] < per_dim[d].len() { + break; + } + idx[d] = 0; + } + } + } + SerializedSelection::Blocks { + rank: _, + starts, + ends, + } => { + // libhdf5 serializes the union as disjoint blocks, so their volumes + // never add up to more than the dataspace. + let mut volume = 0u64; + for (s, e) in starts.chunks_exact(rank).zip(ends.chunks_exact(rank)) { + for d in 0..rank { + if e[d] >= dims[d] { + return Err(vds_err("VDS selection exceeds the dataspace extent")); + } + } + let v = s + .iter() + .zip(e) + .try_fold(1u64, |acc, (&s, &e)| acc.checked_mul(e - s + 1)) + .ok_or_else(overflow)?; + volume = volume.checked_add(v).ok_or_else(overflow)?; + if volume > total { + return Err(vds_err("VDS selection blocks overlap")); + } + } + let mut out = Vec::with_capacity(volume as usize); + for (s, e) in starts.chunks_exact(rank).zip(ends.chunks_exact(rank)) { + let mut cur = s.to_vec(); + 'block: loop { + if (0..rank).all(|d| cur[d] < limit(d)) { + out.push((0..rank).map(|d| cur[d] * row_stride[d]).sum()); + } + for d in (0..rank).rev() { + if cur[d] < e[d] { + cur[d] += 1; + continue 'block; + } + cur[d] = s[d]; + } + break; + } + } + out.sort_unstable(); + out.dedup(); + Ok(out) + } + } +} + +/// A source dataset's decoded contents. +struct SourceData { + dims: Vec, + raw: Vec, +} + +/// Source files and datasets, fetched on demand. The most recently used +/// external file is kept, since consecutive mappings usually share one. +struct Sources<'a, 'r> { + file_data: &'a [u8], + resolver: Option<&'r VdsFileResolver<'r>>, + cached_file: Option<(String, Option>)>, +} + +impl<'a, 'r> Sources<'a, 'r> { + fn new(file_data: &'a [u8], resolver: Option<&'r VdsFileResolver<'r>>) -> Self { + Sources { + file_data, + resolver, + cached_file: None, + } + } + + /// The bytes of source file `name`, or `None` if it does not exist. + fn file(&mut self, name: &str) -> Result, FormatError> { + if name == "." { + return Ok(Some(self.file_data)); + } + if self.cached_file.as_ref().is_none_or(|(n, _)| n != name) { + let resolver = self.resolver.ok_or_else(|| { + vds_err("external-file virtual dataset sources require a file resolver") + })?; + self.cached_file = Some((String::from(name), resolver(name)?)); + } + Ok(self.cached_file.as_ref().and_then(|(_, b)| b.as_deref())) + } + + /// Read source dataset `path` from file `file`, or `None` when either + /// does not exist. Its datatype must be the virtual dataset's: libhdf5 + /// converts between types here, which is not supported. + fn dataset( + &mut self, + file: &str, + path: &str, + datatype: &Datatype, + ) -> Result, FormatError> { + let Some(bytes) = self.file(file)? else { + return Ok(None); + }; + read_source(bytes, path, datatype) + } +} + +/// Read source dataset `path` of the file in `file_data` in full (its own +/// fill value applied to unallocated chunks), or `None` if it does not exist. +fn read_source( + file_data: &[u8], + path: &str, + datatype: &Datatype, +) -> Result, FormatError> { + use crate::filter_pipeline::FilterPipeline; + use crate::message_type::MessageType; + use crate::object_header::ObjectHeader; + use crate::shared_message::message_data_with_sohm; + + let sig = crate::signature::find_signature(file_data)?; + let sb = crate::superblock::Superblock::parse(file_data, sig)?; + let (os, ls) = (sb.offset_size, sb.length_size); + let addr = match crate::group_v2::resolve_path_any(file_data, &sb, path) { + Ok(a) => a, + Err(FormatError::PathNotFound(_)) => return Ok(None), + Err(e) => return Err(e), + }; + let hdr = ObjectHeader::parse(file_data, addr as usize, os, ls)?; + let msg = |t: MessageType| { + hdr.messages + .iter() + .find(|m| m.msg_type == t) + .ok_or_else(|| vds_err(format!("VDS source {path} has no {t:?} message"))) + }; + let dataspace = Dataspace::parse( + &message_data_with_sohm(file_data, msg(MessageType::Dataspace)?, os, ls)?, + ls, + )?; + let (src_type, _) = Datatype::parse(&message_data_with_sohm( + file_data, + msg(MessageType::Datatype)?, + os, + ls, + )?)?; + if &src_type != datatype { + return Err(vds_err(format!( + "VDS source {path} has a different datatype from the virtual dataset \ + (type conversion is not supported)" + ))); + } + let layout = DataLayout::parse(&msg(MessageType::DataLayout)?.data, os, ls)?; + // A source that is itself virtual could form a cycle (A -> B -> A) and + // recurse without bound. Nested virtual sources are not supported. + if matches!(layout, DataLayout::Virtual { .. }) { + return Err(vds_err( + "virtual dataset source is itself virtual (unsupported)", + )); + } + let pipeline = hdr + .messages + .iter() + .find(|m| m.msg_type == MessageType::FilterPipeline) + .map(|m| { + message_data_with_sohm(file_data, m, os, ls).and_then(|d| FilterPipeline::parse(&d)) + }) + .transpose()?; + let raw = crate::fill_value::read_full_with_fill( + &hdr.messages, + file_data, + &layout, + &dataspace, + src_type.type_size() as usize, + os, + ls, + || { + crate::data_read::read_raw_data_full( + file_data, + &layout, + &dataspace, + &src_type, + pipeline.as_ref(), + os, + ls, + ) + }, + )?; + Ok(Some(SourceData { + dims: dataspace.dimensions, + raw, + })) +} diff --git a/crates/clawhdf5-format/tests/integration_test.rs b/crates/clawhdf5-format/tests/integration_test.rs index 7aab9f0..a43ecda 100644 --- a/crates/clawhdf5-format/tests/integration_test.rs +++ b/crates/clawhdf5-format/tests/integration_test.rs @@ -83,6 +83,45 @@ fn read_chunked_dataset(file_data: &[u8], dataset_path: &str) -> (Vec, Datat (raw, datatype, dataspace) } +/// Helper: read a virtual dataset with `vds::read_virtual_dataset`, giving it +/// the dataset's own fill value (same-file sources only). +fn read_virtual_fixture(file_data: &[u8], path: &str) -> (Vec, Datatype) { + let sig = find_signature(file_data).unwrap(); + let sb = Superblock::parse(file_data, sig).unwrap(); + let addr = resolve_path_any(file_data, &sb, path).unwrap(); + let hdr = + ObjectHeader::parse(file_data, addr as usize, sb.offset_size, sb.length_size).unwrap(); + let msg = |t: MessageType| hdr.messages.iter().find(|m| m.msg_type == t).unwrap(); + let ds = Dataspace::parse(&msg(MessageType::Dataspace).data, sb.length_size).unwrap(); + let (dt, _) = Datatype::parse(&msg(MessageType::Datatype).data).unwrap(); + let layout = DataLayout::parse( + &msg(MessageType::DataLayout).data, + sb.offset_size, + sb.length_size, + ) + .unwrap(); + let fill = clawhdf5_format::fill_value::dataset_fill_value_in( + file_data, + &hdr.messages, + sb.offset_size, + sb.length_size, + ) + .unwrap(); + let v = clawhdf5_format::vds::read_virtual_dataset( + file_data, + &layout, + &ds, + &dt, + fill.as_deref(), + sb.offset_size, + sb.length_size, + None, + ) + .unwrap(); + assert_eq!(v.dims, ds.dimensions); + (v.data, dt) +} + /// Helper: read any dataset (contiguous or chunked) as f64. fn read_dataset_f64_any(bytes: &[u8], path: &str) -> Vec { let sig = find_signature(bytes).unwrap(); @@ -672,7 +711,7 @@ fn v4_virtual_dataset_same_file_read() { // virt[4:8] <- (unmapped) => fill 0 // virt[8:12] <- src_b[0:4] (ALL) => 20,21,22,23 let file_data = include_bytes!("fixtures/vds_same_file.h5"); - let (raw, datatype, _) = read_chunked_dataset(file_data, "virt"); + let (raw, datatype) = read_virtual_fixture(file_data, "virt"); let values = read_as_i32(&raw, &datatype).unwrap(); assert_eq!( values, @@ -681,6 +720,38 @@ fn v4_virtual_dataset_same_file_read() { ); } +#[test] +fn v4_virtual_dataset_raw_api_refuses_to_guess_the_fill_value() { + // The raw read API has no fill value message, so a virtual dataset with an + // unmapped region is an error there instead of zeros that may be wrong. + let file_data = include_bytes!("fixtures/vds_same_file.h5"); + let sig = find_signature(file_data).unwrap(); + let sb = Superblock::parse(file_data, sig).unwrap(); + let addr = resolve_path_any(file_data, &sb, "virt").unwrap(); + let hdr = + ObjectHeader::parse(file_data, addr as usize, sb.offset_size, sb.length_size).unwrap(); + let msg = |t: MessageType| hdr.messages.iter().find(|m| m.msg_type == t).unwrap(); + let ds = Dataspace::parse(&msg(MessageType::Dataspace).data, sb.length_size).unwrap(); + let (dt, _) = Datatype::parse(&msg(MessageType::Datatype).data).unwrap(); + let layout = DataLayout::parse( + &msg(MessageType::DataLayout).data, + sb.offset_size, + sb.length_size, + ) + .unwrap(); + let err = read_raw_data_full( + file_data, + &layout, + &ds, + &dt, + None, + sb.offset_size, + sb.length_size, + ) + .unwrap_err(); + assert!(err.to_string().contains("fill value"), "{err}"); +} + #[test] fn v4_virtual_dataset_2d_same_file_read() { // A 4x4 virtual dataset assembled from two 2x2 same-file sources placed as @@ -689,7 +760,7 @@ fn v4_virtual_dataset_2d_same_file_read() { // virt[2:4,2:4] <- src_b = [[5,6],[7,8]] // everything else -> fill 0 let file_data = include_bytes!("fixtures/vds_2d_same_file.h5"); - let (raw, datatype, _) = read_chunked_dataset(file_data, "virt"); + let (raw, datatype) = read_virtual_fixture(file_data, "virt"); let values = read_as_i32(&raw, &datatype).unwrap(); assert_eq!( values, diff --git a/crates/clawhdf5/src/reader.rs b/crates/clawhdf5/src/reader.rs index bb6a8fd..2e7bd3d 100644 --- a/crates/clawhdf5/src/reader.rs +++ b/crates/clawhdf5/src/reader.rs @@ -485,6 +485,7 @@ impl<'f> Dataset<'f> { self.file.length_size(), )?; let fill_matters = !clawhdf5_format::fill_value::has_storage(&dl) + || matches!(dl, DataLayout::Virtual { .. }) || (matches!(dl, DataLayout::Chunked { .. }) && !clawhdf5_format::fill_value::is_default(fill.as_deref())); if fill_matters { @@ -837,24 +838,9 @@ impl<'f> Dataset<'f> { let pipeline = self.filter_pipeline()?; // Virtual datasets are assembled from source datasets; the per-file - // chunk cache does not apply. Route them through the resolver path so - // external sibling files resolve relative to this file's directory. + // chunk cache does not apply. if matches!(dl, DataLayout::Virtual { .. }) { - let base_dir = self.file.base_dir.clone(); - let resolver = move |name: &str| -> Option> { - let dir = base_dir.as_ref()?; - std::fs::read(dir.join(sibling_file_name(name)?)).ok() - }; - return Ok(data_read::read_raw_data_full_with_resolver( - self.file.data.as_bytes(), - &dl, - &ds, - &dt, - pipeline.as_ref(), - self.file.offset_size(), - self.file.length_size(), - Some(&resolver), - )?); + return self.read_virtual(&dl, &ds, &dt); } // Unallocated storage reads as the dataset's fill value. @@ -880,6 +866,62 @@ impl<'f> Dataset<'f> { }, ) } + + /// Resolver for external Virtual Dataset source files: names are + /// resolved against the directory of the file that holds the virtual + /// dataset, as libhdf5 does. A missing file is `Ok(None)` (its mappings + /// read as the fill value); a name that would leave that directory is + /// refused with an error rather than read as fill. + fn vds_resolver(&self) -> impl Fn(&str) -> Result>, FormatError> + use<> { + let base_dir = self.file.base_dir.clone(); + move |name: &str| { + let Some(dir) = base_dir.as_ref() else { + return Err(FormatError::ChunkedReadError(format!( + "virtual dataset source file {name:?} cannot be resolved for an in-memory file" + ))); + }; + let rel = sibling_file_name(name).ok_or_else(|| { + FormatError::ChunkedReadError(format!( + "virtual dataset source file {name:?} is outside the virtual file's \ + directory and is not followed" + )) + })?; + match std::fs::read(dir.join(rel)) { + Ok(bytes) => Ok(Some(bytes)), + Err(e) if e.kind() == std::io::ErrorKind::NotFound => Ok(None), + Err(e) => Err(FormatError::ChunkedReadError(format!( + "cannot read virtual dataset source file {name:?}: {e}" + ))), + } + } + } + + /// Read a whole virtual dataset; unmapped elements hold its fill value. + fn read_virtual( + &self, + dl: &DataLayout, + ds: &Dataspace, + dt: &Datatype, + ) -> Result, Error> { + let fill = clawhdf5_format::fill_value::dataset_fill_value_in( + self.file.data.as_bytes(), + &self.header.messages, + self.file.offset_size(), + self.file.length_size(), + )?; + let resolver = self.vds_resolver(); + let v = clawhdf5_format::vds::read_virtual_dataset( + self.file.data.as_bytes(), + dl, + ds, + dt, + fill.as_deref(), + self.file.offset_size(), + self.file.length_size(), + Some(&resolver), + )?; + Ok(v.data) + } } // --------------------------------------------------------------------------- diff --git a/crates/clawhdf5/tests/vds_interop.rs b/crates/clawhdf5/tests/vds_interop.rs index f3ce7e3..032d77c 100644 --- a/crates/clawhdf5/tests/vds_interop.rs +++ b/crates/clawhdf5/tests/vds_interop.rs @@ -190,3 +190,102 @@ expect("shared.h5", "v", "shared") ); assert_matches_libhdf5(dir.path(), "shared.h5", "v", "shared"); } + +// --------------------------------------------------------------------------- +// Fill value +// --------------------------------------------------------------------------- + +/// Elements no mapping supplies read as the virtual dataset's fill value, not +/// as 0: unmapped regions, a missing source file, a missing source dataset. +/// A source's own unallocated chunks read as *its* fill value. +#[test] +fn vds_unmapped_regions_read_as_fill_value() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + generate( + dir.path(), + r#" +for i in range(3): + with h5py.File(f"src_{i}.h5", "w") as s: + s.create_dataset("data", data=np.arange(10.0) + i * 100) +with h5py.File("sparse_src.h5", "w") as s: + d = s.create_dataset("data", shape=(10,), chunks=(5,), dtype="f8", fillvalue=42.0) + d[0:5] = np.arange(5.0) + 1000 # the second chunk is never written +for libver in ["earliest", "latest"]: + with h5py.File(f"fill_{libver}.h5", "w", libver=libver) as f: + f.create_dataset("local", data=np.arange(10.0) * -1) + lay = h5py.VirtualLayout(shape=(6, 10), dtype="f8") + for i in range(3): + lay[i] = h5py.VirtualSource(f"src_{i}.h5", "data", shape=(10,)) + lay[3] = h5py.VirtualSource("no_such_file.h5", "data", shape=(10,)) + lay[4] = h5py.VirtualSource("src_0.h5", "no_such_dataset", shape=(10,)) + # row 5 is not mapped at all + f.create_virtual_dataset("files", lay, fillvalue=-1.0) + lay = h5py.VirtualLayout(shape=(20,), dtype="f8") + lay[0:10] = h5py.VirtualSource(".", "local", shape=(10,)) + f.create_virtual_dataset("same_file", lay, fillvalue=7.0) + lay = h5py.VirtualLayout(shape=(12,), dtype="f8") + lay[1:11] = h5py.VirtualSource("sparse_src.h5", "data", shape=(10,)) + f.create_virtual_dataset("sparse_source", lay, fillvalue=-3.5) + lay = h5py.VirtualLayout(shape=(3, 4), dtype="i4") + lay[1, :] = h5py.VirtualSource(".", "ints", shape=(4,)) + f.create_dataset("ints", data=np.arange(4, dtype="i4") + 1) + f.create_virtual_dataset("int_fill", lay, fillvalue=-99) + for name in ["files", "same_file", "sparse_source", "int_fill"]: + expect(f"fill_{libver}.h5", name, f"{name}_{libver}") +"#, + ); + for libver in ["earliest", "latest"] { + let file = format!("fill_{libver}.h5"); + for name in ["files", "same_file", "sparse_source", "int_fill"] { + assert_matches_libhdf5(dir.path(), &file, name, &format!("{name}_{libver}")); + } + } + + // A selection read goes through the same fill-aware assembly. + let f = File::open(dir.path().join("fill_latest.h5")).unwrap(); + let sel = clawhdf5::Selection::slice(std::slice::from_ref(&(8..14))); + let got = f + .dataset("same_file") + .unwrap() + .read_f64_selection(&sel) + .unwrap(); + assert_eq!(got, vec![-8.0, -9.0, 7.0, 7.0, 7.0, 7.0]); +} + +/// A source name that would leave the virtual file's directory is refused +/// with an error; it used to be skipped and read silently as fill. +#[test] +fn vds_source_outside_directory_is_an_error_not_fill() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + std::fs::create_dir(dir.path().join("sub")).unwrap(); + generate( + dir.path(), + r#" +with h5py.File("src.h5", "w") as s: + s.create_dataset("data", data=np.arange(4.0)) +with h5py.File("sub/up.h5", "w", libver="latest") as f: + lay = h5py.VirtualLayout(shape=(4,), dtype="f8") + lay[:] = h5py.VirtualSource("../src.h5", "data", shape=(4,)) + f.create_virtual_dataset("v", lay, fillvalue=-1.0) +with h5py.File("nested.h5", "w", libver="latest") as f: + lay = h5py.VirtualLayout(shape=(4,), dtype="f8") + lay[:] = h5py.VirtualSource("sub/inner.h5", "data", shape=(4,)) + f.create_virtual_dataset("v", lay, fillvalue=-1.0) +with h5py.File("sub/inner.h5", "w") as s: + s.create_dataset("data", data=np.arange(4.0) + 10) +expect("nested.h5", "v", "nested") +"#, + ); + // libhdf5 resolves "../src.h5" (and would read [0, 1, 2, 3]); we refuse + // to leave the directory, and say so. + let f = File::open(dir.path().join("sub/up.h5")).unwrap(); + let err = f.dataset("v").unwrap().read_f64().unwrap_err(); + assert!( + err.to_string().contains("not followed"), + "unexpected error: {err}" + ); + // A relative name below the virtual file's directory resolves there. + assert_matches_libhdf5(dir.path(), "nested.h5", "v", "nested"); +} diff --git a/docs/known-issues.md b/docs/known-issues.md index eb7d80f..b1de913 100644 --- a/docs/known-issues.md +++ b/docs/known-issues.md @@ -63,13 +63,15 @@ segfault or abort. ## Gaps found by the 2026-09-25 HDF5 audit (open) -**Status:** open. These fail with an error; none returns wrong data, except -the VDS item, which is marked. +**Status:** open. These fail with an error; none returns wrong data (the VDS +fill-value item that did is fixed). - **Layout message versions 1 and 2** (HDF5 1.6-era files): 84 of the 686 sweep files, `InvalidLayoutVersion`. This is the largest single gap. - **Virtual datasets:** - - **Wrong data:** unmapped regions read as 0 instead of the fill value. + - ~~**Wrong data:** unmapped regions read as 0 instead of the fill value.~~ + Fixed 2026-09-25: unmapped elements and missing sources read as the + virtual dataset's fill value. - `%b` printf-style source names are not expanded. - ~~Hyperslab selection versions 1 and 2 are refused.~~ Fixed 2026-09-25: versions 1-3 and irregular hyperslabs are decoded. From b4a44a2e66adb86b87c45d8e1337a7d7f37fe546 Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 22:10:09 -0500 Subject: [PATCH 4/4] feat(format): read unlimited and printf-style VDS mappings like libhdf5 Unlimited VDS mappings were refused, and printf-style source names ("f-%b.h5") were not expanded, so those regions read as fill (read-matrix case 0470: 29 of 30 values wrong). All 7 virtual datasets in the libhdf5 test set use such mappings. Implement H5Dvirtual.c's semantics in the vds module: - %b is the block number, %% a literal %, other specifiers are an error; block j of the virtual selection comes from the source named with j, probing from 0 to the first missing source (printf gap 0); - unlimited source/virtual selections are clipped to what the source's current extent fills (H5S_hyper_get_clip_extent_match, partial last block included); - the extent is recomputed as H5Dget_space does (view "last available": the largest clip, never below what limited mappings need), exposed as vds::virtual_dataset_extent and used by Dataset::shape(); - a source in the other byte order is byte-swapped; other conversions stay an error. Tests: vds_interop::vds_printf_source_names, vds_unlimited_mappings_follow_source_extents (h5py low-level API, earliest and latest format) and vds_libhdf5_test_files (vds-eiger, 4_vds and vds-percival-unlim-maxmin from HDF5's tools/test/testfiles/vds, committed as fixtures) all compare shape and values with h5py; unit tests for the clip arithmetic, name parsing and mapping rules. Co-Authored-By: Claude Opus 5.5 (1M context) --- CHANGELOG.md | 13 + crates/clawhdf5-format/src/vds.rs | 762 ++++++++++++++++-- crates/clawhdf5/src/reader.rs | 17 +- crates/clawhdf5/tests/fixtures/vds/4_0.h5 | Bin 0 -> 4581 bytes crates/clawhdf5/tests/fixtures/vds/4_1.h5 | Bin 0 -> 4581 bytes crates/clawhdf5/tests/fixtures/vds/4_2.h5 | Bin 0 -> 4581 bytes crates/clawhdf5/tests/fixtures/vds/4_vds.h5 | Bin 0 -> 5496 bytes crates/clawhdf5/tests/fixtures/vds/README.md | 14 + crates/clawhdf5/tests/fixtures/vds/a.h5 | Bin 0 -> 7736 bytes crates/clawhdf5/tests/fixtures/vds/b.h5 | Bin 0 -> 7736 bytes crates/clawhdf5/tests/fixtures/vds/c.h5 | Bin 0 -> 7736 bytes crates/clawhdf5/tests/fixtures/vds/d.h5 | Bin 0 -> 7736 bytes crates/clawhdf5/tests/fixtures/vds/f-0.h5 | Bin 0 -> 4144 bytes crates/clawhdf5/tests/fixtures/vds/f-3.h5 | Bin 0 -> 4144 bytes .../clawhdf5/tests/fixtures/vds/vds-eiger.h5 | Bin 0 -> 5496 bytes .../fixtures/vds/vds-percival-unlim-maxmin.h5 | Bin 0 -> 5496 bytes crates/clawhdf5/tests/vds_interop.rs | 159 ++++ docs/known-issues.md | 8 +- 18 files changed, 889 insertions(+), 84 deletions(-) create mode 100644 crates/clawhdf5/tests/fixtures/vds/4_0.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/4_1.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/4_2.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/4_vds.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/README.md create mode 100644 crates/clawhdf5/tests/fixtures/vds/a.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/b.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/c.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/d.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/f-0.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/f-3.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/vds-eiger.h5 create mode 100644 crates/clawhdf5/tests/fixtures/vds/vds-percival-unlim-maxmin.h5 diff --git a/CHANGELOG.md b/CHANGELOG.md index de54b91..aabbce6 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -290,6 +290,19 @@ **Behaviour change:** the raw-read API (`read_raw_data_full*`), which has no fill value, now returns an error for a virtual dataset with unmapped elements instead of zeros. + - Unlimited and printf-style mappings are supported (all 7 VDS files in the + libhdf5 test set are such mappings, e.g. Eiger/Percival detector layouts). + `%b` in a source file or dataset name is the block number and `%%` a + literal `%` (other `%` sequences are an error, as in libhdf5); block `j` + is read from the source named with `j`, probing from 0 up to the first + missing source. Unlimited source/virtual selections cover as much as the + source's current extent fills, including a partial last block. As + libhdf5 does on `H5Dget_space`, the extent is recomputed from the sources + present (default "last available" view, printf gap 0) — + `vds::virtual_dataset_extent`, used by `Dataset::shape()` — so e.g. + `vds-eiger.h5` is `[5, 10, 10]`, not its stored `[20, 10, 10]`. A source + stored in the other byte order is byte-swapped (libhdf5 converts); + other type conversions remain an error. - Hyperslab selection versions 1 and 2 were refused ("only version-3 hyperslab selections are supported"). Version 1 is what libhdf5 writes for every VDS created with the default format bounds (h5py's default), so diff --git a/crates/clawhdf5-format/src/vds.rs b/crates/clawhdf5-format/src/vds.rs index 3ac0232..099c14a 100644 --- a/crates/clawhdf5-format/src/vds.rs +++ b/crates/clawhdf5-format/src/vds.rs @@ -45,12 +45,149 @@ fn vds_err(msg: impl Into) -> FormatError { FormatError::ChunkedReadError(msg.into()) } +/// Upper bound on the printf-style source datasets probed for one mapping. +const MAX_PRINTF_BLOCKS: u64 = 1 << 20; + +/// A source file or dataset name, parsed for printf-style `%b` block-number +/// substitutions (`H5D_virtual_parse_source_name`): `%b` is the block +/// number, `%%` a literal `%`, and any other `%` sequence is invalid. +#[derive(Debug, Clone, PartialEq)] +struct SourceName { + /// Literal text around the substitutions: `segments.len() == subs + 1`. + segments: Vec, +} + +impl SourceName { + fn parse(name: &str) -> Result { + let mut segments = vec![String::new()]; + let mut chars = name.chars(); + while let Some(c) = chars.next() { + if c != '%' { + segments.last_mut().expect("never empty").push(c); + continue; + } + match chars.next() { + Some('b') => segments.push(String::new()), + Some('%') => segments.last_mut().expect("never empty").push('%'), + _ => { + return Err(vds_err(format!( + "invalid format specifier in VDS source name {name:?}" + ))); + } + } + } + Ok(SourceName { segments }) + } + + /// Number of `%b` substitutions. + fn subs(&self) -> usize { + self.segments.len() - 1 + } + + /// The name with every `%b` replaced by `block`. + fn build(&self, block: u64) -> String { + let mut out = String::new(); + for (i, seg) in self.segments.iter().enumerate() { + if i > 0 { + out.push_str(&format!("{block}")); + } + out.push_str(seg); + } + out + } +} + +/// How a mapping's selections relate, as libhdf5 classifies them. +#[derive(Debug, Clone, Copy, PartialEq)] +enum Kind { + /// Both selections have a fixed size. + Fixed, + /// Both are unlimited (in `vdim` / `sdim`): the mapping grows with the + /// source dataset's extent. + Unlimited { vdim: usize, sdim: usize }, + /// The virtual selection repeats a block without limit in `vdim`; block + /// `j` comes from the source named by substituting `j` for `%b`. + Printf { vdim: usize }, +} + /// One mapping with its selections decoded. struct Mapping { - file: String, - dataset: String, + file: SourceName, + dataset: SourceName, vsel: SerializedSelection, ssel: SerializedSelection, + kind: Kind, +} + +impl Mapping { + fn new(m: VdsMapping) -> Result { + let (vsel, _) = SerializedSelection::decode(&m.virtual_selection)?; + let (ssel, _) = SerializedSelection::decode(&m.source_selection)?; + let file = SourceName::parse(&m.source_file)?; + let dataset = SourceName::parse(&m.source_dataset)?; + let subs = file.subs() + dataset.subs(); + // The checks of H5D_virtual_check_mapping_pre/_post. + let kind = match (vsel.unlimited_dim(), ssel.unlimited_dim()) { + (Some(vdim), None) => { + if subs == 0 { + return Err(vds_err( + "unlimited virtual selection with a limited source selection \ + and no %b in the source names", + )); + } + match &vsel { + SerializedSelection::Regular { count, block, .. } + if count[vdim] == UNLIMITED && block[vdim] != UNLIMITED => {} + _ => { + return Err(vds_err( + "printf VDS mapping needs a virtual selection with an unlimited count", + )); + } + } + Kind::Printf { vdim } + } + (Some(vdim), Some(sdim)) => { + if non_unlimited_elements(&vsel, vdim) != non_unlimited_elements(&ssel, sdim) { + return Err(vds_err( + "unlimited VDS mapping: virtual and source selections differ \ + outside the unlimited dimension", + )); + } + Kind::Unlimited { vdim, sdim } + } + (None, Some(_)) => { + return Err(vds_err( + "VDS mapping with an unlimited source selection and a limited \ + virtual selection is not supported", + )); + } + (None, None) => Kind::Fixed, + }; + if subs > 0 && !matches!(kind, Kind::Printf { .. }) { + return Err(vds_err( + "%b in a VDS source name without an unlimited virtual selection", + )); + } + Ok(Mapping { + file, + dataset, + vsel, + ssel, + kind, + }) + } +} + +/// Elements a regular selection selects outside dimension `skip`. +fn non_unlimited_elements(sel: &SerializedSelection, skip: usize) -> Option { + let SerializedSelection::Regular { count, block, .. } = sel else { + return None; + }; + (0..count.len()) + .filter(|&d| d != skip) + .try_fold(1u64, |acc, d| { + acc.checked_mul(count[d].checked_mul(block[d])?) + }) } /// Load and decode the mapping list of a virtual layout. @@ -82,49 +219,269 @@ fn load_mappings( })?; parse_vds_mappings(&obj.data, length_size)? .into_iter() - .map(|m: VdsMapping| { - let (vsel, _) = SerializedSelection::decode(&m.virtual_selection)?; - let (ssel, _) = SerializedSelection::decode(&m.source_selection)?; - Ok(Mapping { - file: m.source_file, - dataset: m.source_dataset, - vsel, - ssel, - }) - }) + .map(Mapping::new) .collect() } +/// `H5S__hyper_get_clip_diminfo`: the count and block a regular selection has +/// in its unlimited dimension once clipped to `clip`. +fn clip_diminfo(start: u64, stride: u64, count: u64, block: u64, clip: u64) -> (u64, u64) { + if start >= clip { + if block == UNLIMITED { + (count, 0) + } else { + (0, block) + } + } else if block == UNLIMITED || block == stride { + (1, clip - start) + } else { + ((clip - start).div_ceil(stride.max(1)), block) + } +} + +/// The unlimited-dimension parameters (start, stride, count, block) of a +/// regular selection. +fn unlim_diminfo(sel: &SerializedSelection, dim: usize) -> Result<[u64; 4], FormatError> { + match sel { + SerializedSelection::Regular { + start, + stride, + count, + block, + } => Ok([start[dim], stride[dim], count[dim], block[dim]]), + _ => Err(vds_err( + "unlimited VDS selection is not a regular hyperslab", + )), + } +} + +/// `H5S_hyper_get_clip_extent_match` with `incl_trail = false` (the +/// "last available" view): the extent to clip `clip_sel` (unlimited in +/// `clip_dim`) to so that it holds as many slices as `match_sel` (unlimited +/// in `match_dim`) holds when clipped to `match_clip`. +fn clip_extent_match( + clip_sel: &SerializedSelection, + clip_dim: usize, + match_sel: &SerializedSelection, + match_dim: usize, + match_clip: u64, +) -> Result { + let overflow = || FormatError::Overflow("VDS clip extent overflow".into()); + let [mstart, mstride, mcount, mblock] = unlim_diminfo(match_sel, match_dim)?; + let (count, block) = clip_diminfo(mstart, mstride, mcount, mblock, match_clip); + let slices = if block == 0 || count == 0 { + 0 + } else if count == 1 { + block + } else { + let mut n = block.checked_mul(count).ok_or_else(overflow)?; + let span = mstride + .checked_mul(count - 1) + .and_then(|s| s.checked_add(block)) + .ok_or_else(overflow)?; + let room = match_clip - mstart; + if span > room { + n -= span - room; + } + n + }; + + // H5S__hyper_get_clip_extent_real + let [start, stride, _, block] = unlim_diminfo(clip_sel, clip_dim)?; + if slices == 0 { + return Ok(0); + } + let extent = if block == UNLIMITED || block == stride { + start.checked_add(slices) + } else { + let full = slices / block; + let rem = slices - full * block; + if rem > 0 { + full.checked_mul(stride) + .and_then(|o| start.checked_add(o)) + .and_then(|e| e.checked_add(rem)) + } else { + (full - 1) + .checked_mul(stride) + .and_then(|o| start.checked_add(o)) + .and_then(|e| e.checked_add(block)) + } + }; + extent.ok_or_else(overflow) +} + +/// How each mapping is read, and the resulting extent. +struct Plan { + dims: Vec, + steps: Vec, +} + +#[derive(Clone, Copy)] +enum Step { + Fixed, + /// Read with the virtual selection clipped to `vclip` and the source + /// selection to the source's extent; `None` if the source is missing. + Unlimited(Option), + /// Read blocks `0..blocks`. + Printf(u64), +} + +/// Work out the extent libhdf5 gives the virtual dataset +/// (`H5D__virtual_set_extent_unlim`, default view `H5D_VDS_LAST_AVAILABLE` +/// with a printf gap of 0) and how much of each unlimited mapping is read. +fn plan(mappings: &[Mapping], stored: &[u64], sources: &mut Sources) -> Result { + let overflow = || FormatError::Overflow("VDS extent overflow".into()); + let rank = stored.len(); + let mut new_dims: Vec> = vec![None; rank]; + // Minimum extent needed by the limited parts of every virtual selection + // (H5D_virtual_update_min_dims). + let mut min_dims = vec![0u64; rank]; + let mut steps = Vec::with_capacity(mappings.len()); + + for m in mappings { + let skip = match m.kind { + Kind::Unlimited { vdim, .. } | Kind::Printf { vdim } => Some(vdim), + Kind::Fixed => None, + }; + if let Some(ends) = selection_bounds_end(&m.vsel)? { + if ends.len() != rank { + return Err(vds_err("VDS selection rank does not match dataspace rank")); + } + for d in (0..rank).filter(|&d| Some(d) != skip) { + min_dims[d] = min_dims[d].max(ends[d].checked_add(1).ok_or_else(overflow)?); + } + } + + let (vdim, clip, step) = match m.kind { + Kind::Fixed => { + steps.push(Step::Fixed); + continue; + } + Kind::Unlimited { vdim, sdim } => { + match sources.dims(&m.file.build(0), &m.dataset.build(0))? { + Some(src_dims) => { + let extent = *src_dims.get(sdim).ok_or_else(|| { + vds_err("VDS source rank does not match its selection") + })?; + let clip = clip_extent_match(&m.vsel, vdim, &m.ssel, sdim, extent)?; + (vdim, clip, Step::Unlimited(Some(clip))) + } + None => (vdim, 0, Step::Unlimited(None)), + } + } + Kind::Printf { vdim } => { + // With a gap of 0 the search stops at the first missing + // source dataset. + let mut found = 0u64; + while sources + .dims(&m.file.build(found), &m.dataset.build(found))? + .is_some() + { + found += 1; + if found > MAX_PRINTF_BLOCKS { + return Err(vds_err("too many printf-style VDS source datasets")); + } + } + let clip = if found == 0 { + 0 + } else { + // End of block `found - 1` in the unlimited dimension. + let [start, stride, _, block] = unlim_diminfo(&m.vsel, vdim)?; + (found - 1) + .checked_mul(stride) + .and_then(|o| start.checked_add(o)) + .and_then(|e| e.checked_add(block)) + .ok_or_else(overflow)? + }; + (vdim, clip, Step::Printf(found)) + } + }; + if vdim >= rank { + return Err(vds_err("VDS selection rank does not match dataspace rank")); + } + new_dims[vdim] = Some(new_dims[vdim].map_or(clip, |n| n.max(clip))); + steps.push(step); + } + + let dims = (0..rank) + .map(|d| match new_dims[d] { + None => stored[d], + Some(n) => n.max(min_dims[d]), + }) + .collect(); + Ok(Plan { dims, steps }) +} + +/// The last selected coordinate in each dimension (`H5S_SELECT_BOUNDS`), +/// ignoring any unlimited dimension; `None` for ALL/NONE and empty selections. +fn selection_bounds_end(sel: &SerializedSelection) -> Result>, FormatError> { + let overflow = || FormatError::Overflow("VDS selection bounds overflow".into()); + match sel { + SerializedSelection::All | SerializedSelection::None => Ok(None), + SerializedSelection::Regular { + start, + stride, + count, + block, + } => { + if count.contains(&0) || block.contains(&0) { + return Ok(None); + } + let mut ends = Vec::with_capacity(start.len()); + for d in 0..start.len() { + if count[d] == UNLIMITED || block[d] == UNLIMITED { + ends.push(0); + continue; + } + let end = (count[d] - 1) + .checked_mul(stride[d]) + .and_then(|o| start[d].checked_add(o)) + .and_then(|e| e.checked_add(block[d] - 1)) + .ok_or_else(overflow)?; + ends.push(end); + } + Ok(Some(ends)) + } + SerializedSelection::Blocks { rank, ends, .. } => { + if ends.is_empty() { + return Ok(None); + } + let mut max = vec![0u64; *rank]; + for e in ends.chunks_exact(*rank) { + for d in 0..*rank { + max[d] = max[d].max(e[d]); + } + } + Ok(Some(max)) + } + } +} + /// The virtual dataset's extent as libhdf5 reports it (`H5Dget_space`). /// /// For a virtual dataset whose mappings are all of fixed size this is the -/// stored dataspace. +/// stored dataspace. With unlimited or printf-style mappings libhdf5 +/// recomputes the unlimited dimension from the sources present (the default +/// "last available" view: the largest extent any mapping can fill), which +/// needs the source files, read through `resolver`. pub fn virtual_dataset_extent( file_data: &[u8], layout: &DataLayout, dataspace: &Dataspace, _offset_size: u8, length_size: u8, - _resolver: Option<&VdsFileResolver>, + resolver: Option<&VdsFileResolver>, ) -> Result, FormatError> { let mappings = load_mappings(file_data, layout, length_size)?; - check_fixed(&mappings)?; - Ok(dataspace.dimensions.clone()) -} - -fn check_fixed(mappings: &[Mapping]) -> Result<(), FormatError> { - if mappings - .iter() - .any(|m| m.vsel.unlimited_dim().is_some() || m.ssel.unlimited_dim().is_some()) - { - return Err(vds_err( - "unlimited virtual dataset mappings are not supported", - )); + if mappings.iter().all(|m| m.kind == Kind::Fixed) { + return Ok(dataspace.dimensions.clone()); } - Ok(()) + let mut sources = Sources::new(file_data, resolver); + Ok(plan(&mappings, &dataspace.dimensions, &mut sources)?.dims) } -/// Read a whole virtual dataset. +/// Read a whole virtual dataset, at the extent +/// [`virtual_dataset_extent`] reports. /// /// `fill` is the virtual dataset's fill value (from its fill value message; /// `None` for the default of zeros); every element no mapping supplies holds @@ -137,20 +494,13 @@ pub fn read_virtual_dataset( dataspace: &Dataspace, datatype: &Datatype, fill: Option<&[u8]>, - offset_size: u8, + _offset_size: u8, length_size: u8, resolver: Option<&VdsFileResolver>, ) -> Result { let mappings = load_mappings(file_data, layout, length_size)?; - check_fixed(&mappings)?; - let dims = virtual_dataset_extent( - file_data, - layout, - dataspace, - offset_size, - length_size, - resolver, - )?; + let mut sources = Sources::new(file_data, resolver); + let Plan { dims, steps } = plan(&mappings, &dataspace.dimensions, &mut sources)?; let elem_size = datatype.type_size() as usize; let total = dims @@ -167,14 +517,46 @@ pub fn read_virtual_dataset( } let mut mapped = vec![false; usize::try_from(total).map_err(|_| vds_err("VDS too large"))?]; - let mut sources = Sources::new(file_data, resolver); - for m in &mappings { - let Some(src) = sources.dataset(&m.file, &m.dataset, datatype)? else { - continue; // missing source file or dataset: fill - }; - let vidx = selection_indices(&m.vsel, &dims, None)?; - let sidx = selection_indices(&m.ssel, &src.dims, None)?; - scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + for (m, step) in mappings.iter().zip(&steps) { + match (*step, m.kind) { + (Step::Fixed, _) => { + let (file, dset) = (m.file.build(0), m.dataset.build(0)); + let Some(src) = sources.dataset(&file, &dset, datatype)? else { + continue; // missing source file or dataset: fill + }; + let vidx = selection_indices(&m.vsel, &dims, None)?; + let sidx = selection_indices(&m.ssel, &src.dims, None)?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + (Step::Unlimited(Some(vclip)), Kind::Unlimited { vdim, sdim }) => { + let (file, dset) = (m.file.build(0), m.dataset.build(0)); + let Some(src) = sources.dataset(&file, &dset, datatype)? else { + continue; + }; + let vidx = selection_indices(&m.vsel, &dims, Some((vdim, vclip)))?; + let extent = *src + .dims + .get(sdim) + .ok_or_else(|| vds_err("VDS source rank does not match its selection"))?; + let sidx = selection_indices(&m.ssel, &src.dims, Some((sdim, extent)))?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + (Step::Unlimited(None), _) => {} + (Step::Printf(blocks), Kind::Printf { vdim }) => { + for j in 0..blocks { + let Some(src) = + sources.dataset(&m.file.build(j), &m.dataset.build(j), datatype)? + else { + continue; + }; + let vblock = unlim_block(&m.vsel, vdim, j)?; + let vidx = selection_indices(&vblock, &dims, None)?; + let sidx = selection_indices(&m.ssel, &src.dims, None)?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + } + _ => return Err(vds_err("internal error: VDS plan does not match mapping")), + } } let unmapped = mapped.iter().filter(|&&m| !m).count() as u64; @@ -185,6 +567,37 @@ pub fn read_virtual_dataset( }) } +/// `H5S_hyper_get_unlim_block`: block `j` of a selection whose count is +/// unlimited in `dim`. +fn unlim_block( + sel: &SerializedSelection, + dim: usize, + j: u64, +) -> Result { + let SerializedSelection::Regular { + start, + stride, + count, + block, + } = sel + else { + return Err(vds_err("printf VDS selection is not a regular hyperslab")); + }; + let mut start = start.clone(); + let mut count = count.clone(); + start[dim] = j + .checked_mul(stride[dim]) + .and_then(|o| start[dim].checked_add(o)) + .ok_or_else(|| FormatError::Overflow("VDS block start overflow".into()))?; + count[dim] = 1; + Ok(SerializedSelection::Regular { + start, + stride: stride.clone(), + count, + block: block.clone(), + }) +} + /// Copy source element `sidx[i]` to virtual element `vidx[i]` for every `i`. fn scatter( out: &mut [u8], @@ -395,6 +808,15 @@ impl<'a, 'r> Sources<'a, 'r> { Ok(self.cached_file.as_ref().and_then(|(_, b)| b.as_deref())) } + /// The extent of source dataset `path` in file `file`, or `None` when + /// either does not exist. + fn dims(&mut self, file: &str, path: &str) -> Result>, FormatError> { + let Some(bytes) = self.file(file)? else { + return Ok(None); + }; + Ok(open_source(bytes, path)?.map(|s| s.dataspace.dimensions)) + } + /// Read source dataset `path` from file `file`, or `None` when either /// does not exist. Its datatype must be the virtual dataset's: libhdf5 /// converts between types here, which is not supported. @@ -407,20 +829,37 @@ impl<'a, 'r> Sources<'a, 'r> { let Some(bytes) = self.file(file)? else { return Ok(None); }; - read_source(bytes, path, datatype) + let Some(src) = open_source(bytes, path)? else { + return Ok(None); + }; + read_source(bytes, src, path, datatype).map(Some) } } -/// Read source dataset `path` of the file in `file_data` in full (its own -/// fill value applied to unallocated chunks), or `None` if it does not exist. -fn read_source( - file_data: &[u8], +/// An opened source dataset's object header. +struct OpenSource { + offset_size: u8, + length_size: u8, + header: crate::object_header::ObjectHeader, + dataspace: Dataspace, +} + +fn source_message<'h>( + src: &'h OpenSource, path: &str, - datatype: &Datatype, -) -> Result, FormatError> { - use crate::filter_pipeline::FilterPipeline; + t: crate::message_type::MessageType, +) -> Result<&'h crate::object_header::HeaderMessage, FormatError> { + src.header + .messages + .iter() + .find(|m| m.msg_type == t) + .ok_or_else(|| vds_err(format!("VDS source {path} has no {t:?} message"))) +} + +/// Open source dataset `path` of the file in `file_data`, or `None` if there +/// is no such object (libhdf5 reads a missing source as fill). +fn open_source(file_data: &[u8], path: &str) -> Result, FormatError> { use crate::message_type::MessageType; - use crate::object_header::ObjectHeader; use crate::shared_message::message_data_with_sohm; let sig = crate::signature::find_signature(file_data)?; @@ -431,30 +870,55 @@ fn read_source( Err(FormatError::PathNotFound(_)) => return Ok(None), Err(e) => return Err(e), }; - let hdr = ObjectHeader::parse(file_data, addr as usize, os, ls)?; - let msg = |t: MessageType| { - hdr.messages - .iter() - .find(|m| m.msg_type == t) - .ok_or_else(|| vds_err(format!("VDS source {path} has no {t:?} message"))) + let header = crate::object_header::ObjectHeader::parse(file_data, addr as usize, os, ls)?; + let mut src = OpenSource { + offset_size: os, + length_size: ls, + header, + dataspace: Dataspace { + space_type: crate::dataspace::DataspaceType::Null, + rank: 0, + dimensions: Vec::new(), + max_dimensions: None, + }, }; - let dataspace = Dataspace::parse( - &message_data_with_sohm(file_data, msg(MessageType::Dataspace)?, os, ls)?, - ls, - )?; - let (src_type, _) = Datatype::parse(&message_data_with_sohm( - file_data, - msg(MessageType::Datatype)?, - os, - ls, - )?)?; - if &src_type != datatype { + let ds_msg = source_message(&src, path, MessageType::Dataspace)?; + src.dataspace = Dataspace::parse(&message_data_with_sohm(file_data, ds_msg, os, ls)?, ls)?; + Ok(Some(src)) +} + +/// Read an opened source dataset in full (its own fill value applied to +/// unallocated chunks). +fn read_source( + file_data: &[u8], + src: OpenSource, + path: &str, + datatype: &Datatype, +) -> Result { + use crate::filter_pipeline::FilterPipeline; + use crate::message_type::MessageType; + use crate::shared_message::message_data_with_sohm; + + let (os, ls) = (src.offset_size, src.length_size); + let dt_msg = source_message(&src, path, MessageType::Datatype)?; + let (src_type, _) = Datatype::parse(&message_data_with_sohm(file_data, dt_msg, os, ls)?)?; + // libhdf5 converts each source to the virtual dataset's type. Only the + // conversion that is a pure byte swap is done here. + let swap = if &src_type == datatype { + false + } else if differs_only_in_byte_order(&src_type, datatype) { + true + } else { return Err(vds_err(format!( "VDS source {path} has a different datatype from the virtual dataset \ - (type conversion is not supported)" + (only a byte-order conversion is supported)" ))); - } - let layout = DataLayout::parse(&msg(MessageType::DataLayout)?.data, os, ls)?; + }; + let layout = DataLayout::parse( + &source_message(&src, path, MessageType::DataLayout)?.data, + os, + ls, + )?; // A source that is itself virtual could form a cycle (A -> B -> A) and // recurse without bound. Nested virtual sources are not supported. if matches!(layout, DataLayout::Virtual { .. }) { @@ -462,7 +926,8 @@ fn read_source( "virtual dataset source is itself virtual (unsupported)", )); } - let pipeline = hdr + let pipeline = src + .header .messages .iter() .find(|m| m.msg_type == MessageType::FilterPipeline) @@ -471,10 +936,10 @@ fn read_source( }) .transpose()?; let raw = crate::fill_value::read_full_with_fill( - &hdr.messages, + &src.header.messages, file_data, &layout, - &dataspace, + &src.dataspace, src_type.type_size() as usize, os, ls, @@ -482,7 +947,7 @@ fn read_source( crate::data_read::read_raw_data_full( file_data, &layout, - &dataspace, + &src.dataspace, &src_type, pipeline.as_ref(), os, @@ -490,8 +955,141 @@ fn read_source( ) }, )?; - Ok(Some(SourceData { - dims: dataspace.dimensions, + let mut raw = raw; + if swap { + let size = datatype.type_size() as usize; + for element in raw.chunks_exact_mut(size) { + element.reverse(); + } + } + Ok(SourceData { + dims: src.dataspace.dimensions, raw, - })) + }) +} + +/// Whether two numeric types are identical apart from their byte order (both +/// little- or big-endian), so converting one to the other is a byte swap. +fn differs_only_in_byte_order(a: &Datatype, b: &Datatype) -> bool { + use crate::datatype::DatatypeByteOrder::{BigEndian, LittleEndian}; + let mut a = a.clone(); + match &mut a { + Datatype::FixedPoint { byte_order, .. } + | Datatype::FloatingPoint { byte_order, .. } + | Datatype::BitField { byte_order, .. } => { + *byte_order = match byte_order { + LittleEndian => BigEndian, + BigEndian => LittleEndian, + _ => return false, + } + } + _ => return false, + } + &a == b +} + +#[cfg(test)] +mod tests { + use super::*; + + fn regular(start: u64, stride: u64, count: u64, block: u64) -> SerializedSelection { + SerializedSelection::Regular { + start: vec![start], + stride: vec![stride], + count: vec![count], + block: vec![block], + } + } + + #[test] + fn printf_names_follow_libhdf5() { + let n = SourceName::parse("f-%b.h5").unwrap(); + assert_eq!(n.subs(), 1); + assert_eq!(n.build(12), "f-12.h5"); + let n = SourceName::parse("100%%_%b_%b").unwrap(); + assert_eq!(n.build(3), "100%_3_3"); + let n = SourceName::parse("plain%%name").unwrap(); + assert_eq!((n.subs(), n.build(7)), (0, "plain%name".to_string())); + // Anything else after '%' (or a trailing '%') is invalid. + assert!(SourceName::parse("a%d").is_err()); + assert!(SourceName::parse("a%").is_err()); + } + + #[test] + fn clip_extent_matches_libhdf5_arithmetic() { + // 7 source slices (contiguous unlimited source) into blocks of 3 + // every 4: two full blocks and one slice of a third -> extent 9. + let src = regular(0, 1, UNLIMITED, 1); + let v = regular(0, 4, UNLIMITED, 3); + assert_eq!(clip_extent_match(&v, 0, &src, 0, 7).unwrap(), 9); + // Exactly two blocks: the extent ends at the end of the last block. + assert_eq!(clip_extent_match(&v, 0, &src, 0, 6).unwrap(), 7); + // An empty source gives an empty mapping. + assert_eq!(clip_extent_match(&v, 0, &src, 0, 0).unwrap(), 0); + // Unlimited block: the extent is start + slices. + let vb = regular(2, 1, 1, UNLIMITED); + assert_eq!(clip_extent_match(&vb, 0, &src, 0, 5).unwrap(), 7); + // A strided source clipped mid-block counts only the selected slices. + let src2 = regular(1, 4, UNLIMITED, 2); // 1,2, 5,6, 9,10 ... + let dense = regular(0, 1, UNLIMITED, 1); + assert_eq!(clip_extent_match(&dense, 0, &src2, 0, 6).unwrap(), 3); + } + + #[test] + fn clipped_selection_drops_the_partial_tail() { + let v = regular(0, 4, UNLIMITED, 3); + assert_eq!( + selection_indices(&v, &[20], Some((0, 9))).unwrap(), + vec![0, 1, 2, 4, 5, 6, 8] + ); + // Unclipped unlimited selections cannot be enumerated. + assert!(selection_indices(&v, &[20], None).is_err()); + } + + #[test] + fn unlimited_mapping_rules() { + let sel = |s: &SerializedSelection| -> Vec { + // Serialize as version 2 (8-byte regular). + let SerializedSelection::Regular { + start, + stride, + count, + block, + } = s + else { + unreachable!() + }; + let mut b = Vec::new(); + b.extend_from_slice(&2u32.to_le_bytes()); + b.extend_from_slice(&2u32.to_le_bytes()); + b.push(1); + b.extend_from_slice(&0u32.to_le_bytes()); + b.extend_from_slice(&(start.len() as u32).to_le_bytes()); + for d in 0..start.len() { + for v in [start[d], stride[d], count[d], block[d]] { + b.extend_from_slice(&v.to_le_bytes()); + } + } + b + }; + let mapping = |file: &str, v: &SerializedSelection, s: &SerializedSelection| VdsMapping { + source_file: file.into(), + source_dataset: "d".into(), + source_selection: sel(s), + virtual_selection: sel(v), + }; + let unlim = regular(0, 10, UNLIMITED, 10); + let fixed = regular(0, 1, 1, 10); + // Unlimited virtual + limited source needs %b in a name ... + assert!(Mapping::new(mapping("f.h5", &unlim, &fixed)).is_err()); + let m = Mapping::new(mapping("f%b.h5", &unlim, &fixed)).unwrap(); + assert_eq!(m.kind, Kind::Printf { vdim: 0 }); + // ... and %b is refused anywhere else. + assert!(Mapping::new(mapping("f%b.h5", &fixed, &fixed)).is_err()); + let src_unlim = regular(0, 1, UNLIMITED, 1); + let m = Mapping::new(mapping("f.h5", &unlim, &src_unlim)).unwrap(); + assert_eq!(m.kind, Kind::Unlimited { vdim: 0, sdim: 0 }); + // An unlimited source into a limited virtual selection is refused. + assert!(Mapping::new(mapping("f.h5", &fixed, &src_unlim)).is_err()); + } } diff --git a/crates/clawhdf5/src/reader.rs b/crates/clawhdf5/src/reader.rs index 2e7bd3d..20e2c33 100644 --- a/crates/clawhdf5/src/reader.rs +++ b/crates/clawhdf5/src/reader.rs @@ -809,7 +809,22 @@ impl<'f> Dataset<'f> { fn dataspace(&self) -> Result { let data = self.required_payload(MessageType::Dataspace)?; - Ok(Dataspace::parse(&data, self.file.length_size())?) + let mut ds = Dataspace::parse(&data, self.file.length_size())?; + // libhdf5 reports a virtual dataset with unlimited or printf-style + // mappings at the extent its sources currently fill, not the stored + // one (`H5Dget_space`). + if let Ok(dl @ DataLayout::Virtual { .. }) = self.data_layout() { + let resolver = self.vds_resolver(); + ds.dimensions = clawhdf5_format::vds::virtual_dataset_extent( + self.file.data.as_bytes(), + &dl, + &ds, + self.file.offset_size(), + self.file.length_size(), + Some(&resolver), + )?; + } + Ok(ds) } fn data_layout(&self) -> Result { diff --git a/crates/clawhdf5/tests/fixtures/vds/4_0.h5 b/crates/clawhdf5/tests/fixtures/vds/4_0.h5 new file mode 100644 index 0000000000000000000000000000000000000000..5e71d20ca6a499ef7a2d9dc5ed9b4f7023ef4e27 GIT binary patch literal 4581 zcmeHKy-EW?5T3og;1Llog7{P2BUoB0W@X7X3^N7TkonUY^TZ-a%<&msZT|0@L=>GC=qc!y6d_4LGbCSr6B zkB|2A-Bu^>I!NkrNeI961|nC#exM3n^m(n%Y&jp7C-#rGpYhnGim{MAH*p&(*cX&B zuz?|*>Yo-ue!pPjkT{;G(* z+uxn);L@l!-mOf16pe3dPqZFX~&VVi5cax(w4XoYreJp$~hXKa$^s&G%w9POX zqi5#Jg8a->JY+@hQ5{tQRX`O`1yli5Kow90Q~_1sFBIs#yzOP#>m~k$m_C6%hHVV5 SFr2x78w{Hm_A#8hhlei-8D*^i literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/4_1.h5 b/crates/clawhdf5/tests/fixtures/vds/4_1.h5 new file mode 100644 index 0000000000000000000000000000000000000000..edad46e9fa9c7f16d74a51abaeda99560661c468 GIT binary patch literal 4581 zcmeHKy-EW?5T3oe;1Urp+KAQFB4}x8ASbE(7wo(vIV(Xh5z^Y}V~7Yof}M8p2`sEE z^flC(-5E^oFkoTj%#od)Z@-<(ewkUwRlBt@Q(vqDnY9o>?4~k&vb{GvBC)JzJTdGq z+1~|zVt^ii1;!hkzs%d;O;bRogW(F8{0nN%U#~rxd4_fQ4Y`(gCJ`M4FL-1Dt*#1So}CK+I)^n9K|{51KrI`3C}Uiwl5_5&$M- zh*|}(23BCc0GY_d%)|&1{|n3jAPE+z{V*|Z1_iJ>MySUi?qg&~NlnX1EJ2<*^kNVI literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/4_vds.h5 b/crates/clawhdf5/tests/fixtures/vds/4_vds.h5 new file mode 100644 index 0000000000000000000000000000000000000000..64c2288f7e6138ce0ba6f643fa22e33cb6e1547f GIT binary patch literal 5496 zcmeHHJxc>Y5S_if;1Up{AOx#xY}G=+QZc6~#E4(mMKl3B!Ba_Vr+>oQe_(j@@yutA@$a~R3poYp z2Nv}H+KfZmP5(!;`P4hTD=ymu;;~o}%n%nEMYj3BwFACDXTd3(B z6aZN7FJG=)Z_sLWbANj^OyXwIqoIf~{|LxqD8@aE^P6LtVHq(1LkzhhM;%MylVnl- t`_Rey$JvX>cih>Y8;2Gc1IBX8(fH7bU7z4Fo;1}6}Rr>${ literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/README.md b/crates/clawhdf5/tests/fixtures/vds/README.md new file mode 100644 index 0000000..25a1371 --- /dev/null +++ b/crates/clawhdf5/tests/fixtures/vds/README.md @@ -0,0 +1,14 @@ +# VDS test files from libhdf5 + +Copied unchanged from the HDF Group's HDF5 repository, +`tools/test/testfiles/vds/` (the h5dump/h5ls VDS test data). HDF5 is +distributed under a BSD-style license (see `COPYING` in the HDF5 source). + +| File | What it exercises | +|---|---| +| `vds-eiger.h5` + `f-0.h5`, `f-3.h5` | printf-style source name `f-%b.h5`; `f-3.h5` lies past the first missing source and must be ignored (extent 5, not 20) | +| `4_vds.h5` + `4_0.h5`..`4_2.h5` | printf-style `4_%b.h5` with version-2 (1.10 format) hyperslab selections | +| `vds-percival-unlim-maxmin.h5` + `a.h5`..`d.h5` | four interleaved unlimited mappings whose sources have different lengths | + +Used by `crates/clawhdf5/tests/vds_interop.rs::vds_libhdf5_test_files`, +which compares our reads with h5py's. diff --git a/crates/clawhdf5/tests/fixtures/vds/a.h5 b/crates/clawhdf5/tests/fixtures/vds/a.h5 new file mode 100644 index 0000000000000000000000000000000000000000..fa1953587aa29feb90477ff684cb3e34311c7cfd GIT binary patch literal 7736 zcmeI0!Ab)$5QZnowuFMvdQp0iK7z*{Ep%m1Dk{>WkKrTeOL+Fudv9KSByYNtncq@H zuu|AV|16tHCX)&I__OTpQ!_t57>oy^C1f7||QWO0c962vdl z-;Xv25D#6#!oL~K9tn^D36KB@kN^pg011#lmjv{Ck2e`m1JnRDKn-lG0sKg$UH??= Sl0}#A5?q@ENPq;kA@B)9LQp6G literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/b.h5 b/crates/clawhdf5/tests/fixtures/vds/b.h5 new file mode 100644 index 0000000000000000000000000000000000000000..08449ca3909de3755669ea173f0e637a533cb794 GIT binary patch literal 7736 zcmeH~%}PTt5QQhnZ3zXT6+v+!eFT?Xs`ScTsZm^T{COzKYp1&(p-3uR1kB+B$qj3+%Y2in#L_vMu1&N>DcY6Tp2NfPh7J&FSAil|a z-kje+Jh+5~e?FQkDnJFO02QDDRDcRl0V+TRmZ^aL?PaG1WPl8i0W$Dg1~$s~2Rk(& S179<+%-?oAn+p7Q1zrIkUtMwl literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/c.h5 b/crates/clawhdf5/tests/fixtures/vds/c.h5 new file mode 100644 index 0000000000000000000000000000000000000000..ba9af3051ac5f5d045dad907559abe237f0a68da GIT binary patch literal 7736 zcmeI0!Ab)$5QZnowuXYxdQm(`AECz{6}z%0l`7)V$M6x_m+0AxcW+*OByYMinO~_z zP$?AYpCy@OGMSK{FH5`gW^%UMJM4*;oTZY!b=#94|z#p=Fe5f-j#wEArj z*=yHn(M8XXiilA&HyZbFnAX0%^$N(8217zU04Cv>4R%$>7 Q$iRPKV3WVKx%a99??BLUi~s-t literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/d.h5 b/crates/clawhdf5/tests/fixtures/vds/d.h5 new file mode 100644 index 0000000000000000000000000000000000000000..8eceb4af89bff0b062fa2c3040f0bcd661dacd3e GIT binary patch literal 7736 zcmeH}u}T9$5Qb-WFJTD?iH%|*DJ?BiD%{CcDk@@`$H*h-OIX|3+u8a^cDge=zY;|- z60pd>VRv?CXXey=uc$V z{S&YMHSX%7>!d`YKCGw7VN+)r4%&@Upzn^)^}siNYEb&(oI#A%C*ZNi#vg5vq48S> zC3I$$LaI*ykM#nHz=B__VN>B>Gsh*dB0aTw1Igx(DO%wy)t_Rs7+8dbYZ>i+8$`BC zo;FQ-e!Rc_RLIc} z2I9dbEd0ySTu}ikKn17(6`%rCfC^9nD)65Q=>K$fYCs0a02v?yt1@s78{ZWRwgtt)2y;I*C(1p`&Y|l};)u;&`8$|QQ017iQG18q}n@WMyBHB!4btr88#Gqr8BjlgMkVtCk^HDJ8ez_h19Da5#yVbhJ liRL*hOy@_+0n!2KKz}=ce+YS{2f~0bAPfit!a$E1cmwT^HAesd literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/vds-eiger.h5 b/crates/clawhdf5/tests/fixtures/vds/vds-eiger.h5 new file mode 100644 index 0000000000000000000000000000000000000000..23d1fd32d6b3dae1a25aca7b84be19a130ded793 GIT binary patch literal 5496 zcmeH{y-EW?5XWclF0w?#C<#r1B_ug=u~&N*Y}oy zF85LJOB-lT^s5sqbWb{jTsyB=KfORiu+T4J*gWZ<-)toztEAh~p+QmdS(ZuTX!oSZ z#|>RV*IQ0?JUF_Te;E>uB%akYIPRBxqj3-0`K(s{Tk9Mh^m}rm_PN#5C3?m-k+|t> zTy2eZsY^^6S*~K`^FOeTrjNk{^Z6}BoXTXN2`sev)AP;uod?v^i}2h!0Vm)DoPZN> V0#3jQH~}Z%1e|~qa01^Em;ePNd literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/fixtures/vds/vds-percival-unlim-maxmin.h5 b/crates/clawhdf5/tests/fixtures/vds/vds-percival-unlim-maxmin.h5 new file mode 100644 index 0000000000000000000000000000000000000000..b7f8827c74ca881b765574d74a236d322ec5eecc GIT binary patch literal 5496 zcmeHIJxc>Y5S{yo96=!pf>@ z(wp6R7D7mM;V`oyJF~O1v-@7&Zl3zl@m77mPTYckd}WVrXmydVPzS0eCq`Ws9h@t`f+#2fG-7gOYE2EOolX$M&a=G zZjwGEv+#aCOS3Rb@+_T80(9EJf+g0W4WnP%u|f~c3dggRa9ivJB7%i}F^2UD{qxF3 zA!IA*rm0t>z{FkGHH8DYC%%{0OoAL>r^w^RNuO`eW1^`Uuc{lJ4I;kL_7p@-Z&f@0 zuYFE>gH!WEi-3$@Kzs+E^4?L>o*qz_Smyk7IF6jc6jx2+__io?o~j z4hDsz@Unu-=C2v(Q3ZekP+34H4767Ppm;jV0(!$hqY41UeJKm*$7>t^FHm`F6P{XS YKp9X5lmTTx8Bhk40cAiLSa$|~0Ko{0Bme*a literal 0 HcmV?d00001 diff --git a/crates/clawhdf5/tests/vds_interop.rs b/crates/clawhdf5/tests/vds_interop.rs index 032d77c..2346a5e 100644 --- a/crates/clawhdf5/tests/vds_interop.rs +++ b/crates/clawhdf5/tests/vds_interop.rs @@ -289,3 +289,162 @@ expect("nested.h5", "v", "nested") // A relative name below the virtual file's directory resolves there. assert_matches_libhdf5(dir.path(), "nested.h5", "v", "nested"); } + +// --------------------------------------------------------------------------- +// Unlimited and printf-style mappings +// --------------------------------------------------------------------------- + +/// Helpers for building unlimited VDS mappings through h5py's low-level API. +const UNLIMITED_HELPERS: &str = r#" +U = h5py.h5s.UNLIMITED +def space(dims, maxdims, start=None, count=None, stride=None, block=None): + s = h5py.h5s.create_simple(dims, maxdims) + if start is not None: + s.select_hyperslab(start, count, stride, block) + return s +def make_vds(fn, name, dims, maxdims, maps, fill, libver="latest", mode="w"): + # maps: [(vsel_kwargs, source_file, source_dataset, source_space)] + with h5py.File(fn, mode, libver=libver) as f: + dcpl = h5py.h5p.create(h5py.h5p.DATASET_CREATE) + for vsel, sfile, sdset, sspace in maps: + dcpl.set_virtual(space(dims, maxdims, **vsel), sfile.encode(), sdset.encode(), sspace) + dcpl.set_fill_value(np.array(fill, dtype="f8")) + h5py.h5d.create(f.id, name.encode(), h5py.h5t.IEEE_F64LE, + h5py.h5s.create_simple(dims, maxdims), dcpl=dcpl) +"#; + +/// printf-style names: block `j` of the virtual selection comes from the +/// source named with `j` in place of `%b` (`%%` is a literal `%`), probing +/// j = 0, 1, ... until the first missing source. libhdf5 also recomputes the +/// extent from what it finds, so the stored dataspace is not the shape. +#[test] +fn vds_printf_source_names() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let body = format!( + "{UNLIMITED_HELPERS}{}", + r#" +for i in [0, 1, 2, 4]: # 3 is missing: 4 is past the first gap and unused + with h5py.File(f"vds_src_{i}.h5", "w") as s: + s.create_dataset("data", data=np.arange(10.0) + i * 100) + with h5py.File(f"p%c_{i}.h5", "w") as s: + s.create_dataset("data", data=np.arange(10.0) - i * 100) +for libver in ["earliest", "latest"]: + fn = f"printf_{libver}.h5" + make_vds(fn, "files", (10,), (U,), + [(dict(start=(0,), count=(U,), stride=(10,), block=(10,)), "vds_src_%b.h5", "data", + space((10,), (10,), (0,), (1,), (1,), (10,)))], -1.0, libver) + # interleaved blocks with gaps between them, and an escaped percent sign + make_vds(fn, "escaped", (4,), (U,), + [(dict(start=(1,), count=(U,), stride=(6,), block=(4,)), "p%%c_%b.h5", "data", + space((10,), (10,), (2,), (1,), (1,), (4,)))], -5.0, libver, "a") + # printf in the dataset name, same file, 2-D frames + with h5py.File(fn, "a") as f: + for i in range(3): + f.create_dataset(f"frame_{i}", data=np.arange(6.0).reshape(2, 3) + 10 * i) + with h5py.File(fn, "a", libver=libver) as f: + dcpl = h5py.h5p.create(h5py.h5p.DATASET_CREATE) + dcpl.set_virtual(space((1, 2, 3), (U, 2, 3), (0, 0, 0), (U, 1, 1), (1, 1, 1), (1, 2, 3)), + b".", b"frame_%b", space((2, 3), (2, 3))) + h5py.h5d.create(f.id, b"frames", h5py.h5t.IEEE_F64LE, + h5py.h5s.create_simple((1, 2, 3), (U, 2, 3)), dcpl=dcpl) + for name in ["files", "escaped", "frames"]: + expect(fn, name, f"{name}_{libver}") +"# + ); + generate(dir.path(), &body); + for libver in ["earliest", "latest"] { + for name in ["files", "escaped", "frames"] { + let file = format!("printf_{libver}.h5"); + assert_matches_libhdf5(dir.path(), &file, name, &format!("{name}_{libver}")); + } + } + // The shape libhdf5 reports: three 10-element blocks. + let f = File::open(dir.path().join("printf_latest.h5")).unwrap(); + assert_eq!(f.dataset("files").unwrap().shape().unwrap(), vec![30]); +} + +/// Unlimited source and virtual selections: each mapping covers as much as +/// its source's current extent fills (a partial last block included), the +/// extent is the largest of them but never smaller than the limited +/// mappings need, and a missing source contributes nothing. +#[test] +fn vds_unlimited_mappings_follow_source_extents() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let body = format!( + "{UNLIMITED_HELPERS}{}", + r#" +with h5py.File("grow.h5", "w") as s: + s.create_dataset("a", data=np.arange(7.0) + 1, maxshape=(None,)) + s.create_dataset("b", data=np.arange(5.0) + 100, maxshape=(None,)) + s.create_dataset("rows", data=np.arange(12.0).reshape(4, 3) + 50, maxshape=(None, 3)) +unlim_src = lambda: space((1,), (U,), (0,), (U,), (1,), (1,)) +for libver in ["earliest", "latest"]: + fn = f"unlim_{libver}.h5" + make_vds(fn, "interleaved", (1,), (U,), [ + # blocks of 3 every 4: 7 source elements end mid-block + (dict(start=(0,), count=(U,), stride=(4,), block=(3,)), "grow.h5", "a", unlim_src()), + (dict(start=(3,), count=(U,), stride=(4,), block=(1,)), "grow.h5", "b", unlim_src()), + (dict(start=(0,), count=(U,), stride=(1,), block=(1,)), "missing.h5", "a", unlim_src()), + ], -2.0, libver) + make_vds(fn, "rows", (6, 3), (U, 3), [ + # an unlimited *block*, plus a limited mapping reaching row 5 + (dict(start=(0, 0), count=(1, 1), stride=(1, 1), block=(U, 3)), "grow.h5", "rows", + space((1, 3), (U, 3), (0, 0), (1, 1), (1, 1), (U, 3))), + (dict(start=(5, 0), count=(1, 1), stride=(1, 1), block=(1, 3)), "grow.h5", "rows", + space((4, 3), (U, 3), (1, 0), (1, 1), (1, 1), (1, 3))), + ], -4.0, libver, "a") + for name in ["interleaved", "rows"]: + expect(fn, name, f"{name}_{libver}") +"# + ); + generate(dir.path(), &body); + for libver in ["earliest", "latest"] { + for name in ["interleaved", "rows"] { + let file = format!("unlim_{libver}.h5"); + assert_matches_libhdf5(dir.path(), &file, name, &format!("{name}_{libver}")); + } + } + // "a" (7 elements, blocks of 3 every 4) ends at 9; "b" (5 elements from + // 3, every 4) at 20. "rows" fills 4 rows but a limited mapping needs 6. + let f = File::open(dir.path().join("unlim_latest.h5")).unwrap(); + assert_eq!(f.dataset("interleaved").unwrap().shape().unwrap(), vec![20]); + assert_eq!(f.dataset("rows").unwrap().shape().unwrap(), vec![6, 3]); +} + +/// libhdf5's own VDS test files (HDF5 `tools/test/testfiles/vds`): +/// printf-style Eiger frames (with a source past the first gap that must be +/// ignored), a printf mapping in the 1.10 format, and Percival's four +/// interleaved unlimited sources of different lengths. +#[test] +fn vds_libhdf5_test_files() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let fixtures = Path::new(env!("CARGO_MANIFEST_DIR")).join("tests/fixtures/vds"); + for entry in std::fs::read_dir(&fixtures).unwrap() { + let path = entry.unwrap().path(); + if path.extension().is_some_and(|e| e == "h5") { + std::fs::copy(&path, dir.path().join(path.file_name().unwrap())).unwrap(); + } + } + let cases = [ + ("vds-eiger.h5", "/VDS-Eiger"), + ("4_vds.h5", "/vds_dset"), + ("vds-percival-unlim-maxmin.h5", "/VDS-Percival-unlim-maxmin"), + ]; + let mut body = String::new(); + for (i, (file, dset)) in cases.iter().enumerate() { + body.push_str(&format!("expect({file:?}, {dset:?}, \"case{i}\")\n")); + } + generate(dir.path(), &body); + for (i, (file, dset)) in cases.iter().enumerate() { + assert_matches_libhdf5(dir.path(), file, dset, &format!("case{i}")); + } + // Stored as 20 frames; only f-0.h5 is found before the first gap. + let f = File::open(dir.path().join("vds-eiger.h5")).unwrap(); + assert_eq!( + f.dataset("VDS-Eiger").unwrap().shape().unwrap(), + vec![5, 10, 10] + ); +} diff --git a/docs/known-issues.md b/docs/known-issues.md index b1de913..dc5a005 100644 --- a/docs/known-issues.md +++ b/docs/known-issues.md @@ -72,7 +72,13 @@ fill-value item that did is fixed). - ~~**Wrong data:** unmapped regions read as 0 instead of the fill value.~~ Fixed 2026-09-25: unmapped elements and missing sources read as the virtual dataset's fill value. - - `%b` printf-style source names are not expanded. + - ~~`%b` printf-style source names are not expanded.~~ Fixed 2026-09-25: + printf-style and unlimited mappings are read, and the extent is + recomputed from the sources as libhdf5 does. Still open: the + "first missing" view and a printf gap other than 0 (libhdf5 access + properties we always read at their defaults), source-to-virtual type + conversion other than a byte swap, nested virtual sources, and source + files outside the virtual file's directory (refused with an error). - ~~Hyperslab selection versions 1 and 2 are refused.~~ Fixed 2026-09-25: versions 1-3 and irregular hyperslabs are decoded. - ~~The version-1 mapping list written with a 2.0 low bound (flags byte,