gwseq_io/bbi/

header.rs

//! bbi headers: the common header, the zoom ladder, the total summary, autoSql.
//!
//! Field offsets are Supp. Tables
//! 5–7 of the bigWig paper; they are spelled out here rather than derived from
//! a struct layout, because the format's fields are not aligned and a
//! `#[repr(C)]` view of them would be undefined behaviour on half of them.

use crate::error::{Error, Result};
use crate::source::ByteSource;

/// Byte sizes of the fixed-width blocks (Supp. Tables 4–8, 14).
pub const BBI_HEADER_SIZE: u64 = 64;
pub const ZOOM_HEADER_SIZE: u64 = 24;
pub const TOTAL_SUMMARY_SIZE: u64 = 40;
pub const CHR_TREE_HEADER_SIZE: u64 = 32;
/// The R-tree index header the tree nodes follow.
pub const DATA_TREE_HEADER_SIZE: u64 = 48;

pub const BIGWIG_MAGIC_SWAPPED: u32 = 0x26FC_8F88;
pub const BIGBED_MAGIC_SWAPPED: u32 = 0xEBF2_8987;
pub const CHR_TREE_MAGIC: u32 = 0x78CA_8C91;
pub const CHR_TREE_MAGIC_SWAPPED: u32 = 0x91CA_8C78;
pub const DATA_TREE_MAGIC: u32 = 0x2468_ACE0;
pub const DATA_TREE_MAGIC_SWAPPED: u32 = 0xE0AC_6824;

pub const BBI_MIN_VERSION: u16 = 3;
pub const BBI_OUTPUT_VERSION: u16 = 4;

#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum BbiKind {
    BigWig,
    BigBed,
}

impl BbiKind {
    pub fn as_str(self) -> &'static str {
        match self {
            BbiKind::BigWig => "bigwig",
            BbiKind::BigBed => "bigbed",
        }
    }

    pub fn is_bigbed(self) -> bool {
        matches!(self, BbiKind::BigBed)
    }
}

#[derive(Debug, Clone)]
pub struct BbiHeader {
    pub kind: BbiKind,
    pub version: u16,
    pub zoom_levels: u16,
    pub chr_tree_offset: u64,
    pub full_data_offset: u64,
    pub full_index_offset: u64,
    /// bigBed only.
    pub field_count: u16,
    pub defined_field_count: u16,
    pub auto_sql_offset: u64,
    pub total_summary_offset: u64,
    /// 0 when the file stores its blocks uncompressed.
    pub uncompress_buffer_size: u32,
}

#[derive(Debug, Clone, Copy)]
pub struct ZoomHeader {
    pub reduction_level: u32,
    pub data_offset: u64,
    pub index_offset: u64,
}

/// The summary block is optional: a file that carries none sets its offset to
/// 0, and the reader leaves this as it stands rather than parsing 40 bytes of
/// whatever happens to sit at the start of the file. Hence NaN extremes.
#[derive(Debug, Clone, Copy)]
pub struct TotalSummary {
    pub bases_covered: u64,
    pub min_value: f64,
    pub max_value: f64,
    pub sum_data: f64,
    pub sum_squared: f64,
}

impl Default for TotalSummary {
    fn default() -> Self {
        Self {
            bases_covered: 0,
            min_value: f64::NAN,
            max_value: f64::NAN,
            sum_data: 0.0,
            sum_squared: 0.0,
        }
    }
}

#[derive(Debug, Clone, Copy)]
pub struct ChrTreeHeader {
    pub block_size: u32,
    pub key_size: u32,
    pub val_size: u32,
    pub item_count: u64,
}

/// Read the common header.
///
/// A byte-swapped magic is recognised and **refused**: these readers do not
/// swap, and saying so is more useful than reporting the file as unrecognisable.
pub fn read_header(source: &dyn ByteSource) -> Result<BbiHeader> {
    let buf = source.read_exact_at(0, BBI_HEADER_SIZE as usize)?;
    let path = source.path();
    let mut c = crate::bytes::LeCursor::new(&buf, 0, path);

    let magic = c.read_u32()?;
    let kind = match magic {
        super::BIGWIG_MAGIC => BbiKind::BigWig,
        super::BIGBED_MAGIC => BbiKind::BigBed,
        BIGWIG_MAGIC_SWAPPED | BIGBED_MAGIC_SWAPPED => {
            return Err(Error::format(path, "incompatible endianness"))
        }
        _ => return Err(Error::format(path, "not a bigwig or bigbed file")),
    };

    let version = c.read_u16()?;
    if version < BBI_MIN_VERSION {
        return Err(Error::format(
            path,
            format!("bigwig or bigbed version {version} unsupported (>= {BBI_MIN_VERSION})"),
        ));
    }

    Ok(BbiHeader {
        kind,
        version,
        zoom_levels: c.read_u16()?,
        chr_tree_offset: c.read_u64()?,
        full_data_offset: c.read_u64()?,
        full_index_offset: c.read_u64()?,
        field_count: c.read_u16()?,
        defined_field_count: c.read_u16()?,
        auto_sql_offset: c.read_u64()?,
        total_summary_offset: c.read_u64()?,
        uncompress_buffer_size: c.read_u32()?,
        // Eight reserved bytes follow at offset 56. UCSC reads them as an
        // extension offset on version 4 files and follows the pointer when they
        // are not zero, so they cannot be repurposed.
    })
}

