diff --git a/CHANGELOG.md b/CHANGELOG.md index de54b91..aabbce6 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -290,6 +290,19 @@ **Behaviour change:** the raw-read API (`read_raw_data_full*`), which has no fill value, now returns an error for a virtual dataset with unmapped elements instead of zeros. + - Unlimited and printf-style mappings are supported (all 7 VDS files in the + libhdf5 test set are such mappings, e.g. Eiger/Percival detector layouts). + `%b` in a source file or dataset name is the block number and `%%` a + literal `%` (other `%` sequences are an error, as in libhdf5); block `j` + is read from the source named with `j`, probing from 0 up to the first + missing source. Unlimited source/virtual selections cover as much as the + source's current extent fills, including a partial last block. As + libhdf5 does on `H5Dget_space`, the extent is recomputed from the sources + present (default "last available" view, printf gap 0) — + `vds::virtual_dataset_extent`, used by `Dataset::shape()` — so e.g. + `vds-eiger.h5` is `[5, 10, 10]`, not its stored `[20, 10, 10]`. A source + stored in the other byte order is byte-swapped (libhdf5 converts); + other type conversions remain an error. - Hyperslab selection versions 1 and 2 were refused ("only version-3 hyperslab selections are supported"). Version 1 is what libhdf5 writes for every VDS created with the default format bounds (h5py's default), so diff --git a/crates/clawhdf5-format/src/vds.rs b/crates/clawhdf5-format/src/vds.rs index 3ac0232..099c14a 100644 --- a/crates/clawhdf5-format/src/vds.rs +++ b/crates/clawhdf5-format/src/vds.rs @@ -45,12 +45,149 @@ fn vds_err(msg: impl Into) -> FormatError { FormatError::ChunkedReadError(msg.into()) } +/// Upper bound on the printf-style source datasets probed for one mapping. +const MAX_PRINTF_BLOCKS: u64 = 1 << 20; + +/// A source file or dataset name, parsed for printf-style `%b` block-number +/// substitutions (`H5D_virtual_parse_source_name`): `%b` is the block +/// number, `%%` a literal `%`, and any other `%` sequence is invalid. +#[derive(Debug, Clone, PartialEq)] +struct SourceName { + /// Literal text around the substitutions: `segments.len() == subs + 1`. + segments: Vec, +} + +impl SourceName { + fn parse(name: &str) -> Result { + let mut segments = vec![String::new()]; + let mut chars = name.chars(); + while let Some(c) = chars.next() { + if c != '%' { + segments.last_mut().expect("never empty").push(c); + continue; + } + match chars.next() { + Some('b') => segments.push(String::new()), + Some('%') => segments.last_mut().expect("never empty").push('%'), + _ => { + return Err(vds_err(format!( + "invalid format specifier in VDS source name {name:?}" + ))); + } + } + } + Ok(SourceName { segments }) + } + + /// Number of `%b` substitutions. + fn subs(&self) -> usize { + self.segments.len() - 1 + } + + /// The name with every `%b` replaced by `block`. + fn build(&self, block: u64) -> String { + let mut out = String::new(); + for (i, seg) in self.segments.iter().enumerate() { + if i > 0 { + out.push_str(&format!("{block}")); + } + out.push_str(seg); + } + out + } +} + +/// How a mapping's selections relate, as libhdf5 classifies them. +#[derive(Debug, Clone, Copy, PartialEq)] +enum Kind { + /// Both selections have a fixed size. + Fixed, + /// Both are unlimited (in `vdim` / `sdim`): the mapping grows with the + /// source dataset's extent. + Unlimited { vdim: usize, sdim: usize }, + /// The virtual selection repeats a block without limit in `vdim`; block + /// `j` comes from the source named by substituting `j` for `%b`. + Printf { vdim: usize }, +} + /// One mapping with its selections decoded. struct Mapping { - file: String, - dataset: String, + file: SourceName, + dataset: SourceName, vsel: SerializedSelection, ssel: SerializedSelection, + kind: Kind, +} + +impl Mapping { + fn new(m: VdsMapping) -> Result { + let (vsel, _) = SerializedSelection::decode(&m.virtual_selection)?; + let (ssel, _) = SerializedSelection::decode(&m.source_selection)?; + let file = SourceName::parse(&m.source_file)?; + let dataset = SourceName::parse(&m.source_dataset)?; + let subs = file.subs() + dataset.subs(); + // The checks of H5D_virtual_check_mapping_pre/_post. + let kind = match (vsel.unlimited_dim(), ssel.unlimited_dim()) { + (Some(vdim), None) => { + if subs == 0 { + return Err(vds_err( + "unlimited virtual selection with a limited source selection \ + and no %b in the source names", + )); + } + match &vsel { + SerializedSelection::Regular { count, block, .. } + if count[vdim] == UNLIMITED && block[vdim] != UNLIMITED => {} + _ => { + return Err(vds_err( + "printf VDS mapping needs a virtual selection with an unlimited count", + )); + } + } + Kind::Printf { vdim } + } + (Some(vdim), Some(sdim)) => { + if non_unlimited_elements(&vsel, vdim) != non_unlimited_elements(&ssel, sdim) { + return Err(vds_err( + "unlimited VDS mapping: virtual and source selections differ \ + outside the unlimited dimension", + )); + } + Kind::Unlimited { vdim, sdim } + } + (None, Some(_)) => { + return Err(vds_err( + "VDS mapping with an unlimited source selection and a limited \ + virtual selection is not supported", + )); + } + (None, None) => Kind::Fixed, + }; + if subs > 0 && !matches!(kind, Kind::Printf { .. }) { + return Err(vds_err( + "%b in a VDS source name without an unlimited virtual selection", + )); + } + Ok(Mapping { + file, + dataset, + vsel, + ssel, + kind, + }) + } +} + +/// Elements a regular selection selects outside dimension `skip`. +fn non_unlimited_elements(sel: &SerializedSelection, skip: usize) -> Option { + let SerializedSelection::Regular { count, block, .. } = sel else { + return None; + }; + (0..count.len()) + .filter(|&d| d != skip) + .try_fold(1u64, |acc, d| { + acc.checked_mul(count[d].checked_mul(block[d])?) + }) } /// Load and decode the mapping list of a virtual layout. @@ -82,49 +219,269 @@ fn load_mappings( })?; parse_vds_mappings(&obj.data, length_size)? .into_iter() - .map(|m: VdsMapping| { - let (vsel, _) = SerializedSelection::decode(&m.virtual_selection)?; - let (ssel, _) = SerializedSelection::decode(&m.source_selection)?; - Ok(Mapping { - file: m.source_file, - dataset: m.source_dataset, - vsel, - ssel, - }) - }) + .map(Mapping::new) .collect() } +/// `H5S__hyper_get_clip_diminfo`: the count and block a regular selection has +/// in its unlimited dimension once clipped to `clip`. +fn clip_diminfo(start: u64, stride: u64, count: u64, block: u64, clip: u64) -> (u64, u64) { + if start >= clip { + if block == UNLIMITED { + (count, 0) + } else { + (0, block) + } + } else if block == UNLIMITED || block == stride { + (1, clip - start) + } else { + ((clip - start).div_ceil(stride.max(1)), block) + } +} + +/// The unlimited-dimension parameters (start, stride, count, block) of a +/// regular selection. +fn unlim_diminfo(sel: &SerializedSelection, dim: usize) -> Result<[u64; 4], FormatError> { + match sel { + SerializedSelection::Regular { + start, + stride, + count, + block, + } => Ok([start[dim], stride[dim], count[dim], block[dim]]), + _ => Err(vds_err( + "unlimited VDS selection is not a regular hyperslab", + )), + } +} + +/// `H5S_hyper_get_clip_extent_match` with `incl_trail = false` (the +/// "last available" view): the extent to clip `clip_sel` (unlimited in +/// `clip_dim`) to so that it holds as many slices as `match_sel` (unlimited +/// in `match_dim`) holds when clipped to `match_clip`. +fn clip_extent_match( + clip_sel: &SerializedSelection, + clip_dim: usize, + match_sel: &SerializedSelection, + match_dim: usize, + match_clip: u64, +) -> Result { + let overflow = || FormatError::Overflow("VDS clip extent overflow".into()); + let [mstart, mstride, mcount, mblock] = unlim_diminfo(match_sel, match_dim)?; + let (count, block) = clip_diminfo(mstart, mstride, mcount, mblock, match_clip); + let slices = if block == 0 || count == 0 { + 0 + } else if count == 1 { + block + } else { + let mut n = block.checked_mul(count).ok_or_else(overflow)?; + let span = mstride + .checked_mul(count - 1) + .and_then(|s| s.checked_add(block)) + .ok_or_else(overflow)?; + let room = match_clip - mstart; + if span > room { + n -= span - room; + } + n + }; + + // H5S__hyper_get_clip_extent_real + let [start, stride, _, block] = unlim_diminfo(clip_sel, clip_dim)?; + if slices == 0 { + return Ok(0); + } + let extent = if block == UNLIMITED || block == stride { + start.checked_add(slices) + } else { + let full = slices / block; + let rem = slices - full * block; + if rem > 0 { + full.checked_mul(stride) + .and_then(|o| start.checked_add(o)) + .and_then(|e| e.checked_add(rem)) + } else { + (full - 1) + .checked_mul(stride) + .and_then(|o| start.checked_add(o)) + .and_then(|e| e.checked_add(block)) + } + }; + extent.ok_or_else(overflow) +} + +/// How each mapping is read, and the resulting extent. +struct Plan { + dims: Vec, + steps: Vec, +} + +#[derive(Clone, Copy)] +enum Step { + Fixed, + /// Read with the virtual selection clipped to `vclip` and the source + /// selection to the source's extent; `None` if the source is missing. + Unlimited(Option), + /// Read blocks `0..blocks`. + Printf(u64), +} + +/// Work out the extent libhdf5 gives the virtual dataset +/// (`H5D__virtual_set_extent_unlim`, default view `H5D_VDS_LAST_AVAILABLE` +/// with a printf gap of 0) and how much of each unlimited mapping is read. +fn plan(mappings: &[Mapping], stored: &[u64], sources: &mut Sources) -> Result { + let overflow = || FormatError::Overflow("VDS extent overflow".into()); + let rank = stored.len(); + let mut new_dims: Vec> = vec![None; rank]; + // Minimum extent needed by the limited parts of every virtual selection + // (H5D_virtual_update_min_dims). + let mut min_dims = vec![0u64; rank]; + let mut steps = Vec::with_capacity(mappings.len()); + + for m in mappings { + let skip = match m.kind { + Kind::Unlimited { vdim, .. } | Kind::Printf { vdim } => Some(vdim), + Kind::Fixed => None, + }; + if let Some(ends) = selection_bounds_end(&m.vsel)? { + if ends.len() != rank { + return Err(vds_err("VDS selection rank does not match dataspace rank")); + } + for d in (0..rank).filter(|&d| Some(d) != skip) { + min_dims[d] = min_dims[d].max(ends[d].checked_add(1).ok_or_else(overflow)?); + } + } + + let (vdim, clip, step) = match m.kind { + Kind::Fixed => { + steps.push(Step::Fixed); + continue; + } + Kind::Unlimited { vdim, sdim } => { + match sources.dims(&m.file.build(0), &m.dataset.build(0))? { + Some(src_dims) => { + let extent = *src_dims.get(sdim).ok_or_else(|| { + vds_err("VDS source rank does not match its selection") + })?; + let clip = clip_extent_match(&m.vsel, vdim, &m.ssel, sdim, extent)?; + (vdim, clip, Step::Unlimited(Some(clip))) + } + None => (vdim, 0, Step::Unlimited(None)), + } + } + Kind::Printf { vdim } => { + // With a gap of 0 the search stops at the first missing + // source dataset. + let mut found = 0u64; + while sources + .dims(&m.file.build(found), &m.dataset.build(found))? + .is_some() + { + found += 1; + if found > MAX_PRINTF_BLOCKS { + return Err(vds_err("too many printf-style VDS source datasets")); + } + } + let clip = if found == 0 { + 0 + } else { + // End of block `found - 1` in the unlimited dimension. + let [start, stride, _, block] = unlim_diminfo(&m.vsel, vdim)?; + (found - 1) + .checked_mul(stride) + .and_then(|o| start.checked_add(o)) + .and_then(|e| e.checked_add(block)) + .ok_or_else(overflow)? + }; + (vdim, clip, Step::Printf(found)) + } + }; + if vdim >= rank { + return Err(vds_err("VDS selection rank does not match dataspace rank")); + } + new_dims[vdim] = Some(new_dims[vdim].map_or(clip, |n| n.max(clip))); + steps.push(step); + } + + let dims = (0..rank) + .map(|d| match new_dims[d] { + None => stored[d], + Some(n) => n.max(min_dims[d]), + }) + .collect(); + Ok(Plan { dims, steps }) +} + +/// The last selected coordinate in each dimension (`H5S_SELECT_BOUNDS`), +/// ignoring any unlimited dimension; `None` for ALL/NONE and empty selections. +fn selection_bounds_end(sel: &SerializedSelection) -> Result>, FormatError> { + let overflow = || FormatError::Overflow("VDS selection bounds overflow".into()); + match sel { + SerializedSelection::All | SerializedSelection::None => Ok(None), + SerializedSelection::Regular { + start, + stride, + count, + block, + } => { + if count.contains(&0) || block.contains(&0) { + return Ok(None); + } + let mut ends = Vec::with_capacity(start.len()); + for d in 0..start.len() { + if count[d] == UNLIMITED || block[d] == UNLIMITED { + ends.push(0); + continue; + } + let end = (count[d] - 1) + .checked_mul(stride[d]) + .and_then(|o| start[d].checked_add(o)) + .and_then(|e| e.checked_add(block[d] - 1)) + .ok_or_else(overflow)?; + ends.push(end); + } + Ok(Some(ends)) + } + SerializedSelection::Blocks { rank, ends, .. } => { + if ends.is_empty() { + return Ok(None); + } + let mut max = vec![0u64; *rank]; + for e in ends.chunks_exact(*rank) { + for d in 0..*rank { + max[d] = max[d].max(e[d]); + } + } + Ok(Some(max)) + } + } +} + /// The virtual dataset's extent as libhdf5 reports it (`H5Dget_space`). /// /// For a virtual dataset whose mappings are all of fixed size this is the -/// stored dataspace. +/// stored dataspace. With unlimited or printf-style mappings libhdf5 +/// recomputes the unlimited dimension from the sources present (the default +/// "last available" view: the largest extent any mapping can fill), which +/// needs the source files, read through `resolver`. pub fn virtual_dataset_extent( file_data: &[u8], layout: &DataLayout, dataspace: &Dataspace, _offset_size: u8, length_size: u8, - _resolver: Option<&VdsFileResolver>, + resolver: Option<&VdsFileResolver>, ) -> Result, FormatError> { let mappings = load_mappings(file_data, layout, length_size)?; - check_fixed(&mappings)?; - Ok(dataspace.dimensions.clone()) -} - -fn check_fixed(mappings: &[Mapping]) -> Result<(), FormatError> { - if mappings - .iter() - .any(|m| m.vsel.unlimited_dim().is_some() || m.ssel.unlimited_dim().is_some()) - { - return Err(vds_err( - "unlimited virtual dataset mappings are not supported", - )); + if mappings.iter().all(|m| m.kind == Kind::Fixed) { + return Ok(dataspace.dimensions.clone()); } - Ok(()) + let mut sources = Sources::new(file_data, resolver); + Ok(plan(&mappings, &dataspace.dimensions, &mut sources)?.dims) } -/// Read a whole virtual dataset. +/// Read a whole virtual dataset, at the extent +/// [`virtual_dataset_extent`] reports. /// /// `fill` is the virtual dataset's fill value (from its fill value message; /// `None` for the default of zeros); every element no mapping supplies holds @@ -137,20 +494,13 @@ pub fn read_virtual_dataset( dataspace: &Dataspace, datatype: &Datatype, fill: Option<&[u8]>, - offset_size: u8, + _offset_size: u8, length_size: u8, resolver: Option<&VdsFileResolver>, ) -> Result { let mappings = load_mappings(file_data, layout, length_size)?; - check_fixed(&mappings)?; - let dims = virtual_dataset_extent( - file_data, - layout, - dataspace, - offset_size, - length_size, - resolver, - )?; + let mut sources = Sources::new(file_data, resolver); + let Plan { dims, steps } = plan(&mappings, &dataspace.dimensions, &mut sources)?; let elem_size = datatype.type_size() as usize; let total = dims @@ -167,14 +517,46 @@ pub fn read_virtual_dataset( } let mut mapped = vec![false; usize::try_from(total).map_err(|_| vds_err("VDS too large"))?]; - let mut sources = Sources::new(file_data, resolver); - for m in &mappings { - let Some(src) = sources.dataset(&m.file, &m.dataset, datatype)? else { - continue; // missing source file or dataset: fill - }; - let vidx = selection_indices(&m.vsel, &dims, None)?; - let sidx = selection_indices(&m.ssel, &src.dims, None)?; - scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + for (m, step) in mappings.iter().zip(&steps) { + match (*step, m.kind) { + (Step::Fixed, _) => { + let (file, dset) = (m.file.build(0), m.dataset.build(0)); + let Some(src) = sources.dataset(&file, &dset, datatype)? else { + continue; // missing source file or dataset: fill + }; + let vidx = selection_indices(&m.vsel, &dims, None)?; + let sidx = selection_indices(&m.ssel, &src.dims, None)?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + (Step::Unlimited(Some(vclip)), Kind::Unlimited { vdim, sdim }) => { + let (file, dset) = (m.file.build(0), m.dataset.build(0)); + let Some(src) = sources.dataset(&file, &dset, datatype)? else { + continue; + }; + let vidx = selection_indices(&m.vsel, &dims, Some((vdim, vclip)))?; + let extent = *src + .dims + .get(sdim) + .ok_or_else(|| vds_err("VDS source rank does not match its selection"))?; + let sidx = selection_indices(&m.ssel, &src.dims, Some((sdim, extent)))?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + (Step::Unlimited(None), _) => {} + (Step::Printf(blocks), Kind::Printf { vdim }) => { + for j in 0..blocks { + let Some(src) = + sources.dataset(&m.file.build(j), &m.dataset.build(j), datatype)? + else { + continue; + }; + let vblock = unlim_block(&m.vsel, vdim, j)?; + let vidx = selection_indices(&vblock, &dims, None)?; + let sidx = selection_indices(&m.ssel, &src.dims, None)?; + scatter(&mut data, &mut mapped, &src.raw, &vidx, &sidx, elem_size)?; + } + } + _ => return Err(vds_err("internal error: VDS plan does not match mapping")), + } } let unmapped = mapped.iter().filter(|&&m| !m).count() as u64; @@ -185,6 +567,37 @@ pub fn read_virtual_dataset( }) } +/// `H5S_hyper_get_unlim_block`: block `j` of a selection whose count is +/// unlimited in `dim`. +fn unlim_block( + sel: &SerializedSelection, + dim: usize, + j: u64, +) -> Result { + let SerializedSelection::Regular { + start, + stride, + count, + block, + } = sel + else { + return Err(vds_err("printf VDS selection is not a regular hyperslab")); + }; + let mut start = start.clone(); + let mut count = count.clone(); + start[dim] = j + .checked_mul(stride[dim]) + .and_then(|o| start[dim].checked_add(o)) + .ok_or_else(|| FormatError::Overflow("VDS block start overflow".into()))?; + count[dim] = 1; + Ok(SerializedSelection::Regular { + start, + stride: stride.clone(), + count, + block: block.clone(), + }) +} + /// Copy source element `sidx[i]` to virtual element `vidx[i]` for every `i`. fn scatter( out: &mut [u8], @@ -395,6 +808,15 @@ impl<'a, 'r> Sources<'a, 'r> { Ok(self.cached_file.as_ref().and_then(|(_, b)| b.as_deref())) } + /// The extent of source dataset `path` in file `file`, or `None` when + /// either does not exist. + fn dims(&mut self, file: &str, path: &str) -> Result>, FormatError> { + let Some(bytes) = self.file(file)? else { + return Ok(None); + }; + Ok(open_source(bytes, path)?.map(|s| s.dataspace.dimensions)) + } + /// Read source dataset `path` from file `file`, or `None` when either /// does not exist. Its datatype must be the virtual dataset's: libhdf5 /// converts between types here, which is not supported. @@ -407,20 +829,37 @@ impl<'a, 'r> Sources<'a, 'r> { let Some(bytes) = self.file(file)? else { return Ok(None); }; - read_source(bytes, path, datatype) + let Some(src) = open_source(bytes, path)? else { + return Ok(None); + }; + read_source(bytes, src, path, datatype).map(Some) } } -/// Read source dataset `path` of the file in `file_data` in full (its own -/// fill value applied to unallocated chunks), or `None` if it does not exist. -fn read_source( - file_data: &[u8], +/// An opened source dataset's object header. +struct OpenSource { + offset_size: u8, + length_size: u8, + header: crate::object_header::ObjectHeader, + dataspace: Dataspace, +} + +fn source_message<'h>( + src: &'h OpenSource, path: &str, - datatype: &Datatype, -) -> Result, FormatError> { - use crate::filter_pipeline::FilterPipeline; + t: crate::message_type::MessageType, +) -> Result<&'h crate::object_header::HeaderMessage, FormatError> { + src.header + .messages + .iter() + .find(|m| m.msg_type == t) + .ok_or_else(|| vds_err(format!("VDS source {path} has no {t:?} message"))) +} + +/// Open source dataset `path` of the file in `file_data`, or `None` if there +/// is no such object (libhdf5 reads a missing source as fill). +fn open_source(file_data: &[u8], path: &str) -> Result, FormatError> { use crate::message_type::MessageType; - use crate::object_header::ObjectHeader; use crate::shared_message::message_data_with_sohm; let sig = crate::signature::find_signature(file_data)?; @@ -431,30 +870,55 @@ fn read_source( Err(FormatError::PathNotFound(_)) => return Ok(None), Err(e) => return Err(e), }; - let hdr = ObjectHeader::parse(file_data, addr as usize, os, ls)?; - let msg = |t: MessageType| { - hdr.messages - .iter() - .find(|m| m.msg_type == t) - .ok_or_else(|| vds_err(format!("VDS source {path} has no {t:?} message"))) + let header = crate::object_header::ObjectHeader::parse(file_data, addr as usize, os, ls)?; + let mut src = OpenSource { + offset_size: os, + length_size: ls, + header, + dataspace: Dataspace { + space_type: crate::dataspace::DataspaceType::Null, + rank: 0, + dimensions: Vec::new(), + max_dimensions: None, + }, }; - let dataspace = Dataspace::parse( - &message_data_with_sohm(file_data, msg(MessageType::Dataspace)?, os, ls)?, - ls, - )?; - let (src_type, _) = Datatype::parse(&message_data_with_sohm( - file_data, - msg(MessageType::Datatype)?, - os, - ls, - )?)?; - if &src_type != datatype { + let ds_msg = source_message(&src, path, MessageType::Dataspace)?; + src.dataspace = Dataspace::parse(&message_data_with_sohm(file_data, ds_msg, os, ls)?, ls)?; + Ok(Some(src)) +} + +/// Read an opened source dataset in full (its own fill value applied to +/// unallocated chunks). +fn read_source( + file_data: &[u8], + src: OpenSource, + path: &str, + datatype: &Datatype, +) -> Result { + use crate::filter_pipeline::FilterPipeline; + use crate::message_type::MessageType; + use crate::shared_message::message_data_with_sohm; + + let (os, ls) = (src.offset_size, src.length_size); + let dt_msg = source_message(&src, path, MessageType::Datatype)?; + let (src_type, _) = Datatype::parse(&message_data_with_sohm(file_data, dt_msg, os, ls)?)?; + // libhdf5 converts each source to the virtual dataset's type. Only the + // conversion that is a pure byte swap is done here. + let swap = if &src_type == datatype { + false + } else if differs_only_in_byte_order(&src_type, datatype) { + true + } else { return Err(vds_err(format!( "VDS source {path} has a different datatype from the virtual dataset \ - (type conversion is not supported)" + (only a byte-order conversion is supported)" ))); - } - let layout = DataLayout::parse(&msg(MessageType::DataLayout)?.data, os, ls)?; + }; + let layout = DataLayout::parse( + &source_message(&src, path, MessageType::DataLayout)?.data, + os, + ls, + )?; // A source that is itself virtual could form a cycle (A -> B -> A) and // recurse without bound. Nested virtual sources are not supported. if matches!(layout, DataLayout::Virtual { .. }) { @@ -462,7 +926,8 @@ fn read_source( "virtual dataset source is itself virtual (unsupported)", )); } - let pipeline = hdr + let pipeline = src + .header .messages .iter() .find(|m| m.msg_type == MessageType::FilterPipeline) @@ -471,10 +936,10 @@ fn read_source( }) .transpose()?; let raw = crate::fill_value::read_full_with_fill( - &hdr.messages, + &src.header.messages, file_data, &layout, - &dataspace, + &src.dataspace, src_type.type_size() as usize, os, ls, @@ -482,7 +947,7 @@ fn read_source( crate::data_read::read_raw_data_full( file_data, &layout, - &dataspace, + &src.dataspace, &src_type, pipeline.as_ref(), os, @@ -490,8 +955,141 @@ fn read_source( ) }, )?; - Ok(Some(SourceData { - dims: dataspace.dimensions, + let mut raw = raw; + if swap { + let size = datatype.type_size() as usize; + for element in raw.chunks_exact_mut(size) { + element.reverse(); + } + } + Ok(SourceData { + dims: src.dataspace.dimensions, raw, - })) + }) +} + +/// Whether two numeric types are identical apart from their byte order (both +/// little- or big-endian), so converting one to the other is a byte swap. +fn differs_only_in_byte_order(a: &Datatype, b: &Datatype) -> bool { + use crate::datatype::DatatypeByteOrder::{BigEndian, LittleEndian}; + let mut a = a.clone(); + match &mut a { + Datatype::FixedPoint { byte_order, .. } + | Datatype::FloatingPoint { byte_order, .. } + | Datatype::BitField { byte_order, .. } => { + *byte_order = match byte_order { + LittleEndian => BigEndian, + BigEndian => LittleEndian, + _ => return false, + } + } + _ => return false, + } + &a == b +} + +#[cfg(test)] +mod tests { + use super::*; + + fn regular(start: u64, stride: u64, count: u64, block: u64) -> SerializedSelection { + SerializedSelection::Regular { + start: vec![start], + stride: vec![stride], + count: vec![count], + block: vec![block], + } + } + + #[test] + fn printf_names_follow_libhdf5() { + let n = SourceName::parse("f-%b.h5").unwrap(); + assert_eq!(n.subs(), 1); + assert_eq!(n.build(12), "f-12.h5"); + let n = SourceName::parse("100%%_%b_%b").unwrap(); + assert_eq!(n.build(3), "100%_3_3"); + let n = SourceName::parse("plain%%name").unwrap(); + assert_eq!((n.subs(), n.build(7)), (0, "plain%name".to_string())); + // Anything else after '%' (or a trailing '%') is invalid. + assert!(SourceName::parse("a%d").is_err()); + assert!(SourceName::parse("a%").is_err()); + } + + #[test] + fn clip_extent_matches_libhdf5_arithmetic() { + // 7 source slices (contiguous unlimited source) into blocks of 3 + // every 4: two full blocks and one slice of a third -> extent 9. + let src = regular(0, 1, UNLIMITED, 1); + let v = regular(0, 4, UNLIMITED, 3); + assert_eq!(clip_extent_match(&v, 0, &src, 0, 7).unwrap(), 9); + // Exactly two blocks: the extent ends at the end of the last block. + assert_eq!(clip_extent_match(&v, 0, &src, 0, 6).unwrap(), 7); + // An empty source gives an empty mapping. + assert_eq!(clip_extent_match(&v, 0, &src, 0, 0).unwrap(), 0); + // Unlimited block: the extent is start + slices. + let vb = regular(2, 1, 1, UNLIMITED); + assert_eq!(clip_extent_match(&vb, 0, &src, 0, 5).unwrap(), 7); + // A strided source clipped mid-block counts only the selected slices. + let src2 = regular(1, 4, UNLIMITED, 2); // 1,2, 5,6, 9,10 ... + let dense = regular(0, 1, UNLIMITED, 1); + assert_eq!(clip_extent_match(&dense, 0, &src2, 0, 6).unwrap(), 3); + } + + #[test] + fn clipped_selection_drops_the_partial_tail() { + let v = regular(0, 4, UNLIMITED, 3); + assert_eq!( + selection_indices(&v, &[20], Some((0, 9))).unwrap(), + vec![0, 1, 2, 4, 5, 6, 8] + ); + // Unclipped unlimited selections cannot be enumerated. + assert!(selection_indices(&v, &[20], None).is_err()); + } + + #[test] + fn unlimited_mapping_rules() { + let sel = |s: &SerializedSelection| -> Vec { + // Serialize as version 2 (8-byte regular). + let SerializedSelection::Regular { + start, + stride, + count, + block, + } = s + else { + unreachable!() + }; + let mut b = Vec::new(); + b.extend_from_slice(&2u32.to_le_bytes()); + b.extend_from_slice(&2u32.to_le_bytes()); + b.push(1); + b.extend_from_slice(&0u32.to_le_bytes()); + b.extend_from_slice(&(start.len() as u32).to_le_bytes()); + for d in 0..start.len() { + for v in [start[d], stride[d], count[d], block[d]] { + b.extend_from_slice(&v.to_le_bytes()); + } + } + b + }; + let mapping = |file: &str, v: &SerializedSelection, s: &SerializedSelection| VdsMapping { + source_file: file.into(), + source_dataset: "d".into(), + source_selection: sel(s), + virtual_selection: sel(v), + }; + let unlim = regular(0, 10, UNLIMITED, 10); + let fixed = regular(0, 1, 1, 10); + // Unlimited virtual + limited source needs %b in a name ... + assert!(Mapping::new(mapping("f.h5", &unlim, &fixed)).is_err()); + let m = Mapping::new(mapping("f%b.h5", &unlim, &fixed)).unwrap(); + assert_eq!(m.kind, Kind::Printf { vdim: 0 }); + // ... and %b is refused anywhere else. + assert!(Mapping::new(mapping("f%b.h5", &fixed, &fixed)).is_err()); + let src_unlim = regular(0, 1, UNLIMITED, 1); + let m = Mapping::new(mapping("f.h5", &unlim, &src_unlim)).unwrap(); + assert_eq!(m.kind, Kind::Unlimited { vdim: 0, sdim: 0 }); + // An unlimited source into a limited virtual selection is refused. + assert!(Mapping::new(mapping("f.h5", &fixed, &src_unlim)).is_err()); + } } diff --git a/crates/clawhdf5/src/reader.rs b/crates/clawhdf5/src/reader.rs index 2e7bd3d..20e2c33 100644 --- a/crates/clawhdf5/src/reader.rs +++ b/crates/clawhdf5/src/reader.rs @@ -809,7 +809,22 @@ impl<'f> Dataset<'f> { fn dataspace(&self) -> Result { let data = self.required_payload(MessageType::Dataspace)?; - Ok(Dataspace::parse(&data, self.file.length_size())?) + let mut ds = Dataspace::parse(&data, self.file.length_size())?; + // libhdf5 reports a virtual dataset with unlimited or printf-style + // mappings at the extent its sources currently fill, not the stored + // one (`H5Dget_space`). + if let Ok(dl @ DataLayout::Virtual { .. }) = self.data_layout() { + let resolver = self.vds_resolver(); + ds.dimensions = clawhdf5_format::vds::virtual_dataset_extent( + self.file.data.as_bytes(), + &dl, + &ds, + self.file.offset_size(), + self.file.length_size(), + Some(&resolver), + )?; + } + Ok(ds) } fn data_layout(&self) -> Result { diff --git a/crates/clawhdf5/tests/fixtures/vds/4_0.h5 b/crates/clawhdf5/tests/fixtures/vds/4_0.h5 new file mode 100644 index 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 index 032d77c..2346a5e 100644 --- a/crates/clawhdf5/tests/vds_interop.rs +++ b/crates/clawhdf5/tests/vds_interop.rs @@ -289,3 +289,162 @@ expect("nested.h5", "v", "nested") // A relative name below the virtual file's directory resolves there. assert_matches_libhdf5(dir.path(), "nested.h5", "v", "nested"); } + +// --------------------------------------------------------------------------- +// Unlimited and printf-style mappings +// --------------------------------------------------------------------------- + +/// Helpers for building unlimited VDS mappings through h5py's low-level API. +const UNLIMITED_HELPERS: &str = r#" +U = h5py.h5s.UNLIMITED +def space(dims, maxdims, start=None, count=None, stride=None, block=None): + s = h5py.h5s.create_simple(dims, maxdims) + if start is not None: + s.select_hyperslab(start, count, stride, block) + return s +def make_vds(fn, name, dims, maxdims, maps, fill, libver="latest", mode="w"): + # maps: [(vsel_kwargs, source_file, source_dataset, source_space)] + with h5py.File(fn, mode, libver=libver) as f: + dcpl = h5py.h5p.create(h5py.h5p.DATASET_CREATE) + for vsel, sfile, sdset, sspace in maps: + dcpl.set_virtual(space(dims, maxdims, **vsel), sfile.encode(), sdset.encode(), sspace) + dcpl.set_fill_value(np.array(fill, dtype="f8")) + h5py.h5d.create(f.id, name.encode(), h5py.h5t.IEEE_F64LE, + h5py.h5s.create_simple(dims, maxdims), dcpl=dcpl) +"#; + +/// printf-style names: block `j` of the virtual selection comes from the +/// source named with `j` in place of `%b` (`%%` is a literal `%`), probing +/// j = 0, 1, ... until the first missing source. libhdf5 also recomputes the +/// extent from what it finds, so the stored dataspace is not the shape. +#[test] +fn vds_printf_source_names() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let body = format!( + "{UNLIMITED_HELPERS}{}", + r#" +for i in [0, 1, 2, 4]: # 3 is missing: 4 is past the first gap and unused + with h5py.File(f"vds_src_{i}.h5", "w") as s: + s.create_dataset("data", data=np.arange(10.0) + i * 100) + with h5py.File(f"p%c_{i}.h5", "w") as s: + s.create_dataset("data", data=np.arange(10.0) - i * 100) +for libver in ["earliest", "latest"]: + fn = f"printf_{libver}.h5" + make_vds(fn, "files", (10,), (U,), + [(dict(start=(0,), count=(U,), stride=(10,), block=(10,)), "vds_src_%b.h5", "data", + space((10,), (10,), (0,), (1,), (1,), (10,)))], -1.0, libver) + # interleaved blocks with gaps between them, and an escaped percent sign + make_vds(fn, "escaped", (4,), (U,), + [(dict(start=(1,), count=(U,), stride=(6,), block=(4,)), "p%%c_%b.h5", "data", + space((10,), (10,), (2,), (1,), (1,), (4,)))], -5.0, libver, "a") + # printf in the dataset name, same file, 2-D frames + with h5py.File(fn, "a") as f: + for i in range(3): + f.create_dataset(f"frame_{i}", data=np.arange(6.0).reshape(2, 3) + 10 * i) + with h5py.File(fn, "a", libver=libver) as f: + dcpl = h5py.h5p.create(h5py.h5p.DATASET_CREATE) + dcpl.set_virtual(space((1, 2, 3), (U, 2, 3), (0, 0, 0), (U, 1, 1), (1, 1, 1), (1, 2, 3)), + b".", b"frame_%b", space((2, 3), (2, 3))) + h5py.h5d.create(f.id, b"frames", h5py.h5t.IEEE_F64LE, + h5py.h5s.create_simple((1, 2, 3), (U, 2, 3)), dcpl=dcpl) + for name in ["files", "escaped", "frames"]: + expect(fn, name, f"{name}_{libver}") +"# + ); + generate(dir.path(), &body); + for libver in ["earliest", "latest"] { + for name in ["files", "escaped", "frames"] { + let file = format!("printf_{libver}.h5"); + assert_matches_libhdf5(dir.path(), &file, name, &format!("{name}_{libver}")); + } + } + // The shape libhdf5 reports: three 10-element blocks. + let f = File::open(dir.path().join("printf_latest.h5")).unwrap(); + assert_eq!(f.dataset("files").unwrap().shape().unwrap(), vec![30]); +} + +/// Unlimited source and virtual selections: each mapping covers as much as +/// its source's current extent fills (a partial last block included), the +/// extent is the largest of them but never smaller than the limited +/// mappings need, and a missing source contributes nothing. +#[test] +fn vds_unlimited_mappings_follow_source_extents() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let body = format!( + "{UNLIMITED_HELPERS}{}", + r#" +with h5py.File("grow.h5", "w") as s: + s.create_dataset("a", data=np.arange(7.0) + 1, maxshape=(None,)) + s.create_dataset("b", data=np.arange(5.0) + 100, maxshape=(None,)) + s.create_dataset("rows", data=np.arange(12.0).reshape(4, 3) + 50, maxshape=(None, 3)) +unlim_src = lambda: space((1,), (U,), (0,), (U,), (1,), (1,)) +for libver in ["earliest", "latest"]: + fn = f"unlim_{libver}.h5" + make_vds(fn, "interleaved", (1,), (U,), [ + # blocks of 3 every 4: 7 source elements end mid-block + (dict(start=(0,), count=(U,), stride=(4,), block=(3,)), "grow.h5", "a", unlim_src()), + (dict(start=(3,), count=(U,), stride=(4,), block=(1,)), "grow.h5", "b", unlim_src()), + (dict(start=(0,), count=(U,), stride=(1,), block=(1,)), "missing.h5", "a", unlim_src()), + ], -2.0, libver) + make_vds(fn, "rows", (6, 3), (U, 3), [ + # an unlimited *block*, plus a limited mapping reaching row 5 + (dict(start=(0, 0), count=(1, 1), stride=(1, 1), block=(U, 3)), "grow.h5", "rows", + space((1, 3), (U, 3), (0, 0), (1, 1), (1, 1), (U, 3))), + (dict(start=(5, 0), count=(1, 1), stride=(1, 1), block=(1, 3)), "grow.h5", "rows", + space((4, 3), (U, 3), (1, 0), (1, 1), (1, 1), (1, 3))), + ], -4.0, libver, "a") + for name in ["interleaved", "rows"]: + expect(fn, name, f"{name}_{libver}") +"# + ); + generate(dir.path(), &body); + for libver in ["earliest", "latest"] { + for name in ["interleaved", "rows"] { + let file = format!("unlim_{libver}.h5"); + assert_matches_libhdf5(dir.path(), &file, name, &format!("{name}_{libver}")); + } + } + // "a" (7 elements, blocks of 3 every 4) ends at 9; "b" (5 elements from + // 3, every 4) at 20. "rows" fills 4 rows but a limited mapping needs 6. + let f = File::open(dir.path().join("unlim_latest.h5")).unwrap(); + assert_eq!(f.dataset("interleaved").unwrap().shape().unwrap(), vec![20]); + assert_eq!(f.dataset("rows").unwrap().shape().unwrap(), vec![6, 3]); +} + +/// libhdf5's own VDS test files (HDF5 `tools/test/testfiles/vds`): +/// printf-style Eiger frames (with a source past the first gap that must be +/// ignored), a printf mapping in the 1.10 format, and Percival's four +/// interleaved unlimited sources of different lengths. +#[test] +fn vds_libhdf5_test_files() { + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let fixtures = Path::new(env!("CARGO_MANIFEST_DIR")).join("tests/fixtures/vds"); + for entry in std::fs::read_dir(&fixtures).unwrap() { + let path = entry.unwrap().path(); + if path.extension().is_some_and(|e| e == "h5") { + std::fs::copy(&path, dir.path().join(path.file_name().unwrap())).unwrap(); + } + } + let cases = [ + ("vds-eiger.h5", "/VDS-Eiger"), + ("4_vds.h5", "/vds_dset"), + ("vds-percival-unlim-maxmin.h5", "/VDS-Percival-unlim-maxmin"), + ]; + let mut body = String::new(); + for (i, (file, dset)) in cases.iter().enumerate() { + body.push_str(&format!("expect({file:?}, {dset:?}, \"case{i}\")\n")); + } + generate(dir.path(), &body); + for (i, (file, dset)) in cases.iter().enumerate() { + assert_matches_libhdf5(dir.path(), file, dset, &format!("case{i}")); + } + // Stored as 20 frames; only f-0.h5 is found before the first gap. + let f = File::open(dir.path().join("vds-eiger.h5")).unwrap(); + assert_eq!( + f.dataset("VDS-Eiger").unwrap().shape().unwrap(), + vec![5, 10, 10] + ); +} diff --git a/docs/known-issues.md b/docs/known-issues.md index b1de913..dc5a005 100644 --- a/docs/known-issues.md +++ b/docs/known-issues.md @@ -72,7 +72,13 @@ fill-value item that did is fixed). - ~~**Wrong data:** unmapped regions read as 0 instead of the fill value.~~ Fixed 2026-09-25: unmapped elements and missing sources read as the virtual dataset's fill value. - - `%b` printf-style source names are not expanded. + - ~~`%b` printf-style source names are not expanded.~~ Fixed 2026-09-25: + printf-style and unlimited mappings are read, and the extent is + recomputed from the sources as libhdf5 does. Still open: the + "first missing" view and a printf gap other than 0 (libhdf5 access + properties we always read at their defaults), source-to-virtual type + conversion other than a byte swap, nested virtual sources, and source + files outside the virtual file's directory (refused with an error). - ~~Hyperslab selection versions 1 and 2 are refused.~~ Fixed 2026-09-25: versions 1-3 and irregular hyperslabs are decoded. - ~~The version-1 mapping list written with a 2.0 low bound (flags byte,