diff --git a/crates/clawhdf5-format/src/data_read.rs b/crates/clawhdf5-format/src/data_read.rs index e3c1bbd..82e6544 100644 --- a/crates/clawhdf5-format/src/data_read.rs +++ b/crates/clawhdf5-format/src/data_read.rs @@ -773,14 +773,7 @@ pub fn read_as_f64_zerocopy<'a>(raw: &'a [u8], datatype: &Datatype) -> Option<&' // Only native LE f64 is eligible #[cfg(target_endian = "little")] { - if !matches!( - datatype, - Datatype::FloatingPoint { - size: 8, - byte_order: DatatypeByteOrder::LittleEndian, - .. - } - ) { + if !is_native_le_float(datatype, FloatFormat::Double) { return None; } if !raw.len().is_multiple_of(8) { @@ -809,14 +802,7 @@ pub fn read_as_f64_zerocopy<'a>(raw: &'a [u8], datatype: &Datatype) -> Option<&' pub fn read_as_f32_zerocopy<'a>(raw: &'a [u8], datatype: &Datatype) -> Option<&'a [f32]> { #[cfg(target_endian = "little")] { - if !matches!( - datatype, - Datatype::FloatingPoint { - size: 4, - byte_order: DatatypeByteOrder::LittleEndian, - .. - } - ) { + if !is_native_le_float(datatype, FloatFormat::Single) { return None; } if !raw.len().is_multiple_of(4) { @@ -902,9 +888,9 @@ fn native_le_to_vec(raw: &[u8], count: usize) -> Vec { /// Convert raw bytes to `f64` values. pub fn read_as_f64(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { - // Array datatypes (e.g. an array-typed compound member) are read as a flat - // sequence of their base elements. - if let Datatype::Array { base_type, .. } = datatype { + // Array datatypes read as a flat sequence of their base elements, and + // enumerations (h5py's bool among them) as their integer values. + if let Datatype::Array { base_type, .. } | Datatype::Enumeration { base_type, .. } = datatype { return read_as_f64(raw, base_type); } ensure_numeric(datatype, "FloatingPoint or FixedPoint")?; @@ -919,20 +905,19 @@ pub fn read_as_f64(raw: &[u8], datatype: &Datatype) -> Result, FormatEr // Fast path: native-endian f64 — single bulk memcpy #[cfg(target_endian = "little")] - if matches!( - datatype, - Datatype::FloatingPoint { - size: 8, - byte_order: DatatypeByteOrder::LittleEndian, - .. - } - ) { + if is_native_le_float(datatype, FloatFormat::Double) { return Ok(native_le_to_vec::(raw, count)); } let order = get_byte_order(datatype); let mut result = Vec::with_capacity(count); - + if let Datatype::FloatingPoint { .. } = datatype { + let format = FloatFormat::of(datatype)?; + for chunk in raw.chunks_exact(elem_size) { + result.push(format.decode(chunk, &order)); + } + return Ok(result); + } for i in 0..count { let chunk = &raw[i * elem_size..(i + 1) * elem_size]; let val = convert_to_f64(chunk, datatype, &order)?; @@ -947,18 +932,7 @@ fn convert_to_f64( order: &DatatypeByteOrder, ) -> Result { match dt { - Datatype::FloatingPoint { size, .. } => match size { - 4 => { - let v = read_f32_bytes(bytes, order); - Ok(v as f64) - } - 8 => Ok(read_f64_bytes(bytes, order)), - 2 => Ok(read_f16_bytes(bytes, order) as f64), - _ => Err(FormatError::DataSizeMismatch { - expected: 8, - actual: *size as usize, - }), - }, + Datatype::FloatingPoint { .. } => Ok(FloatFormat::of(dt)?.decode(bytes, order)), Datatype::FixedPoint { size, signed, @@ -982,9 +956,83 @@ 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 { + // Every conversion follows libhdf5's default (hard) conversions: a value + // outside the target type's range saturates to its minimum or maximum — + // including a negative value read as unsigned, which reads as 0 — rather + // than being truncated to its low bits. Floats truncate toward zero; 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) => i64::try_from(v).unwrap_or(i64::MAX), + Scalar::Float(v) => v as i64, + } + } + + fn to_u64(self) -> u64 { + match self { + Scalar::Signed(v) => u64::try_from(v).unwrap_or(0), + Scalar::Unsigned(v) => v, + Scalar::Float(v) => v as u64, + } + } + + fn to_i32(self) -> i32 { + match self { + Scalar::Signed(v) => v.clamp(i32::MIN.into(), i32::MAX.into()) as i32, + Scalar::Unsigned(v) => i32::try_from(v).unwrap_or(i32::MAX), + 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. +/// +/// Values are converted the way libhdf5 converts them: integers outside the +/// target range saturate at its minimum or maximum (a negative value read as +/// unsigned is 0), and floating-point data is truncated toward zero and +/// saturated, with NaN read as 0. pub fn read_as_i64(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { - if let Datatype::Array { base_type, .. } = datatype { + // Array datatypes read as a flat sequence of their base elements, and + // enumerations (h5py's bool among them) as their integer values. + if let Datatype::Array { base_type, .. } | Datatype::Enumeration { base_type, .. } = datatype { return read_as_i64(raw, base_type); } ensure_numeric(datatype, "FixedPoint (signed)")?; @@ -1014,19 +1062,24 @@ 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. +/// +/// Values are converted the way libhdf5 converts them: integers outside the +/// target range saturate at its minimum or maximum (a negative value read as +/// unsigned is 0), and floating-point data is truncated toward zero and +/// saturated, with NaN read as 0. pub fn read_as_u64(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { - if let Datatype::Array { base_type, .. } = datatype { + // Array datatypes read as a flat sequence of their base elements, and + // enumerations (h5py's bool among them) as their integer values. + if let Datatype::Array { base_type, .. } | Datatype::Enumeration { base_type, .. } = datatype { return read_as_u64(raw, base_type); } ensure_numeric(datatype, "FixedPoint (unsigned)")?; @@ -1039,19 +1092,19 @@ 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) } /// Convert raw bytes to `f32` values. pub fn read_as_f32(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { - if let Datatype::Array { base_type, .. } = datatype { + // Array datatypes read as a flat sequence of their base elements, and + // enumerations (h5py's bool among them) as their integer values. + if let Datatype::Array { base_type, .. } | Datatype::Enumeration { base_type, .. } = datatype { return read_as_f32(raw, base_type); } ensure_numeric(datatype, "FloatingPoint")?; @@ -1066,25 +1119,11 @@ pub fn read_as_f32(raw: &[u8], datatype: &Datatype) -> Result, FormatEr // Fast path: native-endian f32 — single bulk memcpy #[cfg(target_endian = "little")] - if matches!( - datatype, - Datatype::FloatingPoint { - size: 4, - byte_order: DatatypeByteOrder::LittleEndian, - .. - } - ) { + if is_native_le_float(datatype, FloatFormat::Single) { return Ok(native_le_to_vec::(raw, count)); } - // Little-endian half precision (numpy float16): widen directly. - if matches!( - datatype, - Datatype::FloatingPoint { - size: 2, - byte_order: DatatypeByteOrder::LittleEndian, - .. - } - ) { + // Little-endian IEEE half precision (numpy float16): widen directly. + if is_native_le_float(datatype, FloatFormat::Half) { let (halves, _) = raw[..count * 2].as_chunks::<2>(); return Ok(halves .iter() @@ -1094,18 +1133,22 @@ pub fn read_as_f32(raw: &[u8], datatype: &Datatype) -> Result, FormatEr let order = get_byte_order(datatype); let mut result = Vec::with_capacity(count); + if let Datatype::FloatingPoint { .. } = datatype { + let format = FloatFormat::of(datatype)?; + for chunk in raw.chunks_exact(elem_size) { + result.push(match format { + FloatFormat::Single => read_f32_bytes(chunk, &order), + FloatFormat::Half => read_f16_bytes(chunk, &order), + // Double rounds; every other supported layout (bfloat16, FP8) + // is exact in f32. + _ => format.decode(chunk, &order) as f32, + }); + } + return Ok(result); + } for i in 0..count { let chunk = &raw[i * elem_size..(i + 1) * elem_size]; match datatype { - Datatype::FloatingPoint { size: 4, .. } => { - result.push(read_f32_bytes(chunk, &order)); - } - Datatype::FloatingPoint { size: 8, .. } => { - result.push(read_f64_bytes(chunk, &order) as f32); - } - Datatype::FloatingPoint { size: 2, .. } => { - result.push(read_f16_bytes(chunk, &order)); - } Datatype::FixedPoint { signed: true, size, @@ -1140,8 +1183,15 @@ pub fn read_as_f32(raw: &[u8], datatype: &Datatype) -> Result, FormatEr } /// Convert raw bytes to `i32` values. +/// +/// Values are converted the way libhdf5 converts them: integers outside the +/// target range saturate at its minimum or maximum (a negative value read as +/// unsigned is 0), and floating-point data is truncated toward zero and +/// saturated, with NaN read as 0. pub fn read_as_i32(raw: &[u8], datatype: &Datatype) -> Result, FormatError> { - if let Datatype::Array { base_type, .. } = datatype { + // Array datatypes read as a flat sequence of their base elements, and + // enumerations (h5py's bool among them) as their integer values. + if let Datatype::Array { base_type, .. } | Datatype::Enumeration { base_type, .. } = datatype { return read_as_i32(raw, base_type); } ensure_numeric(datatype, "FixedPoint")?; @@ -1162,6 +1212,7 @@ pub fn read_as_i32(raw: &[u8], datatype: &Datatype) -> Result, FormatEr datatype, Datatype::FixedPoint { byte_order: DatatypeByteOrder::LittleEndian, + signed: true, .. } ) @@ -1170,12 +1221,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) } @@ -1616,6 +1665,174 @@ fn reorder_bytes(bytes: &[u8], order: &DatatypeByteOrder) -> [u8; 8] { buf } +/// How the bits of a floating-point datatype are laid out, read from the +/// datatype message's fields rather than assumed from its size (a 2-byte +/// float may be IEEE half or bfloat16). +#[derive(Debug, Clone, Copy, PartialEq)] +enum FloatFormat { + /// IEEE-754 binary16. + Half, + /// IEEE-754 binary32. + Single, + /// IEEE-754 binary64. + Double, + /// Any other IEEE-style layout (implied leading mantissa bit, all-ones + /// exponent for infinity/NaN) whose values are all exact in `f64`: + /// bfloat16, the FP8 formats, and similar. + Other(FloatLayout), +} + +#[derive(Debug, Clone, Copy, PartialEq)] +struct FloatLayout { + exponent_location: u32, + exponent_size: u32, + mantissa_location: u32, + mantissa_size: u32, + exponent_bias: u32, +} + +impl FloatFormat { + fn of(dt: &Datatype) -> Result { + let Datatype::FloatingPoint { + size, + exponent_location, + exponent_size, + mantissa_location, + mantissa_size, + exponent_bias, + .. + } = dt + else { + return Err(FormatError::TypeMismatch { + expected: "FloatingPoint", + actual: datatype_name(dt), + }); + }; + let layout = FloatLayout { + exponent_location: u32::from(*exponent_location), + exponent_size: u32::from(*exponent_size), + mantissa_location: u32::from(*mantissa_location), + mantissa_size: u32::from(*mantissa_size), + exponent_bias: *exponent_bias, + }; + let fields = ( + layout.exponent_location, + layout.exponent_size, + layout.mantissa_location, + layout.mantissa_size, + layout.exponent_bias, + ); + let bits = size.saturating_mul(8); + // The sign bit is not kept in `Datatype`; every standard layout has it + // directly above the exponent, with the mantissa below. + let well_formed = layout.exponent_size > 0 + && layout.mantissa_size > 0 + && layout.mantissa_location + layout.mantissa_size <= layout.exponent_location + && layout.exponent_location + layout.exponent_size < bits; + match (size, fields) { + (2, (10, 5, 0, 10, 15)) => Ok(FloatFormat::Half), + (4, (23, 8, 0, 23, 127)) => Ok(FloatFormat::Single), + (8, (52, 11, 0, 52, 1023)) => Ok(FloatFormat::Double), + _ if well_formed + && *size <= 8 + && layout.exponent_size <= 11 + && layout.mantissa_size <= 52 => + { + Ok(FloatFormat::Other(layout)) + } + // Fields that cannot describe any float (e.g. left zeroed by a + // hand-built datatype): fall back to the IEEE type of that size. + (2, _) if !well_formed => Ok(FloatFormat::Half), + (4, _) if !well_formed => Ok(FloatFormat::Single), + (8, _) if !well_formed => Ok(FloatFormat::Double), + // x87 80-bit extended, binary128, ...: not representable in f64. + _ => Err(FormatError::TypeMismatch { + expected: "floating point of at most 64 bits (IEEE-style layout)", + actual: "FloatingPoint", + }), + } + } + + fn decode(self, bytes: &[u8], order: &DatatypeByteOrder) -> f64 { + match self { + FloatFormat::Half => f64::from(read_f16_bytes(bytes, order)), + FloatFormat::Single => f64::from(read_f32_bytes(bytes, order)), + FloatFormat::Double => read_f64_bytes(bytes, order), + FloatFormat::Other(layout) => { + layout.decode(read_unsigned_int(bytes, bytes.len(), order)) + } + } + } +} + +impl FloatLayout { + /// Decode the value held in the low `size * 8` bits of `bits`. + fn decode(self, bits: u64) -> f64 { + let field = |location: u32, size: u32| (bits >> location) & ((1u64 << size) - 1); + let exponent = field(self.exponent_location, self.exponent_size); + let mantissa = field(self.mantissa_location, self.mantissa_size); + let negative = field(self.exponent_location + self.exponent_size, 1) == 1; + let max_exponent = (1u64 << self.exponent_size) - 1; + let magnitude = if exponent == max_exponent { + if mantissa == 0 { + f64::INFINITY + } else { + f64::NAN + } + } else { + let bias = i64::from(self.exponent_bias); + let msize = i64::from(self.mantissa_size); + // value = significand * 2^power, with an implied leading 1 unless + // the number is subnormal (exponent field 0). + let (significand, power) = if exponent == 0 { + (mantissa, 1 - bias - msize) + } else { + ( + mantissa | (1u64 << self.mantissa_size), + exponent as i64 - bias - msize, + ) + }; + scale_by_pow2(significand as f64, power) + }; + if negative { -magnitude } else { magnitude } + } +} + +/// `x * 2^power` without `std` (no `powi`/`libm`). `x` is a non-negative +/// integer below 2^53, so it is exact. +fn scale_by_pow2(x: f64, power: i64) -> f64 { + if x == 0.0 || power < -1200 { + return 0.0; + } + if power > 1100 { + return f64::INFINITY; + } + let pow2 = |p: i64| f64::from_bits(((p + 1023) as u64) << 52); + let mut x = x; + let mut power = power; + while power > 1023 { + x *= pow2(1023); + power -= 1023; + } + while power < -1022 { + x *= pow2(-1022); + power += 1022; + } + x * pow2(power) +} + +/// Whether `datatype` is the little-endian IEEE float `format`, whose bytes +/// can be copied straight into native values on a little-endian target. +fn is_native_le_float(datatype: &Datatype, format: FloatFormat) -> bool { + matches!( + datatype, + Datatype::FloatingPoint { + byte_order: DatatypeByteOrder::LittleEndian, + .. + } + ) && FloatFormat::of(datatype).is_ok_and(|f| f == format) +} + fn read_f64_bytes(bytes: &[u8], order: &DatatypeByteOrder) -> f64 { let buf = reorder_bytes(bytes, order); f64::from_le_bytes(buf) @@ -1666,20 +1883,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 +2095,67 @@ 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 bfloat16_and_fp8_decode_by_fields() { + // bfloat16 is a 2-byte float that is not IEEE half. + let bf16 = Datatype::FloatingPoint { + size: 2, + byte_order: DatatypeByteOrder::LittleEndian, + bit_offset: 0, + bit_precision: 16, + exponent_location: 7, + exponent_size: 8, + mantissa_location: 0, + mantissa_size: 7, + exponent_bias: 127, + }; + let raw: Vec = [0x3FC0u16, 0xC010, 0x7F80, 0x0001] + .iter() + .flat_map(|v| v.to_le_bytes()) + .collect(); + let got = read_as_f64(&raw, &bf16).unwrap(); + assert_eq!(&got[..3], &[1.5, -2.25, f64::INFINITY]); + assert_eq!(got[3], 2f64.powi(-133)); // smallest subnormal + assert_eq!(read_as_f32(&raw, &bf16).unwrap()[..2], [1.5, -2.25]); + + // FP8 E4M3: 1, -1, 2, 0, NaN (IEEE-style, as libhdf5 treats it). + let e4m3 = Datatype::FloatingPoint { + size: 1, + byte_order: DatatypeByteOrder::LittleEndian, + bit_offset: 0, + bit_precision: 8, + exponent_location: 3, + exponent_size: 4, + mantissa_location: 0, + mantissa_size: 3, + exponent_bias: 7, + }; + let got = read_as_f64(&[0x38, 0xB8, 0x40, 0x00, 0x7E], &e4m3).unwrap(); + assert_eq!(&got[..4], &[1.0, -1.0, 2.0, 0.0]); + assert!(got[4].is_nan()); + } + #[test] fn full_width_signed_unchanged() { // Regression: full-width 32-bit signed must be unaffected. diff --git a/crates/clawhdf5-format/src/vl_data.rs b/crates/clawhdf5-format/src/vl_data.rs index 6b4dcfc..9a50d51 100644 --- a/crates/clawhdf5-format/src/vl_data.rs +++ b/crates/clawhdf5-format/src/vl_data.rs @@ -148,7 +148,12 @@ pub fn read_vl_strings( Ok(result) } -/// Resolve VL byte sequences from raw data. +/// Resolve VL sequences from raw data, returning each element's bytes. +/// +/// Each element is the sequence's full encoding — element count × base type +/// size bytes, in the base type's byte order — so a sequence of `i32` yields +/// four bytes per value. Decode it with the base type (e.g. +/// [`crate::data_read::read_as_i64`]). pub fn read_vl_bytes( file_data: &[u8], raw_data: &[u8], @@ -177,8 +182,10 @@ pub fn read_vl_bytes( }, )?; - let len = (vl.length as usize).min(obj.data.len()); - result.push(obj.data[..len].to_vec()); + // The heap object holds the whole sequence. `vl.length` counts + // elements, not bytes, so it is only the byte length when the base + // type is one byte wide. + result.push(obj.data.clone()); } Ok(result) diff --git a/crates/clawhdf5/tests/numeric_conversion_interop.rs b/crates/clawhdf5/tests/numeric_conversion_interop.rs new file mode 100644 index 0000000..357a37d --- /dev/null +++ b/crates/clawhdf5/tests/numeric_conversion_interop.rs @@ -0,0 +1,374 @@ +//! 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]); +} + +#[test] +fn integer_reads_saturate_out_of_range_values_like_libhdf5() { + // Narrowing reads used to keep the low bits (i64 2^40+5 read as i32 was + // 5, u64::MAX read as i64 was -1) and signed-to-unsigned reads wrapped + // (-1 read as u64 was 4294967295). + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("int_narrowing.h5"); + let script = format!( + r#"{PRELUDE} +data = {{ + "i8": np.array([2**40 + 5, -(2**35), 7, -1, 2**63 - 1, -(2**63)], "i8"), + "u8": np.array([2**64 - 1, 2**63, 5, 0], "i2"), + "u1": np.array([255, 0, 128], "u1"), +}} +with h5py.File("{path}", "w") as f: + for name, arr in data.items(): + f.create_dataset(name, data=arr) +with h5py.File("{path}", "r") as f: + for name in data: + d = f[name] + emit(name + ":i32", d.astype("(&expected[&format!("{name}:i32")]), + "{name} as i32" + ); + assert_eq!( + ds.read_i64().unwrap(), + parse::(&expected[&format!("{name}:i64")]), + "{name} as i64" + ); + // libhdf5 wraps a negative big-endian i64 read as little-endian u64 + // (it only byte-swaps when the sizes match and the order differs); + // every other signed-to-unsigned read saturates at 0, so do that. + if name != "i8be" { + assert_eq!( + ds.read_u64().unwrap(), + parse::(&expected[&format!("{name}:u64")]), + "{name} as u64" + ); + } + } + let i8be = file.dataset("i8be").unwrap(); + assert_eq!(i8be.read_u64().unwrap(), vec![(1 << 40) + 5, 0, 7, 0]); +} + +#[test] +fn floats_decode_by_their_datatype_fields() { + // Every 2-byte float used to decode as IEEE half, so bfloat16 1.5 read as + // 1.9375 and +inf as NaN; 1-byte FP8 floats were refused. + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("float_layouts.h5"); + let script = format!( + r#"{PRELUDE} +def custom(base, fields, bias, size): + # fields: (sign pos, exponent pos, exponent size, mantissa pos, mantissa size) + t = base.copy() + t.set_fields(*fields) + t.set_ebias(bias) + t.set_precision(size * 8) + t.set_size(size) + return t + +def write(f, name, ftype, raw): + raw = np.ascontiguousarray(raw) + space = h5py.h5s.create_simple(raw.shape) + ds = h5py.h5d.create(f.id, name.encode(), ftype, space) + ds.write(h5py.h5s.ALL, h5py.h5s.ALL, raw, mtype=ftype) + +# bfloat16: 1.5, -2.25, +inf, 0, 3.140625, 1, -0, smallest subnormal, +# largest finite, NaN +bf16 = np.array([0x3FC0, 0xC010, 0x7F80, 0x0000, 0x4049, 0x3F80, 0x8000, 0x0001, + 0x7F7F, 0x7FC1], "f2")) + f.create_dataset("f4_be", data=vals.astype(">f4")) + f.create_dataset("f8_le", data=vals.astype(" = parse(&expected[name]); + let got: Vec = ds + .read_f64() + .unwrap() + .into_iter() + .map(|v| canonical(v).to_bits()) + .collect(); + assert_eq!(got, want, "{name} as f64"); + // f32 reads agree too (every value here is exact in f32 except the + // f64 dataset, which rounds like `as f32`). + let got32: Vec = ds + .read_f32() + .unwrap() + .into_iter() + .map(|v| if v.is_nan() { f32::NAN } else { v }.to_bits()) + .collect(); + let want32: Vec = want + .iter() + .map(|&b| canonical(f64::from_bits(b)) as f32) + .map(|v| if v.is_nan() { f32::NAN } else { v }.to_bits()) + .collect(); + assert_eq!(got32, want32, "{name} as f32"); + } +} + +#[test] +fn enum_and_bool_datasets_read_as_their_integer_values() { + // Enumerations (h5py stores bool as an enum of int8) were refused by the + // numeric readers with a type mismatch. + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("enums.h5"); + let script = format!( + r#"{PRELUDE} +with h5py.File("{path}", "w") as f: + f.create_dataset("bool", data=np.array([True, False, True])) + e = h5py.enum_dtype({{"RED": 0, "GREEN": 7, "BLUE": -3}}, basetype=">i2") + f.create_dataset("enum_i2be", data=np.array([0, 7, -3, 7], ">i2"), dtype=e) + e = h5py.enum_dtype({{"LOW": 0, "HIGH": 200}}, basetype="u1") + f.create_dataset("enum_u1", data=np.array([200, 0, 200], "u1"), dtype=e) + e = h5py.enum_dtype({{"A": -(2**40), "B": 2**40}}, basetype=" = parse(&expected[name]); + assert_eq!(ds.read_i64().unwrap(), want, "{name} as i64"); + let want_f64: Vec = want.iter().map(|&v| v as f64).collect(); + assert_eq!(ds.read_f64().unwrap(), want_f64, "{name} as f64"); + } + let bools = file.dataset("bool").unwrap(); + assert_eq!(bools.read_u64().unwrap(), vec![1, 0, 1]); + assert_eq!(bools.read_i32().unwrap(), vec![1, 0, 1]); +} + +#[test] +fn vl_sequences_of_wide_base_types_read_whole() { + // read_vl_bytes took the sequence's element count as its byte length, so + // [1, 2, 3] as VL int32 came back as 3 bytes instead of 12. + use clawhdf5::Selection; + use clawhdf5_format::datatype::{Datatype, DatatypeByteOrder}; + use clawhdf5_format::vl_data::read_vl_bytes; + + skip_if_no_python!(); + let dir = tempfile::tempdir().unwrap(); + let path = dir.path().join("vlen.h5"); + let script = format!( + r#"{PRELUDE} +data = {{ + "i4": (h5py.vlen_dtype("> = (0..count) + .map(|i| parse(expected.get(&format!("{name}:{i}")).unwrap())) + .collect(); + assert_eq!(got, want, "{name}"); + if name == "i4" { + let i32_le = Datatype::FixedPoint { + size: 4, + byte_order: DatatypeByteOrder::LittleEndian, + signed: true, + bit_offset: 0, + bit_precision: 32, + }; + let values = clawhdf5_format::data_read::read_as_i64(&got[0], &i32_le).unwrap(); + assert_eq!(values, vec![1, 2, 3]); + assert_eq!(got[3].len(), 40 * 4); + } + } +}