pub fn read_zoom_headers(source: &dyn ByteSource, count: u16) -> Result<Vec<ZoomHeader>> {
    if count == 0 {
        return Ok(Vec::new());
    }
    let len = count as u64 * ZOOM_HEADER_SIZE;
    let buf = source.read_exact_at(BBI_HEADER_SIZE, len as usize)?;
    let mut c = crate::bytes::LeCursor::new(&buf, BBI_HEADER_SIZE, source.path());
    let mut headers = Vec::with_capacity(count as usize);
    for _ in 0..count {
        let reduction_level = c.read_u32()?;
        c.skip(4)?; // reserved
        headers.push(ZoomHeader {
            reduction_level,
            data_offset: c.read_u64()?,
            index_offset: c.read_u64()?,
        });
    }
    Ok(headers)
}

/// Read the total summary, or the default when the file declares none.
pub fn read_total_summary(source: &dyn ByteSource, offset: u64) -> Result<TotalSummary> {
    if offset == 0 {
        return Ok(TotalSummary::default());
    }
    let buf = source.read_exact_at(offset, TOTAL_SUMMARY_SIZE as usize)?;
    let mut c = crate::bytes::LeCursor::new(&buf, offset, source.path());
    Ok(TotalSummary {
        bases_covered: c.read_u64()?,
        min_value: c.read_f64()?,
        max_value: c.read_f64()?,
        sum_data: c.read_f64()?,
        sum_squared: c.read_f64()?,
    })
}

/// The R-tree index header that precedes the tree nodes.
///
/// Only the magic is checked here: everything else the walk needs, it reads off
/// the nodes themselves.
pub fn check_data_tree_magic(source: &dyn ByteSource, offset: u64) -> Result<()> {
    let buf = source.read_exact_at(offset, 4)?;
    let magic = u32::from_le_bytes([buf[0], buf[1], buf[2], buf[3]]);
    match magic {
        DATA_TREE_MAGIC => Ok(()),
        DATA_TREE_MAGIC_SWAPPED => Err(Error::format(
            source.path(),
            "incompatible endianness (data tree)",
        )),
        _ => Err(Error::format(
            source.path(),
            "invalid data tree magic number",
        )),
    }
}

/// The standard predefined BED field names, in order (Supp. Table 1), used to
/// name the columns of a file carrying no autoSql of its own.
pub const BED_FIELD_NAMES: &[&str] = &[
    "chrom",
    "chromStart",
    "chromEnd",
    "name",
    "score",
    "strand",
    "thickStart",
    "thickEnd",
    "itemRgb",
    "blockCount",
    "blockSizes",
    "blockStarts",
];

/// The bigBed column declaration. Field order matters — the first three are the
/// coordinates — so this is an `IndexMap`, not a `HashMap`.
///
/// An offset of 0 means the file carries no autoSql, and the standard BED names
/// stand in for as many columns as its header declares.
///
/// The text is one field per line, `type name1, name2;`. The parser takes
/// everything up to the first semicolon on a line as the field list, which is
/// why the writer must keep semicolons out of comments: one there would be read
/// as another field and the count would stop matching the header.
pub fn read_auto_sql(
    source: &dyn ByteSource,
    offset: u64,
    field_count: u16,
) -> Result<indexmap::IndexMap<String, String>> {
    if offset == 0 {
        return Ok(default_bed_fields(field_count as usize));
    }
    let text = read_nul_terminated(source, offset)?;

    let mut fields = indexmap::IndexMap::new();
    for line in text.lines() {
        // `<type> <names>;` — the first whitespace-delimited token is the type,
        // everything up to the first `;` is the comma-separated name list.
        let Some((head, rest)) = split_once_whitespace(line.trim()) else {
            continue;
        };
        let Some(list) = rest.split(';').next().filter(|l| !l.trim().is_empty()) else {
            continue;
        };
        if !rest.contains(';') {
            continue;
        }
        for name in list.split(',') {
            // The declaration may carry an array size — `int[blockCount] sizes`
            // — and a trailing comment; the name is the first token either way.
            let name = name.split_whitespace().next().unwrap_or("");
            if !name.is_empty() {
                fields.insert(name.to_string(), head.to_string());
            }
        }
    }

    validate_bed_field_names(&fields, source.path())?;
    if fields.len() != field_count as usize {
        return Err(Error::format(
            source.path(),
            format!(
                "field count {field_count} does not match autosql field count {}",
                fields.len()
            ),
        ));
    }
    Ok(fields)
}

