Skip to content

Commit

Permalink
moving average and base-quality stubs
Browse files Browse the repository at this point in the history
  • Loading branch information
brentp committed Aug 18, 2023
1 parent fe8af74 commit aa2b612
Show file tree
Hide file tree
Showing 5 changed files with 129 additions and 2 deletions.
8 changes: 8 additions & 0 deletions src/bin/commands/trimmer.rs
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,14 @@ pub(crate) struct TrimmerOpts {
#[clap(long, short = 't', default_value = "5")]
threads: usize,

/// Minimum base-quality to keep a base when trimming tails.
#[clap(long, short = 'q', default_value = "20")]
trim_tail_quality: u8,

/// Window size for moving average when trimming tails.
#[clap(long, short = 'w', default_value = "20")]
trim_tail_window: u8,

/// Level of compression to use to compress outputs.
#[clap(long, short = 'c', default_value = "5")]
compression_level: usize,
Expand Down
4 changes: 2 additions & 2 deletions src/bin/main.rs
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ pub mod commands;
use anyhow::Result;
use clap::Parser;
use commands::command::Command;
use commands::{demux::Demux, trimmer::Trimmer};
use commands::{demux::Demux, trimmer::TrimmerOpts};
use enum_dispatch::enum_dispatch;
use env_logger::Env;

Expand All @@ -23,7 +23,7 @@ struct Args {
#[command(version)]
enum Subcommand {
Demux(Demux),
Trimmer(Trimmer),
Trimmer(TrimmerOpts),
}

fn main() -> Result<()> {
Expand Down
78 changes: 78 additions & 0 deletions src/lib/base_quality.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,78 @@
use std::ops::Range;

pub(crate) fn find_oscillating_quals(bqs: &[u8]) -> Range<usize> {
return 0..0;
}

pub(crate) enum Tail {
Left,
Right,
Both,
}

/// A simple moving average calculator.
/// Only requires that T is convertable to f64.
/// Uses space of window * size_of(T) bytes.
struct MovingAverage<T> {
window: usize,
values: Vec<T>,
sum: f64,
idx: usize,
}

impl<T: Copy + Default + std::convert::Into<f64>> MovingAverage<T> {
/// create a new moving average calculator with a window of `window` values.
fn new(window: usize) -> Self {
Self { window, values: vec![T::default(); window], sum: 0.0, idx: 0 }
}

/// push a new value into the moving average calculator and get the new mean.
fn push(&mut self, value: T) -> f64 {
let old_value = self.values[self.idx];
self.values[self.idx] = value;
self.sum = self.sum + value.into() - old_value.into();
self.idx = (self.idx + 1) % self.window;
self.mean()
}

/// get the current mean.
#[inline]
fn mean(&self) -> f64 {
self.sum / (self.window as f64)
}
}

pub(crate) fn find_low_quality_bases(
bqs: &[u8],
min_quality: u8,
window: u8,
tail: Tail,
) -> Range<usize> {
if matches!(tail, Tail::Left | Tail::Both) {
let mut i = 0;
while i < bqs.len() && bqs[i] < min_quality {
i += 1;
}
return 0..i;
}
if matches!(tail, Tail::Right | Tail::Both) {
let mut i = bqs.len() - 1;
while i > 0 && bqs[i] < min_quality {
i -= 1;
}
return i..bqs.len();
}
0..0
}

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

#[test]
fn test_find_oscillating_quals() {
let bqs = b"!\"#$%&'()*+,-./0123456789:;<=>?@ABCDEFGHIJ";
let range = find_oscillating_quals(bqs);
assert_eq!(range, 0..0);
}
}
1 change: 1 addition & 0 deletions src/lib/mod.rs
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
pub mod barcode_matching;
pub mod base_quality;
pub mod pair_overlap;
pub mod samples;

Expand Down
40 changes: 40 additions & 0 deletions src/lib/moving_average.rs
Original file line number Diff line number Diff line change
@@ -0,0 +1,40 @@
/// A simple moving average calculator.
/// Only requires that T is convertable to f64.
/// Uses space of window * size_of(T) bytes.
struct MovingAverage<T> {
window: usize,
values: Vec<T>,
sum: f64,
idx: usize,
}

impl<T: Copy + Default + std::convert::Into<f64>> MovingAverage<T> {
/// create a new moving average calculator with a window of `window` values.
fn new(window: usize) -> Self {
Self { window, values: vec![T::default(); window], sum: 0.0, idx: 0 }
}

/// push a new value into the moving average calculator and get the new mean.
fn push(&mut self, value: T) -> f64 {
let old_value = self.values[self.idx];
self.values[self.idx] = value;
self.sum = self.sum + value.into() - old_value.into();
self.idx = (self.idx + 1) % self.window;
self.mean()
}

/// get the current mean.
#[inline]
fn mean(&self) -> f64 {
self.sum / (self.window as f64)
}
}

// write some tests for the calculator
#[cfg(test)]
mod tests {
use super::*;

#[test]
fn test_moving_average() {}
}

0 comments on commit aa2b612

Please sign in to comment.