binseq 0.11.0

A high efficiency binary format for sequencing data
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145

//! Shared nucleotide encoder used by the BQ and VBQ writers.

use bitnuc_deprec::BitSize;
use rand::{SeedableRng, rngs::SmallRng};

use crate::{Policy, error::Result, policy::RNG_SEED};

/// Encodes nucleotide sequences into a compact multi-bit binary format.
///
/// The encoder maintains internal buffers to avoid repeated allocations and
/// handles invalid nucleotides according to a configurable [`Policy`].
#[derive(Clone)]
pub(crate) struct Encoder {
    /// Bitsize of the nucleotides
    bitsize: BitSize,

    /// Reusable buffers for encoded nucleotides
    sbuffer: Vec<u64>,
    xbuffer: Vec<u64>,

    /// Reusable buffers for invalid nucleotide sequences
    s_ibuf: Vec<u8>,
    x_ibuf: Vec<u8>,

    /// Invalid Nucleotide Policy
    policy: Policy,

    /// Random Number Generator (seeded with `RNG_SEED` for reproducibility)
    rng: SmallRng,
}

impl Encoder {
    /// Initialize a new encoder with the given policy.
    pub fn with_policy(bitsize: BitSize, policy: Policy) -> Self {
        Self {
            bitsize,
            policy,
            sbuffer: Vec::default(),
            xbuffer: Vec::default(),
            s_ibuf: Vec::default(),
            x_ibuf: Vec::default(),
            rng: SmallRng::seed_from_u64(RNG_SEED),
        }
    }

    /// Returns the invalid-nucleotide policy of this encoder.
    pub fn policy(&self) -> Policy {
        self.policy
    }

    /// Encodes a single sequence.
    ///
    /// Will return `None` if the sequence is invalid and the policy does not allow correction.
    pub fn encode_single(&mut self, primary: &[u8]) -> Result<Option<&[u64]>> {
        self.clear();
        if self.bitsize.encode(primary, &mut self.sbuffer).is_err() {
            self.clear();
            if self
                .policy
                .handle(primary, &mut self.s_ibuf, &mut self.rng)?
            {
                self.bitsize.encode(&self.s_ibuf, &mut self.sbuffer)?;
            } else {
                return Ok(None);
            }
        }
        Ok(Some(&self.sbuffer))
    }

    /// Encodes a pair of sequences.
    ///
    /// Will return `None` if either sequence is invalid and the policy does not allow correction.
    pub fn encode_paired(
        &mut self,
        primary: &[u8],
        extended: &[u8],
    ) -> Result<Option<(&[u64], &[u64])>> {
        self.clear();
        if self.bitsize.encode(primary, &mut self.sbuffer).is_err()
            || self.bitsize.encode(extended, &mut self.xbuffer).is_err()
        {
            self.clear();
            if self
                .policy
                .handle(primary, &mut self.s_ibuf, &mut self.rng)?
                && self
                    .policy
                    .handle(extended, &mut self.x_ibuf, &mut self.rng)?
            {
                self.bitsize.encode(&self.s_ibuf, &mut self.sbuffer)?;
                self.bitsize.encode(&self.x_ibuf, &mut self.xbuffer)?;
            } else {
                return Ok(None);
            }
        }
        Ok(Some((&self.sbuffer, &self.xbuffer)))
    }

    /// Clear all buffers and reset the encoder.
    pub fn clear(&mut self) {
        self.sbuffer.clear();
        self.xbuffer.clear();
        self.s_ibuf.clear();
        self.x_ibuf.clear();
    }
}

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

    #[test]
    fn test_encode_single_invalid_ignored() {
        let mut encoder = Encoder::with_policy(BitSize::Two, Policy::IgnoreSequence);
        assert!(encoder.encode_single(b"ACGTNNNN").unwrap().is_none());
    }

    #[test]
    fn test_encode_single_invalid_corrected() {
        let mut encoder = Encoder::with_policy(BitSize::Two, Policy::SetToA);
        assert!(encoder.encode_single(b"ACGTNNNN").unwrap().is_some());
    }

    #[test]
    fn test_encode_paired_invalid_ignored() {
        let mut encoder = Encoder::with_policy(BitSize::Two, Policy::IgnoreSequence);
        assert!(
            encoder
                .encode_paired(b"ACGTNNNN", b"ACGTACGT")
                .unwrap()
                .is_none()
        );
    }

    #[test]
    fn test_encode_paired_invalid_corrected() {
        let mut encoder = Encoder::with_policy(BitSize::Two, Policy::SetToA);
        assert!(
            encoder
                .encode_paired(b"ACGTNNNN", b"NNNNACGT")
                .unwrap()
                .is_some()
        );
    }
}