diff --git a/CHANGELOG.md b/CHANGELOG.md index f3f46dc..f8568f7 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -331,6 +331,47 @@ The legacy per-member dimension fields are now decoded into an array type, as libhdf5 does; more than four dimensions, or a zero-sized one, is an error. +- `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. + - 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 + 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. + - 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 ece868c..60a8243 100644 --- a/crates/clawhdf5-format/src/data_layout.rs +++ b/crates/clawhdf5-format/src/data_layout.rs @@ -74,21 +74,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). @@ -114,7 +126,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> { @@ -134,17 +146,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)?; @@ -1035,6 +1087,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-format/src/data_read.rs b/crates/clawhdf5-format/src/data_read.rs index 8d3b16e..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,160 +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::split_user_block; - use crate::superblock::Superblock; - - // An external source file is handed over whole, user block included; - // its addresses are relative to its superblock. - let (_, file_data) = split_user_block(file_data)?; - let sb = Superblock::parse(file_data, 0)?; - 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 ffed10f..ef09c1b 100644 --- a/crates/clawhdf5-format/src/lib.rs +++ b/crates/clawhdf5-format/src/lib.rs @@ -101,6 +101,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/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-format/src/vds.rs b/crates/clawhdf5-format/src/vds.rs new file mode 100644 index 0000000..23a3a52 --- /dev/null +++ b/crates/clawhdf5-format/src/vds.rs @@ -0,0 +1,1100 @@ +//! 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()) +} + +/// 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: 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. +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(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. 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>, +) -> Result, FormatError> { + let mappings = load_mappings(file_data, layout, length_size)?; + if mappings.iter().all(|m| m.kind == Kind::Fixed) { + return Ok(dataspace.dimensions.clone()); + } + let mut sources = Sources::new(file_data, resolver); + Ok(plan(&mappings, &dataspace.dimensions, &mut sources)?.dims) +} + +/// 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 +/// 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)?; + 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 + .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"))?]; + + 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; + Ok(VirtualData { + dims, + data, + unmapped, + }) +} + +/// `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], + 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)?)); + } + // An external file is handed over whole; its addresses are relative + // to its superblock, so skip any user block. + match self.cached_file.as_ref().and_then(|(_, b)| b.as_deref()) { + Some(bytes) => Ok(Some(crate::signature::split_user_block(bytes)?.1)), + None => Ok(None), + } + } + + /// 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. + fn dataset( + &mut self, + file: &str, + path: &str, + datatype: &Datatype, + ) -> Result, FormatError> { + let Some(bytes) = self.file(file)? else { + return Ok(None); + }; + let Some(src) = open_source(bytes, path)? else { + return Ok(None); + }; + read_source(bytes, src, path, datatype).map(Some) + } +} + +/// 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, + 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::shared_message::message_data_with_sohm; + + // `file_data` starts at the superblock (see `Sources::file`). + let sb = crate::superblock::Superblock::parse(file_data, 0)?; + 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 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 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 \ + (only a byte-order conversion is supported)" + ))); + }; + 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 { .. }) { + return Err(vds_err( + "virtual dataset source is itself virtual (unsupported)", + )); + } + let pipeline = src + .header + .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( + &src.header.messages, + file_data, + &layout, + &src.dataspace, + src_type.type_size() as usize, + os, + ls, + || { + crate::data_read::read_raw_data_full( + file_data, + &layout, + &src.dataspace, + &src_type, + pipeline.as_ref(), + os, + ls, + ) + }, + )?; + 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-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 04e3c2b..65c0292 100644 --- a/crates/clawhdf5/src/reader.rs +++ b/crates/clawhdf5/src/reader.rs @@ -522,6 +522,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 { @@ -855,7 +856,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 { @@ -884,24 +900,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. @@ -927,6 +928,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/fixtures/vds/4_0.h5 b/crates/clawhdf5/tests/fixtures/vds/4_0.h5 new file mode 100644 index 0000000..5e71d20 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/4_0.h5 differ 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 0000000..edad46e Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/4_1.h5 differ diff --git a/crates/clawhdf5/tests/fixtures/vds/4_2.h5 b/crates/clawhdf5/tests/fixtures/vds/4_2.h5 new file mode 100644 index 0000000..a82b012 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/4_2.h5 differ 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 0000000..64c2288 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/4_vds.h5 differ 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 0000000..fa19535 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/a.h5 differ diff --git a/crates/clawhdf5/tests/fixtures/vds/b.h5 b/crates/clawhdf5/tests/fixtures/vds/b.h5 new file mode 100644 index 0000000..08449ca Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/b.h5 differ diff --git a/crates/clawhdf5/tests/fixtures/vds/c.h5 b/crates/clawhdf5/tests/fixtures/vds/c.h5 new file mode 100644 index 0000000..ba9af30 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/c.h5 differ diff --git a/crates/clawhdf5/tests/fixtures/vds/d.h5 b/crates/clawhdf5/tests/fixtures/vds/d.h5 new file mode 100644 index 0000000..8eceb4a Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/d.h5 differ diff --git a/crates/clawhdf5/tests/fixtures/vds/f-0.h5 b/crates/clawhdf5/tests/fixtures/vds/f-0.h5 new file mode 100644 index 0000000..242df07 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/f-0.h5 differ diff --git a/crates/clawhdf5/tests/fixtures/vds/f-3.h5 b/crates/clawhdf5/tests/fixtures/vds/f-3.h5 new file mode 100644 index 0000000..539c1c9 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/f-3.h5 differ 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 0000000..23d1fd3 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/vds-eiger.h5 differ 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 0000000..b7f8827 Binary files /dev/null and b/crates/clawhdf5/tests/fixtures/vds/vds-percival-unlim-maxmin.h5 differ diff --git a/crates/clawhdf5/tests/vds_interop.rs b/crates/clawhdf5/tests/vds_interop.rs new file mode 100644 index 0000000..2346a5e --- /dev/null +++ b/crates/clawhdf5/tests/vds_interop.rs @@ -0,0 +1,450 @@ +//! 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"); +} + +// --------------------------------------------------------------------------- +// 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"); +} + +// --------------------------------------------------------------------------- +// 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"); +} + +// --------------------------------------------------------------------------- +// 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 a98c2f1..4ea7518 100644 --- a/docs/known-issues.md +++ b/docs/known-issues.md @@ -73,8 +73,8 @@ 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. @@ -85,9 +85,20 @@ the VDS item, which is marked. per-member dimensions were skipped, so an array member read as one scalar. **Fixed 2026-09-25.** - **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. + - ~~**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.~~ 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, + shared names) was misparsed.~~ Found and fixed 2026-09-25. - **Files with a user block:** the base address is not applied. **Fixed 2026-09-25:** every reader views the file from the superblock on (`twithub.h5`, `twithub513.h5`, `h5clear_fsm_persist_user_*.h5`; the