/// Split on the whitespace that separates an autoSql field's type from its
/// names.
///
/// Not simply the first run of it: a real bigBed declares
/// `enum("+", "-", ".") strand;`, and `set(a, b) flags;`, whose types carry
/// spaces inside their parentheses. Whitespace inside a bracket or a quoted
/// string is part of the type, and the split is the first run outside both — so
/// the line above is the type `enum("+", "-", ".")` and the name `strand`,
/// rather than the type `enum("+",` and a field count one short.
fn split_once_whitespace(line: &str) -> Option<(&str, &str)> {
    let (mut depth, mut quote) = (0i32, None::<char>);
    let mut head_end = None;
    for (at, c) in line.char_indices() {
        match quote {
            Some(q) => {
                if c == q {
                    quote = None;
                }
            }
            None => match c {
                '"' | '\'' => quote = Some(c),
                '(' | '[' | '{' => depth += 1,
                ')' | ']' | '}' => depth -= 1,
                _ if depth <= 0 && c.is_whitespace() => {
                    head_end = Some(at);
                    break;
                }
                _ => {}
            },
        }
    }
    let (head, rest) = line.split_at(head_end?);
    Some((head, rest.trim_start()))
}

/// The autoSql text is NUL-terminated and its length is nowhere in the header,
/// so it is read in chunks until the terminator turns up.
fn read_nul_terminated(source: &dyn ByteSource, offset: u64) -> Result<String> {
    const CHUNK: usize = 4096;
    // An autoSql block describing a few hundred columns is a few kilobytes; a
    // file claiming megabytes of it is corrupt, not verbose.
    const MAX: usize = 1 << 20;
    let mut buf: Vec<u8> = Vec::new();
    loop {
        let chunk = source.read_at(offset + buf.len() as u64, CHUNK)?;
        if chunk.is_empty() {
            return Err(Error::corrupt(
                source.path(),
                offset,
                "autosql block is not NUL-terminated before the end of the file",
            ));
        }
        if let Some(at) = memchr::memchr(0, &chunk) {
            buf.extend_from_slice(&chunk[..at]);
            break;
        }
        buf.extend_from_slice(&chunk);
        if buf.len() > MAX {
            return Err(Error::corrupt(
                source.path(),
                offset,
                format!("autosql block is longer than {MAX} bytes"),
            ));
        }
    }
    String::from_utf8(buf)
        .map_err(|_| Error::corrupt(source.path(), offset, "autosql is not valid UTF-8"))
}

/// Check that the first three fields are the chromosome, the start and the end,
/// in that order.
///
/// Shared by the reader, which validates the autoSql it parses, and by the
/// writer, which validates the fields it was handed before turning them into
/// autoSql: a description this accepts on the way out must be one it accepts on
/// the way back in.
///
/// `chrom` may be spelled `chr`, since a BED written from this library names
/// its first column that, and the start and end may or may not carry the
/// `chrom` prefix. Hand-written matchers rather than regexes — three fixed
/// shapes do not need a regex engine.
pub fn validate_bed_field_names(
    fields: &indexmap::IndexMap<String, String>,
    path: &str,
) -> Result<()> {
    let bad = || Error::format(path, "missing or misplaced chr, start or end in autosql");
    if fields.len() < 3 {
        return Err(bad());
    }
    let name = |i: usize| fields.get_index(i).map(|(k, _)| k.as_str()).unwrap_or("");
    if !is_chr_name(name(0)) || !is_coord_name(name(1), "start") || !is_coord_name(name(2), "end") {
        return Err(bad());
    }
    Ok(())
}

/// `chr`, `chrom`, and either with an `id` or `name` suffix: `chr_id`,
/// `chromName`, `chrom_id`, …
fn is_chr_name(name: &str) -> bool {
    let lower = name.to_ascii_lowercase();
    let rest = match strip_chr_prefix(&lower) {
        Some(rest) => rest,
        None => return false,
    };
    let rest = rest.strip_prefix('_').unwrap_or(rest);
    matches!(rest, "" | "id" | "name")
}

/// `start`/`end`, with an optional `chr`/`chrom` prefix: `chromStart`,
/// `chrom_end`, `start`, …
fn is_coord_name(name: &str, suffix: &str) -> bool {
    let lower = name.to_ascii_lowercase();
    if lower == suffix {
        return true;
    }
    let Some(rest) = strip_chr_prefix(&lower) else {
        return false;
    };
    rest.strip_prefix('_').unwrap_or(rest) == suffix
}

