From 67120cd50d45ebc0693cdfa4f9bce8c235067c47 Mon Sep 17 00:00:00 2001 From: "Mitchell R. Vollger" Date: Wed, 18 Mar 2026 09:51:20 -0600 Subject: [PATCH] fix: --skip-no-m6a no longer skips reads that have m6A but lack MSP annotations The skip condition incorrectly combined m6A emptiness with MSP count, causing all reads to be skipped when the BAM had m6A calls (MM/ML tags) but no pre-computed MSP annotations (as tag). This affected ONT data where dorado provides m6A calls but fiberseq MSP/nucleosome tags are absent. Fixes #100 Co-Authored-By: Claude Opus 4.6 --- src/subcommands/fire.rs | 19 ++++++++++++------- 1 file changed, 12 insertions(+), 7 deletions(-) diff --git a/src/subcommands/fire.rs b/src/subcommands/fire.rs index 446c975a..9cdbfa9d 100644 --- a/src/subcommands/fire.rs +++ b/src/subcommands/fire.rs @@ -72,20 +72,25 @@ pub fn add_fire_to_bam(fire_opts: &mut FireOptions) -> Result<(), anyhow::Error> let n_msps = rec.msp.annotations.len(); if fire_opts.skip_no_m6a || fire_opts.min_msp > 0 || fire_opts.min_ave_msp_size > 0 { - // skip no calls - if rec.m6a.annotations.is_empty() || n_msps == 0 { + // skip reads with no m6a calls + if fire_opts.skip_no_m6a && rec.m6a.annotations.is_empty() { skip_because_no_m6a += 1; continue; } - //let max_msp_len = *rec.msp.lengths.iter().flatten().max().unwrap_or(&0); if n_msps < fire_opts.min_msp { skip_because_num_msp += 1; continue; } - let ave_msp_size = rec.msp.lengths().iter().sum::() / n_msps as i64; - if ave_msp_size < fire_opts.min_ave_msp_size { - skip_because_ave_msp_length += 1; - continue; + if fire_opts.min_ave_msp_size > 0 { + if n_msps == 0 { + skip_because_ave_msp_length += 1; + continue; + } + let ave_msp_size = rec.msp.lengths().iter().sum::() / n_msps as i64; + if ave_msp_size < fire_opts.min_ave_msp_size { + skip_because_ave_msp_length += 1; + continue; + } } } out.write(&rec.record)?;