1 
use std::path::PathBuf;
2 

            

        

            

                3
            

            

                
            

            

                
use color_eyre::{eyre::eyre, Result};

            

        

            

                4
            

            

                
            

            

                
use rayon::prelude::*;

            

        

            

                5
            

            

                
            

            

                
use rust_htslib::bcf::{IndexedReader, Read};

            

        

            

                6
            

            

                
            

            

                
use termion::color;

            

        

            

                7
            

            

                
            

            

                

            

        

            

                8
            

            

                
            

            

                
#[doc(hidden)]

            

        

            

                9
            

            

                12
            

            

                
#[tracing::instrument]

            

        

            

                10
            

            

                6
            

            

                
pub fn run(path: PathBuf, bp_per_snp: i64, npipes: i64) -> Result<()> {

            

        

            

                11
            

            

                
            

            

                
    let contigs = crate::core::get_htslib_bcf_contigs(&path)?;

            

        

            

                12
            

            

                
            

            

                

            

        

            

                13
            

            

                
            

            

                
    let lines = contigs

            

        

            

                14
            

            

                
            

            

                
        .par_iter()

            

        

            

                15
            

            

                
            

            

                
        .map_with(path, |path: &mut PathBuf, (contig, contig_size)| {

            

        

            

                16
            

            

                
            

            

                
            // Fetch readers

            

        

            

                17
            

            

                6
            

            

                
            let mut reader = IndexedReader::from_path(path)?;

            

        

            

                18
            

            

                6
            

            

                
            let header_view = reader.header();

            

        

            

                19
            

            

                
            

            

                

            

        

            

                20
            

            

                6
            

            

                
            let rid = header_view.name2rid(contig.as_bytes())?;

            

        

            

                21
            

            

                6
            

            

                
            reader.fetch(rid, 0, None)?;

            

        

            

                22
            

            

                
            

            

                

            

        

            

                23
            

            

                6
            

            

                
            let records = reader.records();

            

        

            

                24
            

            

                
            

            

                

            

        

            

                25
            

            

                
            

            

                
            // Program logic

            

        

            

                26
            

            

                
            

            

                
            let counts_per_window =

            

        

            

                27
            

            

                6
            

            

                
                counts_per_window(records, *contig_size, npipes)?;

            

        

            

                28
            

            

                
            

            

                
            let line = create_line(counts_per_window, contig, bp_per_snp);

            

        

            

                29
            

            

                
            

            

                

            

        

            

                30
            

            

                
            

            

                
            Ok(line)

            

        

            

                31
            

            

                6
            

            

                
        })

            

        

            

                32
            

            

                
            

            

                
        .collect::<Result<Vec<String>>>()?;

            

        

            

                33
            

            

                
            

            

                

            

        

            

                34
            

            

                
            

            

                
    for line in lines {

            

        

            

                35
            

            

                
            

            

                
        println!("{line}");

            

        

            

                36
            

            

                
            

            

                
    }

            

        

            

                37
            

            

                
            

            

                

            

        

            

                38
            

            

                
            

            

                
    Ok(())

            

        

            

                39
            

            

                
            

            

                
}

            

        

            

                40
            

            

                
            

            

                

            

        

            

                41
            

            

                6
            

            

                
fn counts_per_window<R: Read>(

            

        

            

                42
            

            

                6
            

            

                
    records: rust_htslib::bcf::Records<R>,

            

        

            

                43
            

            

                6
            

            

                
    contig_size: i64,

            

        

            

                44
            

            

                6
            

            

                
    npipes: i64,

            

        

            

                45
            

            

                6
            

            

                
) -> Result<Vec<i64>> {

            

        

            

                46
            

            

                6
            

            

                
    let window_size = contig_size / npipes;

            

        

            

                47
            

            

                6
            

            

                

            

        

            

                48
            

            

                6
            

            

                
    // Count hits per window in a single iteration

            

        

            

                49
            

            

                6
            

            

                
    let mut window = 1;

            

        

            

                50
            

            

                6
            

            

                
    let mut counter = 1;

            

        

            

                51
            

            

                6
            

            

                
    let mut counts_per_window = vec![];

            

        

            

                52
            

            

                6
            

            

                
    let mut record_counter = 0;

            

        

            

                53
            

            

                384
            

            

                
    for record in records {

            

        

            

                54
            

            

                378
            

            

                
        record_counter += 1;

            

        

            

                55
            

            

                378
            

            

                
        let record = record?;

            

        

            

                56
            

            

                
            

            

                

            

        

            

                57
            

            

                378
            

            

                
        if record.pos() < window * window_size {

            

        

            

                58
            

            

                318
            

            

                
            counter += 1;

            

        

            

                59
            

            

                318
            

            

                
        } else {

            

        

            

                60
            

            

                60
            

            

                
            counts_per_window.push(window_size / counter);

            

        

            

                61
            

            

                60
            

            

                
            window += 1;

            

        

            

                62
            

            

                60
            

            

                
            counter = 1;

            

        

            

                63
            

            

                60
            

            

                
        }

            

        

            

                64
            

            

                
            

            

                
    }

            

        

            

                65
            

            

                
            

            

                

            

        

            

                66
            

            

                6
            

            

                
    if counts_per_window.len() != npipes as usize {

            

        

            

                67
            

            

                
            

            

                
        counts_per_window.push(window_size / counter);

            

        

            

                68
            

            

                6
            

            

                
    }

            

        

            

                69
            

            

                6
            

            

                
    while counts_per_window.len() != npipes as usize {

            

        

            

                70
            

            

                
            

            

                
        counts_per_window.push(window_size)

            

        

            

                71
            

            

                
            

            

                
    }

            

        

            

                72
            

            

                
            

            

                

            

        

            

                73
            

            

                6
            

            

                
    if record_counter < 99 {

            

        

            

                74
            

            

                6
            

            

                
        Err(eyre!("Contig has under 100 records: {record_counter}"))

            

        

            

                75
            

            

                
            

            

                
    } else {

            

        

            

                76
            

            

                
            

            

                
        Ok(counts_per_window)

            

        

            

                77
            

            

                
            

            

                
    }

            

        

            

                78
            

            

                6
            

            

                
}

            

        

            

                79
            

            

                
            

            

                

            

        

            

                80
            

            

                
            

            

                
fn create_line(window_counts: Vec<i64>, contig: &str, bp_per_snp: i64) -> String {

            

        

            

                81
            

            

                
            

            

                
    let start = format!("\n{}{contig}\t", color::Fg(color::Reset));

            

        

            

                82
            

            

                
            

            

                

            

        

            

                83
            

            

                
            

            

                
    window_counts.iter().fold(start, |acc, count| {

            

        

            

                84
            

            

                
            

            

                
        if count <= &bp_per_snp {

            

        

            

                85
            

            

                
            

            

                
            format!("{acc}{}|", color::Fg(color::AnsiValue::rgb(0, 5, 1)))

            

        

            

                86
            

            

                
            

            

                
        } else {

            

        

            

                87
            

            

                
            

            

                
            format!("{acc}{}-", color::Fg(color::Magenta))

            

        

            

                88
            

            

                
            

            

                
        }

            

        

            

                89
            

            

                
            

            

                
    })

            

        

            

                90
            

            

                
            

            

                
}