fn strip_chr_prefix(lower: &str) -> Option<&str> {
    let rest = lower.strip_prefix("chr")?;
    Some(rest.strip_prefix("om").unwrap_or(rest))
}

/// The standard BED column names for a given width — `field13` and up past the
/// twelve the format names.
pub fn default_bed_fields(col_count: usize) -> indexmap::IndexMap<String, String> {
    (0..col_count)
        .map(|i| {
            let name = BED_FIELD_NAMES
                .get(i)
                .map(|n| (*n).to_string())
                .unwrap_or_else(|| format!("field{}", i + 1));
            (name, "string".to_string())
        })
        .collect()
}

/// Narrow a count or an offset to the 32 bits a bbi file stores it on.
///
/// Every coordinate, count and index the writer puts in a
/// file goes through this: the format is 32-bit throughout, and a value that
/// does not fit has to be refused where it is written rather than truncated
/// into a file that reads back as something else.
pub fn to_bbi_u32(value: i64, what: &str) -> Result<u32> {
    if !(0..=0xFFFF_FFFF).contains(&value) {
        return Err(Error::invalid(format!(
            "{what} {value} does not fit the 32 bits a bigwig or bigbed file stores it on"
        )));
    }
    Ok(value as u32)
}

/// Narrow an `f64` accumulator to the `f32` a zoom record stores it on,
/// clamping to the range `f32` can represent.
///
/// An `as` cast saturates to infinity, and a reader would carry that value
/// into every mean it computes. Clamping to the extremes keeps a summary
/// finite.
pub fn to_bbi_f32(value: f64) -> f32 {
    const LIMIT: f64 = f32::MAX as f64;
    if value > LIMIT {
        return f32::MAX;
    }
    if value < -LIMIT {
        return f32::MIN;
    }
    value as f32
}

/// The 64 bytes of the common header (Supp. Table 5), **with the magic left as
/// zero**.
///
/// The writer stamps the magic separately and last, so a write that dies partway
/// leaves a file no reader accepts rather than a header pointing at data that
/// was never written.
pub fn write_header(header: &BbiHeader) -> Result<Vec<u8>> {
    let mut out = Vec::with_capacity(BBI_HEADER_SIZE as usize);
    out.extend_from_slice(&0u32.to_le_bytes()); // magic, stamped separately
    out.extend_from_slice(&header.version.to_le_bytes());
    out.extend_from_slice(&header.zoom_levels.to_le_bytes());
    out.extend_from_slice(&header.chr_tree_offset.to_le_bytes());
    out.extend_from_slice(&header.full_data_offset.to_le_bytes());
    out.extend_from_slice(&header.full_index_offset.to_le_bytes());
    out.extend_from_slice(&header.field_count.to_le_bytes());
    out.extend_from_slice(&header.defined_field_count.to_le_bytes());
    out.extend_from_slice(&header.auto_sql_offset.to_le_bytes());
    out.extend_from_slice(&header.total_summary_offset.to_le_bytes());
    out.extend_from_slice(&header.uncompress_buffer_size.to_le_bytes());
    out.extend_from_slice(&[0u8; 8]); // reserved
    debug_assert_eq!(out.len(), BBI_HEADER_SIZE as usize);
    Ok(out)
}

/// The 24 bytes of one zoom header (Supp. Table 6).
pub fn write_zoom_header(header: &ZoomHeader) -> Vec<u8> {
    let mut out = Vec::with_capacity(ZOOM_HEADER_SIZE as usize);
    out.extend_from_slice(&header.reduction_level.to_le_bytes());
    out.extend_from_slice(&0u32.to_le_bytes()); // reserved
    out.extend_from_slice(&header.data_offset.to_le_bytes());
    out.extend_from_slice(&header.index_offset.to_le_bytes());
    out
}

/// The 40 bytes of the total summary (Supp. Table 7).
///
/// A file carrying no data writes zero extremes rather than the NaN the struct
/// defaults to: the block is present, so a reader parses it, and NaN extremes
/// would propagate into everything computed from them.
pub fn write_total_summary(summary: &TotalSummary) -> Vec<u8> {
    let empty = summary.bases_covered == 0;
    let mut out = Vec::with_capacity(TOTAL_SUMMARY_SIZE as usize);
    out.extend_from_slice(&summary.bases_covered.to_le_bytes());
    out.extend_from_slice(&if empty { 0.0 } else { summary.min_value }.to_le_bytes());
    out.extend_from_slice(&if empty { 0.0 } else { summary.max_value }.to_le_bytes());
    out.extend_from_slice(&summary.sum_data.to_le_bytes());
    out.extend_from_slice(&summary.sum_squared.to_le_bytes());
    out
}

/// The field types an autoSql declaration may name.
pub const BED_FIELD_TYPES: &[&str] = &["string", "int", "uint", "float"];

