Add PackedNSeq::from_fastq_with_quality behind file_io feature

Author: RagnarGrootKoerkamp

Cargo.lock

@@ -24,10 +24,16 @@ epserde = { version = "0.8", optional = true }
 # Optional pyclass atrributes for some objects.
 pyo3 = { version = "0.25", features = ["extension-module"], optional = true }
 sux = "0.10.2"
+# Optional, for reading PackedSeq[N]Vec directly from a file
+needletail = { version = "0.6.3", optional = true }
 
 [features]
 # Also needed for tests.
-default = ["rand"]
+default = ["rand", "file_io"]
+
+# Enables functions to read PackedNSeq directly from a file.
+file_io = ["needletail"]
 
 # Hides the `simd` warnings when neither AVX2 nor NEON is detected.
 scalar = ["ensure_simd/scalar"]
+

Cargo.toml

@@ -134,4 +134,36 @@ impl PackedNSeqVec {
             *ambi = (u32::from_ne_bytes(*ambi) | mask).to_ne_bytes();
         }
     }
+
+    #[cfg(feature = "file_io")]
+    pub fn from_fastq_with_quality(path: &std::path::Path, min_qual: u8) -> Self {
+        let mut dna = vec![];
+        let mut qual = vec![];
+
+        let mut seqs = Self::default();
+
+        let mut reader = needletail::parse_fastx_file(&path).unwrap();
+        while let Some(record) = reader.next() {
+            let seqrec = record.expect("Invalid FASTA/Q record");
+            dna.extend_from_slice(&seqrec.seq());
+            qual.extend_from_slice(&seqrec.qual().unwrap());
+            dna.push(b'N');
+            qual.push(0);
+
+            // Convert from ascii dna+qual to packed dna+mask in batches of 16kbp.
+            if dna.len() > 16000 {
+                seqs.push_from_ascii_and_quality(&dna, &qual, min_qual);
+                dna.clear();
+                qual.clear();
+            }
+        }
+
+        if dna.len() > 0 {
+            seqs.push_from_ascii_and_quality(&dna, &qual, min_qual);
+            dna.clear();
+            qual.clear();
+        }
+
+        seqs
+    }
 }