From 081341b433b15319ab7f50cef7a6c1d414f9f07f Mon Sep 17 00:00:00 2001 From: osobh Date: Fri, 25 Sep 2026 21:02:39 -0500 Subject: [PATCH] fix(format): convert float data read as integers instead of returning bit patterns read_i32/read_i64/read_u64 on a floating-point dataset reinterpreted the IEEE bits (1.5 read as i64 was 4609434218613702656). Convert like libhdf5's hard conversions instead: truncate toward zero and saturate at the target range; NaN reads as 0. Co-Authored-By: Claude Opus 5.5 (1M context) --- crates/clawhdf5-format/src/data_read.rs | 120 ++++++++++++---- .../tests/numeric_conversion_interop.rs | 129 ++++++++++++++++++ 2 files changed, 226 insertions(+), 23 deletions(-) create mode 100644 crates/clawhdf5/tests/numeric_conversion_interop.rs diff --git a/crates/clawhdf5-format/src/data_read.rs b/crates/clawhdf5-format/src/data_read.rs index e3c1bbd..c998fad 100644 --- a/crates/clawhdf5-format/src/data_read.rs +++ b/crates/clawhdf5-format/src/data_read.rs @@ -982,7 +982,74 @@ fn convert_to_f64( } } +/// One numeric element as stored, before conversion to the caller's type. +#[derive(Debug, Clone, Copy, PartialEq)] +enum Scalar { + Signed(i64), + Unsigned(u64), + Float(f64), +} + +impl Scalar { + /// Float to integer conversions follow libhdf5's hard conversions: + /// truncate toward zero, and saturate a value outside the target range to + /// its minimum or maximum. NaN converts to 0 (libhdf5 leaves that case to + /// the C cast, whose result is platform-dependent). + fn to_i64(self) -> i64 { + match self { + Scalar::Signed(v) => v, + Scalar::Unsigned(v) => v as i64, + Scalar::Float(v) => v as i64, + } + } + + fn to_u64(self) -> u64 { + match self { + Scalar::Signed(v) => v as u64, + Scalar::Unsigned(v) => v, + Scalar::Float(v) => v as u64, + } + } + + fn to_i32(self) -> i32 { + match self { + Scalar::Signed(v) => v as i32, + Scalar::Unsigned(v) => v as i32, + Scalar::Float(v) => v as i32, + } + } +} + +/// Decode one element of a numeric datatype. +fn decode_scalar( + bytes: &[u8], + dt: &Datatype, + order: &DatatypeByteOrder, +) -> Result { + match dt { + Datatype::FixedPoint { + size, + signed, + bit_offset, + bit_precision, + .. + } => { + let full = read_unsigned_int(bytes, *size as usize, order); + let (off, prec) = effective_bits(*size as usize, *bit_offset, *bit_precision); + Ok(if *signed { + Scalar::Signed(extract_signed(full, off, prec)) + } else { + Scalar::Unsigned(extract_unsigned(full, off, prec)) + }) + } + _ => convert_to_f64(bytes, dt, order).map(Scalar::Float), + } +} + /// Convert raw bytes to `i64` values. +/// +/// Floating-point data is converted the way libhdf5 converts it: truncated +/// toward zero, saturating at the target type's range, with NaN read as 0. pub fn read_as_i64(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { if let Datatype::Array { base_type, .. } = datatype { return read_as_i64(raw, base_type); @@ -1014,17 +1081,18 @@ pub fn read_as_i64(raw: &[u8], datatype: &Datatype) -> Result, FormatEr } let order = get_byte_order(datatype); - let (off, prec) = fixed_bits(datatype); let mut result = Vec::with_capacity(count); for i in 0..count { let chunk = &raw[i * elem_size..(i + 1) * elem_size]; - let full = read_unsigned_int(chunk, elem_size, &order); - result.push(extract_signed(full, off, prec)); + result.push(decode_scalar(chunk, datatype, &order)?.to_i64()); } Ok(result) } /// Convert raw bytes to `u64` values. +/// +/// Floating-point data is converted the way libhdf5 converts it: truncated +/// toward zero, saturating at the target type's range, with NaN read as 0. pub fn read_as_u64(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { if let Datatype::Array { base_type, .. } = datatype { return read_as_u64(raw, base_type); @@ -1039,12 +1107,10 @@ pub fn read_as_u64(raw: &[u8], datatype: &Datatype) -> Result, FormatEr } let count = raw.len() / elem_size; let order = get_byte_order(datatype); - let (off, prec) = fixed_bits(datatype); let mut result = Vec::with_capacity(count); for i in 0..count { let chunk = &raw[i * elem_size..(i + 1) * elem_size]; - let full = read_unsigned_int(chunk, elem_size, &order); - result.push(extract_unsigned(full, off, prec)); + result.push(decode_scalar(chunk, datatype, &order)?.to_u64()); } Ok(result) } @@ -1140,6 +1206,9 @@ pub fn read_as_f32(raw: &[u8], datatype: &Datatype) -> Result, FormatEr } /// Convert raw bytes to `i32` values. +/// +/// Floating-point data is converted the way libhdf5 converts it: truncated +/// toward zero, saturating at the target type's range, with NaN read as 0. pub fn read_as_i32(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { if let Datatype::Array { base_type, .. } = datatype { return read_as_i32(raw, base_type); @@ -1170,12 +1239,10 @@ pub fn read_as_i32(raw: &[u8], datatype: &Datatype) -> Result, FormatEr } let order = get_byte_order(datatype); - let (off, prec) = fixed_bits(datatype); let mut result = Vec::with_capacity(count); for i in 0..count { let chunk = &raw[i * elem_size..(i + 1) * elem_size]; - let full = read_unsigned_int(chunk, elem_size, &order); - result.push(extract_signed(full, off, prec) as i32); + result.push(decode_scalar(chunk, datatype, &order)?.to_i32()); } Ok(result) } @@ -1666,20 +1733,6 @@ fn effective_bits(size: usize, bit_offset: u16, bit_precision: u16) -> (u32, u32 (bit_offset as u32, prec) } -/// `(bit_offset, bit_precision)` for a fixed-point datatype, full width for -/// other types. -fn fixed_bits(datatype: &Datatype) -> (u32, u32) { - match datatype { - Datatype::FixedPoint { - size, - bit_offset, - bit_precision, - .. - } => effective_bits(*size as usize, *bit_offset, *bit_precision), - _ => (0, 0), - } -} - /// Whether a datatype occupies its full storage width (bit offset 0, precision /// == size·8), in which case the bulk-copy fast read paths apply. Non /// fixed-point types are treated as full width. @@ -1892,6 +1945,27 @@ mod tests { assert_eq!(read_as_u64(&raw, &dt).unwrap(), vec![4095, 1, 2048]); } + #[test] + fn float_to_int_truncates_and_saturates() { + // Values libhdf5 hands to an undefined C cast: NaN reads as 0 and + // exactly 2^63 saturates instead of wrapping to i64::MIN. + let dt = make_f64_le_type(); + let vals = [f64::NAN, 2f64.powi(63), -2.5, 2.0f64.powi(64)]; + let raw: Vec = vals.iter().flat_map(|v| v.to_le_bytes()).collect(); + assert_eq!( + read_as_i64(&raw, &dt).unwrap(), + vec![0, i64::MAX, -2, i64::MAX] + ); + assert_eq!( + read_as_u64(&raw, &dt).unwrap(), + vec![0, 1 << 63, 0, u64::MAX] + ); + assert_eq!( + read_as_i32(&raw, &dt).unwrap(), + vec![0, i32::MAX, -2, i32::MAX] + ); + } + #[test] fn full_width_signed_unchanged() { // Regression: full-width 32-bit signed must be unaffected. diff --git a/crates/clawhdf5/tests/numeric_conversion_interop.rs b/crates/clawhdf5/tests/numeric_conversion_interop.rs new file mode 100644 index 0000000..5c0fc3c --- /dev/null +++ b/crates/clawhdf5/tests/numeric_conversion_interop.rs @@ -0,0 +1,129 @@ +//! Numeric conversions on read, checked against h5py/libhdf5. +//! +//! h5py writes each file and prints what libhdf5 converts the data to +//! (`Dataset.astype`); the typed readers must return the same values. +//! Skipped when python3 with h5py is unavailable, unless +//! `CLAWHDF5_REQUIRE_INTEROP=1`. + +use std::collections::HashMap; +use std::process::Command; + +use clawhdf5::File; + +fn python() -> String { + std::env::var("CLAWHDF5_PYTHON").unwrap_or_else(|_| "python3".to_string()) +} + +fn interop_required() -> bool { + std::env::var("CLAWHDF5_REQUIRE_INTEROP").is_ok_and(|v| v == "1") +} + +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; + } + }; +} + +/// Run `script` (which writes the file at `path`) and return its stdout as +/// `key -> values`, one `key v1 v2 ...` line per key. +fn run_python(script: &str) -> HashMap> { + let output = Command::new(python()) + .args(["-c", script]) + .output() + .expect("failed to run python"); + assert!( + output.status.success(), + "python failed:\n{}", + String::from_utf8_lossy(&output.stderr) + ); + String::from_utf8_lossy(&output.stdout) + .lines() + .filter_map(|line| { + let mut words = line.split_whitespace().map(str::to_string); + Some((words.next()?, words.collect())) + }) + .collect() +} + +fn parse(values: &[String]) -> Vec +where + T::Err: std::fmt::Debug, +{ + values.iter().map(|v| v.parse().unwrap()).collect() +} + +/// Python prelude: `emit(key, array)` prints one line of integers. +const PRELUDE: &str = r#" +import h5py, numpy as np +def emit(key, arr): + print(key, *[int(v) for v in np.asarray(arr).ravel()]) +"#; + +#[test] +fn float_dataset_read_as_integers_converts_like_libhdf5() { + // read_i32/read_i64/read_u64 on a float dataset used to return the raw + // IEEE bit patterns (1.5 read as i64 was 4609434218613702656). + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("float_to_int.h5"); + let script = format!( + r#"{PRELUDE} +vals = [1.5, -2.75, 3e9, 1e300, -1e300, -0.5, 0.0, 7.99, np.inf, -np.inf, 1e19, -1e19] +with h5py.File("{path}", "w") as f: + for name, dt in (("f8", "f8"), ("f4", "(&expected[&format!("{name}:i32")]), + "{name} as i32" + ); + assert_eq!( + ds.read_i64().unwrap(), + parse::(&expected[&format!("{name}:i64")]), + "{name} as i64" + ); + if name != "f2" { + assert_eq!( + ds.read_u64().unwrap(), + parse::(&expected[&format!("{name}:u64")]), + "{name} as u64" + ); + } + } + // Negative values saturate at 0 rather than wrapping. + let f2 = file.dataset("f2").unwrap(); + assert_eq!(f2.read_u64().unwrap(), vec![1, 0, 0, 0, 7, 65504, 0]); +}