/// An autoSql declaration and the two counts the header states alongside it.
#[derive(Debug, Clone)]
pub struct BedAutoSql {
    pub text: String,
    pub field_count: u16,
    /// How many leading columns are the standard BED ones, which is what tells
    /// a browser a bigBed is a bed4, a bed6 or a bed12.
    pub defined_field_count: u16,
}

/// What the standard columns are for, as UCSC describes them (Supp. Table 2).
///
/// A browser shows these, so the ones it already knows are named the way it
/// expects and the rest fall back to the name the caller gave.
const BED_FIELD_COMMENTS: &[&str] = &[
    "Reference sequence chromosome or scaffold",
    "Start position in chromosome",
    "End position in chromosome",
    "Name of item",
    "Score from 0-1000",
    "+ or - for strand",
    "Start of where display should be thick",
    "End of where display should be thick",
    "Used as itemRgb as of 2004-11-22",
    "Number of blocks",
    "Comma separated list of block sizes",
    "Start positions relative to chromStart",
];

/// Turn a writer's declared columns into the autoSql its header points at.
///
/// The first three are the coordinates whatever the caller called them, and go
/// out under the standard names and types: a browser reading the file expects
/// `chrom`, `chromStart`, `chromEnd`, and the writer's own names are only ever
/// used to address the values it is handed.
pub fn build_auto_sql(fields: &indexmap::IndexMap<String, String>) -> Result<BedAutoSql> {
    validate_bed_field_names(fields, "fields")?;
    for (name, kind) in fields {
        if !BED_FIELD_TYPES.contains(&kind.as_str()) {
            return Err(Error::invalid(format!(
                "field type {kind} of {name} invalid (string, int, uint, float)"
            )));
        }
        // autoSql separates a type from its names with whitespace, its names
        // from each other with commas, and a declaration from the next with a
        // semicolon. A name carrying any of those writes a block this library's
        // own reader parses into a different number of fields than the header
        // declares — a file written here and then refused here.
        if let Some(bad) = name
            .chars()
            .find(|c| c.is_whitespace() || matches!(c, ';' | ',' | '(' | ')' | '"'))
        {
            return Err(Error::invalid(format!(
                "field name {name:?} contains {bad:?}, which autosql uses as punctuation                  (names may not hold whitespace, ';', ',', parentheses or quotes)"
            )));
        }
        if name.is_empty() {
            return Err(Error::invalid("a field name may not be empty"));
        }
    }

    let named = |index: usize, given: &str| -> String {
        if index < 3 {
            BED_FIELD_NAMES[index].to_string()
        } else {
            given.to_string()
        }
    };

    let mut text = String::from("table gwseqBed\n\"gwseq_io bed entries\"\n(\n");
    for (index, (given, kind)) in fields.iter().enumerate() {
        let name = named(index, given);
        let standard =
            index < 3 || (index < BED_FIELD_NAMES.len() && given == BED_FIELD_NAMES[index]);
        let kind = if index < 3 {
            if index == 0 {
                "string"
            } else {
                "uint"
            }
        } else {
            kind.as_str()
        };
        let comment = if standard {
            BED_FIELD_COMMENTS[index]
        } else {
            &name
        };
        text.push_str(&format!("{kind} {name};\t\"{comment}\"\n"));
    }
    text.push_str(")\n");

    let mut defined_field_count = 0u16;
    for (index, (given, _)) in fields.iter().enumerate() {
        let name = named(index, given);
        if index >= BED_FIELD_NAMES.len() || name != BED_FIELD_NAMES[index] {
            break;
        }
        defined_field_count += 1;
    }

    Ok(BedAutoSql {
        text,
        field_count: to_bbi_u32(fields.len() as i64, "fieldCount")? as u16,
        defined_field_count,
    })
}

#[cfg(test)]
mod tests {
    use super::*;
    use crate::source::testing::MemorySource;

    /// A minimal but valid 64-byte header, plus whatever follows.
    fn header_bytes(magic: u32, version: u16, zoom_levels: u16) -> Vec<u8> {
        let mut b = Vec::new();
        b.extend_from_slice(&magic.to_le_bytes());
        b.extend_from_slice(&version.to_le_bytes());
        b.extend_from_slice(&zoom_levels.to_le_bytes());
        b.extend_from_slice(&1000u64.to_le_bytes()); // chromTreeOffset
        b.extend_from_slice(&2000u64.to_le_bytes()); // fullDataOffset
        b.extend_from_slice(&3000u64.to_le_bytes()); // fullIndexOffset
        b.extend_from_slice(&9u16.to_le_bytes()); // fieldCount
        b.extend_from_slice(&6u16.to_le_bytes()); // definedFieldCount
        b.extend_from_slice(&400u64.to_le_bytes()); // autoSqlOffset
        b.extend_from_slice(&500u64.to_le_bytes()); // totalSummaryOffset
        b.extend_from_slice(&32768u32.to_le_bytes()); // uncompressBufSize
        b.extend_from_slice(&0u64.to_le_bytes()); // reserved
        assert_eq!(b.len(), 64);
        b
    }

