|
|
@@ -181,3 +181,34 @@ pub fn range_depths_qual(
|
|
|
}
|
|
|
Ok(depths)
|
|
|
}
|
|
|
+
|
|
|
+pub fn range_len_qual(
|
|
|
+ bam: &mut rust_htslib::bam::IndexedReader,
|
|
|
+ chr: &str,
|
|
|
+ start: i32,
|
|
|
+ stop: i32,
|
|
|
+ mapq: u8,
|
|
|
+) -> Result<Vec<usize>> {
|
|
|
+ bam.fetch((chr, start, stop))?;
|
|
|
+
|
|
|
+ let mut depths = vec![0];
|
|
|
+ depths.resize((stop - start) as usize, 0);
|
|
|
+ for p in bam.pileup() {
|
|
|
+ let pileup = p.context(format!("eRR"))?;
|
|
|
+ let rstart = pileup.pos() as i32;
|
|
|
+
|
|
|
+ if rstart >= start && rstart < stop {
|
|
|
+ let v = depths
|
|
|
+ .get_mut((rstart - start) as usize)
|
|
|
+ .context(format!("Errrr {}", rstart - start))?;
|
|
|
+ *v = pileup
|
|
|
+ .alignments()
|
|
|
+ .filter(|e| e.record().mapq() >= mapq)
|
|
|
+ .map(|e| e.record().seq_len())
|
|
|
+ .sum();
|
|
|
+ } else if rstart >= stop {
|
|
|
+ break;
|
|
|
+ }
|
|
|
+ }
|
|
|
+ Ok(depths)
|
|
|
+}
|