    #[test]
    fn reads_a_bigwig_header() {
        let source = MemorySource::new(header_bytes(super::super::BIGWIG_MAGIC, 4, 0));
        let h = read_header(&source).unwrap();
        assert_eq!(h.kind, BbiKind::BigWig);
        assert_eq!(h.version, 4);
        assert_eq!(h.chr_tree_offset, 1000);
        assert_eq!(h.full_index_offset, 3000);
        assert_eq!(h.field_count, 9);
        assert_eq!(h.uncompress_buffer_size, 32768);
    }

    #[test]
    fn a_bigbed_magic_gives_the_other_kind() {
        let source = MemorySource::new(header_bytes(super::super::BIGBED_MAGIC, 4, 0));
        assert_eq!(read_header(&source).unwrap().kind, BbiKind::BigBed);
    }

    #[test]
    fn a_byte_swapped_file_is_refused_by_name() {
        for magic in [BIGWIG_MAGIC_SWAPPED, BIGBED_MAGIC_SWAPPED] {
            let source = MemorySource::new(header_bytes(magic, 4, 0));
            let err = read_header(&source).unwrap_err().to_string();
            assert!(err.contains("incompatible endianness"), "{err}");
        }
    }

    #[test]
    fn something_else_entirely_is_not_a_bbi_file() {
        let source = MemorySource::new(header_bytes(0xDEAD_BEEF, 4, 0));
        let err = read_header(&source).unwrap_err().to_string();
        assert!(err.contains("not a bigwig or bigbed file"), "{err}");
    }

    #[test]
    fn an_old_version_names_the_floor() {
        let source = MemorySource::new(header_bytes(super::super::BIGWIG_MAGIC, 2, 0));
        let err = read_header(&source).unwrap_err().to_string();
        assert!(err.contains("version 2 unsupported (>= 3)"), "{err}");
    }

    #[test]
    fn a_truncated_header_is_corrupt_not_a_panic() {
        let source = MemorySource::new(vec![0u8; 40]);
        assert!(matches!(
            read_header(&source),
            Err(crate::error::Error::Corrupt { .. })
        ));
    }

    #[test]
    fn zoom_headers_follow_the_common_one() {
        let mut b = header_bytes(super::super::BIGWIG_MAGIC, 4, 2);
        for (reduction, data, index) in [(10u32, 100u64, 200u64), (40, 300, 400)] {
            b.extend_from_slice(&reduction.to_le_bytes());
            b.extend_from_slice(&0u32.to_le_bytes());
            b.extend_from_slice(&data.to_le_bytes());
            b.extend_from_slice(&index.to_le_bytes());
        }
        let source = MemorySource::new(b);
        let zooms = read_zoom_headers(&source, 2).unwrap();
        assert_eq!(zooms.len(), 2);
        assert_eq!(zooms[0].reduction_level, 10);
        assert_eq!(zooms[1].index_offset, 400);
        assert!(read_zoom_headers(&source, 0).unwrap().is_empty());
    }

    #[test]
    fn a_summary_offset_of_zero_means_the_file_carries_none() {
        let source = MemorySource::new(vec![0u8; 8]);
        let s = read_total_summary(&source, 0).unwrap();
        assert_eq!(s.bases_covered, 0);
        assert!(s.min_value.is_nan() && s.max_value.is_nan());
    }

    #[test]
    fn reads_a_total_summary() {
        let mut b = vec![0u8; 64];
        b.extend_from_slice(&1234u64.to_le_bytes());
        b.extend_from_slice(&(-1.5f64).to_le_bytes());
        b.extend_from_slice(&9.25f64.to_le_bytes());
        b.extend_from_slice(&100.0f64.to_le_bytes());
        b.extend_from_slice(&500.0f64.to_le_bytes());
        let source = MemorySource::new(b);
        let s = read_total_summary(&source, 64).unwrap();
        assert_eq!(s.bases_covered, 1234);
        assert_eq!((s.min_value, s.max_value), (-1.5, 9.25));
        assert_eq!((s.sum_data, s.sum_squared), (100.0, 500.0));
    }

    #[test]
    fn the_data_tree_magic_is_checked_both_ways() {
        let source = MemorySource::new(DATA_TREE_MAGIC.to_le_bytes().to_vec());
        assert!(check_data_tree_magic(&source, 0).is_ok());
        let source = MemorySource::new(DATA_TREE_MAGIC_SWAPPED.to_le_bytes().to_vec());
        let err = check_data_tree_magic(&source, 0).unwrap_err().to_string();
        assert!(err.contains("incompatible endianness"), "{err}");
        let source = MemorySource::new(vec![1, 2, 3, 4]);
        assert!(check_data_tree_magic(&source, 0).is_err());
    }
}

#[cfg(test)]
mod write_tests {
    use super::*;
    use crate::source::testing::MemorySource;

    fn fields(pairs: &[(&str, &str)]) -> indexmap::IndexMap<String, String> {
        pairs
            .iter()
            .map(|(a, b)| ((*a).to_string(), (*b).to_string()))
            .collect()
    }

    #[test]
    fn a_written_header_reads_back_as_itself() {
        let header = BbiHeader {
            kind: BbiKind::BigBed,
            version: BBI_OUTPUT_VERSION,
            zoom_levels: 7,
            chr_tree_offset: 4096,
            full_data_offset: 512,
            full_index_offset: 2048,
            field_count: 6,
            defined_field_count: 4,
            auto_sql_offset: 368,
            total_summary_offset: 400,
            uncompress_buffer_size: 32768,
        };
        let mut bytes = write_header(&header).unwrap();
        // The magic is stamped separately, so the bytes as written are refused.
        let err = read_header(&MemorySource::new(bytes.clone())).unwrap_err();
        assert!(err.to_string().contains("not a bigwig or bigbed"), "{err}");

        bytes[..4].copy_from_slice(&super::super::BIGBED_MAGIC.to_le_bytes());
        let read = read_header(&MemorySource::new(bytes)).unwrap();
        assert_eq!(read.kind, BbiKind::BigBed);
        assert_eq!(read.version, header.version);
        assert_eq!(read.zoom_levels, header.zoom_levels);
        assert_eq!(read.chr_tree_offset, header.chr_tree_offset);
        assert_eq!(read.full_data_offset, header.full_data_offset);
        assert_eq!(read.full_index_offset, header.full_index_offset);
        assert_eq!(read.field_count, header.field_count);
        assert_eq!(read.defined_field_count, header.defined_field_count);
        assert_eq!(read.auto_sql_offset, header.auto_sql_offset);
        assert_eq!(read.total_summary_offset, header.total_summary_offset);
        assert_eq!(read.uncompress_buffer_size, header.uncompress_buffer_size);
    }

    #[test]
    fn a_written_zoom_header_and_summary_read_back_as_themselves() {
        let zooms = [
            ZoomHeader {
                reduction_level: 10,
                data_offset: 100,
                index_offset: 200,
            },
            ZoomHeader {
                reduction_level: 40,
                data_offset: 300,
                index_offset: 400,
            },
        ];
        let mut bytes = vec![0u8; BBI_HEADER_SIZE as usize];
        for z in &zooms {
            bytes.extend_from_slice(&write_zoom_header(z));
        }
        let read = read_zoom_headers(&MemorySource::new(bytes), 2).unwrap();
        assert_eq!(read.len(), 2);
        assert_eq!(read[1].reduction_level, 40);
        assert_eq!(read[1].data_offset, 300);
        assert_eq!(read[1].index_offset, 400);

        let summary = TotalSummary {
            bases_covered: 1234,
            min_value: -1.5,
            max_value: 9.25,
            sum_data: 1000.0,
            sum_squared: 5000.0,
        };
        let read = read_total_summary(&MemorySource::new(write_total_summary(&summary)), 0);
        // Offset 0 means "no summary" to a reader, so it is read at 1 instead.
        assert_eq!(read.unwrap().bases_covered, 0);
        let mut padded = vec![0u8];
        padded.extend_from_slice(&write_total_summary(&summary));
        let read = read_total_summary(&MemorySource::new(padded), 1).unwrap();
        assert_eq!(read.bases_covered, 1234);
        assert_eq!(read.min_value, -1.5);
        assert_eq!(read.sum_squared, 5000.0);
    }

    #[test]
    fn an_empty_summary_writes_zero_extremes_rather_than_nan() {
        let bytes = write_total_summary(&TotalSummary::default());
        let mut padded = vec![0u8];
        padded.extend_from_slice(&bytes);
        let read = read_total_summary(&MemorySource::new(padded), 1).unwrap();
        assert_eq!(read.min_value, 0.0);
        assert_eq!(read.max_value, 0.0);
    }

    #[test]
    fn built_autosql_reads_back_as_the_fields_it_describes() {
        // The round trip is the point: a writer that emitted a declaration this
        // library's own reader refuses would produce an unreadable file.
        let declared = fields(&[
            ("chr", "string"),
            ("start", "uint"),
            ("end", "uint"),
            ("name", "string"),
            ("score", "uint"),
        ]);
        let sql = build_auto_sql(&declared).unwrap();
        assert_eq!(sql.field_count, 5);
        // chr/chromStart/chromEnd/name/score are the first five standard ones.
        assert_eq!(sql.defined_field_count, 5);

        let mut bytes = vec![0u8];
        bytes.extend_from_slice(sql.text.as_bytes());
        bytes.push(0);
        let read = read_auto_sql(&MemorySource::new(bytes), 1, 5).unwrap();
        assert_eq!(
            read.keys().map(String::as_str).collect::<Vec<_>>(),
            ["chrom", "chromStart", "chromEnd", "name", "score"]
        );
        assert_eq!(read["chrom"], "string");
        assert_eq!(read["score"], "uint");
    }

    #[test]
    fn a_non_standard_column_stops_the_defined_field_count() {
        let declared = fields(&[
            ("chr", "string"),
            ("start", "uint"),
            ("end", "uint"),
            ("pvalue", "float"),
            ("score", "uint"),
        ]);
        let sql = build_auto_sql(&declared).unwrap();
        // Three coordinates, then a column BED does not know.
        assert_eq!(sql.defined_field_count, 3);
        assert_eq!(sql.field_count, 5);
        // And it keeps its own name and type, with itself as the comment.
        assert!(
            sql.text.contains("float pvalue;\t\"pvalue\"\n"),
            "{}",
            sql.text
        );
    }

    #[test]
    fn a_bad_field_name_or_type_is_refused_before_anything_is_written() {
        let err = build_auto_sql(&fields(&[
            ("name", "string"),
            ("start", "uint"),
            ("end", "uint"),
        ]))
        .unwrap_err()
        .to_string();
        assert!(err.contains("missing or misplaced"), "{err}");
        let err = build_auto_sql(&fields(&[
            ("chr", "string"),
            ("start", "uint"),
            ("end", "uint"),
            ("x", "double"),
        ]))
        .unwrap_err()
        .to_string();
        assert!(err.contains("field type double"), "{err}");
    }

    #[test]
    fn a_coordinate_past_32_bits_is_refused_rather_than_truncated() {
        assert_eq!(to_bbi_u32(4_294_967_295, "chromEnd").unwrap(), u32::MAX);
        let err = to_bbi_u32(4_294_967_296, "chromEnd")
            .unwrap_err()
            .to_string();
        assert!(err.contains("chromEnd 4294967296"), "{err}");
        assert!(to_bbi_u32(-1, "chromStart").is_err());
    }

    #[test]
    fn narrowing_a_summary_clamps_rather_than_reaching_infinity() {
        assert_eq!(to_bbi_f32(1.0e300), f32::MAX);
        assert_eq!(to_bbi_f32(-1.0e300), f32::MIN);
        assert_eq!(to_bbi_f32(1.5), 1.5f32);
    }
}

#[cfg(test)]
mod autosql_tests {
    use super::*;

    fn fields(pairs: &[(&str, &str)]) -> indexmap::IndexMap<String, String> {
        pairs
            .iter()
            .map(|(k, v)| ((*k).to_string(), (*v).to_string()))
            .collect()
    }

    /// A real bigBed declares `enum("+", "-", ".") strand;`, whose type carries
    /// spaces inside its parentheses. Split on the first whitespace, that line
    /// gave the type `enum("+",` and one field too few — and the whole file was
    /// then refused for a field-count mismatch.
    #[test]
    fn a_type_with_spaces_inside_its_brackets_is_one_token() {
        for (line, kind, rest) in [
            (
                "enum(\"+\", \"-\", \".\") strand;",
                "enum(\"+\", \"-\", \".\")",
                "strand;",
            ),
            ("set(a, b) flags;", "set(a, b)", "flags;"),
            (
                "int[blockCount] blockSizes;",
                "int[blockCount]",
                "blockSizes;",
            ),
            ("uint  score;", "uint", "score;"),
            ("string name; comment", "string", "name; comment"),
        ] {
            assert_eq!(split_once_whitespace(line), Some((kind, rest)), "{line}");
        }
        // Nothing to split on at all.
        assert_eq!(split_once_whitespace("uint"), None);
    }

    /// A name carrying autoSql's own punctuation writes a block this reader
    /// parses into a different field count than the header declares — a file
    /// written here and refused here.
    #[test]
    fn a_field_name_may_not_hold_autosql_punctuation() {
        for bad in [
            "two words",
            "semi;colon",
            "com,ma",
            "paren(s)",
            "quo\"te",
            "",
        ] {
            let f = fields(&[
                ("chr", "string"),
                ("start", "uint"),
                ("end", "uint"),
                (bad, "string"),
            ]);
            assert!(
                build_auto_sql(&f).is_err(),
                "{bad:?} was accepted as a field name"
            );
        }
        let ok = fields(&[
            ("chr", "string"),
            ("start", "uint"),
            ("end", "uint"),
            ("itemRgb", "string"),
        ]);
        assert!(build_auto_sql(&ok).is_ok());
    }
}