desc_seq_lib/<

pub mod breakpoint;
pub mod genomic_graph;
// mod phase;

use anyhow::{anyhow, Ok, Result};
use fasta::record::Sequence;
use log::info;
use minimap2::{Aligner, Mapping};
use noodles_fasta as fasta;
use num_format::{CustomFormat, Grouping, ToFormattedString, WriteFormatted};
use petgraph::{dot::Dot, prelude::*};
use rust_htslib::bam::{self, Record};
use std::{
    collections::HashMap,
    fmt,
    fs::{self, File},
    io::{BufReader, BufWriter, Write},
    path::PathBuf,
    process::{Command, Stdio},
};
use uuid::Uuid;

#[derive(Debug, Clone)]
pub struct Genome {
    pub chromosomes: HashMap<String, Chromosome>,
}

#[derive(Debug, Clone)]
pub struct Chromosome {
    contigs: Vec<Contig>,
}

#[derive(Debug, Clone)]
pub struct Contig {
    pub id: String,
    // contig seq on ref
    pub mappings: Vec<Mapping>,
    // reads on ref
    pub supporting_records: Option<Vec<Record>>,
    pub sequence: String,
    pub contig_ref: ContigRef,
}

#[derive(Debug, Clone)]
pub enum ContigRef {
    Unique(Mapping),
    Chimeric((Mapping, Mapping)),
    ChimericTriple((Mapping, Mapping, Mapping)),
    ChimericMultiple((Mapping, Vec<Mapping>, Mapping)),
    LeftAmbiguity((Vec<Mapping>, Mapping)),
    RightAmbiguity((Mapping, Vec<Mapping>)),
    Ambigous(Vec<Mapping>),
}

impl fmt::Display for ContigRef {
    fn fmt(&self, fmt: &mut fmt::Formatter) -> fmt::Result {
        let str = match self {
            ContigRef::Unique(m) => mapping_to_string(m),
            ContigRef::Chimeric((a, b)) => {
                format!("{}<->{}", mapping_to_string(a), mapping_to_string(b))
            }
            ContigRef::ChimericMultiple((a, v, b)) => format!(
                "{}<->{}<->{}",
                mapping_to_string(a),
                mappings_to_string(v),
                mapping_to_string(b)
            ),
            ContigRef::LeftAmbiguity((v, b)) => {
                format!("{}<->{}", mappings_to_string(v), mapping_to_string(b))
            }
            ContigRef::RightAmbiguity((a, v)) => {
                format!("{}<->{}", mapping_to_string(a), mappings_to_string(v))
            }
            ContigRef::Ambigous(v) => format!("{}", mappings_to_string(v)),
            ContigRef::ChimericTriple((a, b, c)) => format!(
                "{}<->{}<->{}",
                mapping_to_string(a),
                mapping_to_string(b),
                mapping_to_string(c)
            ),
        };
        fmt.write_str(&str).unwrap();

        std::result::Result::Ok(())
    }
}

impl ContigRef {
    pub fn breakpoints(&self) -> Option<Vec<Mapping>> {
        match self {
            ContigRef::Unique(_) => None,
            ContigRef::Chimeric((a, b)) => Some(vec![a.clone(), b.clone()]),
            ContigRef::ChimericTriple((a, b, c)) => Some(vec![a.clone(), b.clone(), c.clone()]),
            ContigRef::ChimericMultiple(_) => None,
            ContigRef::LeftAmbiguity(_) => None,
            ContigRef::RightAmbiguity(_) => None,
            ContigRef::Ambigous(_) => None,
        }
    }
    pub fn breakpoints_repr(&self) -> Option<Vec<String>> {
        let left = "►";
        let right = "◄";
        let bp_right = "▐";
        let bp_left = "▌";
        let format = CustomFormat::builder()
            .grouping(Grouping::Standard)
            .minus_sign("-")
            .separator("_")
            .build()
            .unwrap();

        let get_sign = |m: &Mapping| -> &str {
            match m.strand {
                minimap2::Strand::Forward => left,
                minimap2::Strand::Reverse => right,
            }
        };

        let uk = "UNKNOWN".to_string();
        let mut res = Vec::new();
        if let Some(breakpoints) = self.breakpoints() {
            for v in breakpoints.windows(2) {
                let mut bp_string = format!("{}:", v[0].target_name.clone().unwrap_or(uk.clone()));
                let _ = bp_string
                    .write_formatted(&v[0].target_end, &format)
                    .unwrap();
                bp_string = format!(
                    "{bp_string}{}{bp_left}{bp_right}{}{}:",
                    get_sign(&v[0]),
                    get_sign(&v[1]),
                    v[1].target_name.clone().unwrap_or(uk.clone())
                );
                let _ = bp_string
                    .write_formatted(&v[1].target_start, &format)
                    .unwrap();
                res.push(bp_string);
            }
        }

        if !res.is_empty() {
            Some(res)
        } else {
            None
        }
    }
    pub fn desc(&self) -> Option<String> {
        let uk = "UNKNOWN".to_string();
        let to_desc = |v: &mut Vec<Mapping>| -> String {
            v.sort_by(|a, b| a.query_start.cmp(&b.query_start));
            let v: Vec<String> = v
                .into_iter()
                .map(|e| {
                    let strand = match e.strand {
                        minimap2::Strand::Forward => "",
                        minimap2::Strand::Reverse => "_rev",
                    };
                    format!(
                        "{}:{}_{}{}",
                        e.target_name.clone().unwrap_or(uk.clone()),
                        e.target_start,
                        e.target_end,
                        strand
                    )
                })
                .collect();
            format!("[{}]", v.join(";"))
        };

        match self {
            ContigRef::Unique(a) => Some(format!(
                "{}:{}_{}",
                a.target_name.clone().unwrap_or(uk.clone()),
                a.target_start,
                a.target_end
            )),
            ContigRef::Chimeric((a, b)) => Some(to_desc(&mut vec![a.to_owned(), b.to_owned()])),
            ContigRef::ChimericTriple((a, b, c)) => {
                Some(to_desc(&mut vec![a.to_owned(), b.to_owned(), c.to_owned()]))
            }
            ContigRef::ChimericMultiple(_) => None,
            ContigRef::LeftAmbiguity(_) => None,
            ContigRef::RightAmbiguity(_) => None,
            ContigRef::Ambigous(a) => Some(to_desc(&mut a.to_owned())),
        }
    }

    pub fn hgvs(&self) -> Option<String> {
        let uk = "UNKNOWN".to_string();
        match self {
            ContigRef::Unique(_) => None,
            ContigRef::Chimeric((a, b)) => {
                if a.target_name == b.target_name {
                    let chr = a.target_name.clone().unwrap_or(uk.clone());
                    let del_start = a.target_end;
                    let del_end = b.target_start;
                    let hgvs = format!("{chr}:{del_start}_{del_end}");
                    Some(hgvs)
                } else {
                    let a_chr = a.target_name.clone().unwrap_or(uk.clone());
                    let a_bp = a.target_end;
                    let b_chr = b.target_name.clone().unwrap_or(uk.clone());
                    let b_bp = b.target_end;
                    let hgvs = format!("{a_chr}:{a_bp}delins[{b_chr}:{b_bp}]");
                    Some(hgvs)
                }
            }
            ContigRef::ChimericMultiple(_) => None,
            ContigRef::LeftAmbiguity(_) => None,
            ContigRef::RightAmbiguity(_) => None,
            ContigRef::Ambigous(_) => None,
            ContigRef::ChimericTriple((a, b, c)) => {
                let mut v = [a, b, c];
                v.sort_by(|a, b| a.query_start.cmp(&b.query_start));
                let (a, b, c) = (
                    *v.get(0).clone().unwrap(),
                    *v.get(1).clone().unwrap(),
                    *v.get(2).clone().unwrap(),
                );
                let a_target_name = a.target_name.clone().unwrap_or(uk.clone());
                let b_target_name = b.target_name.clone().unwrap_or(uk.clone());
                let c_target_name = c.target_name.clone().unwrap_or(uk.clone());

                // if a_target_name != b_target_name {}

                // Insertions
                // prioritize first len
                let (bp_a_1, bp_a_2) = if a.query_end <= b.query_end {
                    // TODO add inserted nt
                    (
                        (a.target_name.clone().unwrap_or(uk.clone()), a.target_end),
                        (b.target_name.clone().unwrap_or(uk.clone()), b.target_start),
                    )
                } else {
                    let diff = a.query_end - b.query_start;
                    (
                        (a.target_name.clone().unwrap_or(uk.clone()), a.target_end),
                        (
                            b.target_name.clone().unwrap_or(uk.clone()),
                            b.target_start + diff,
                        ),
                    )
                };
                let (bp_b_1, bp_b_2) = if b.query_end <= c.query_end {
                    // TODO add inserted nt
                    (
                        (b.target_name.clone().unwrap_or(uk.clone()), b.target_end),
                        (c.target_name.clone().unwrap_or(uk.clone()), c.target_start),
                    )
                } else {
                    let diff = b.query_end - c.query_start;
                    (
                        (b.target_name.clone().unwrap_or(uk.clone()), b.target_end),
                        (
                            c.target_name.clone().unwrap_or(uk.clone()),
                            c.target_start + diff,
                        ),
                    )
                };
                if bp_a_1.0 == bp_b_2.0 {
                    let hgvs = format!(
                        "{}:{}_{}ins[{}:{}_{}]",
                        bp_a_1.0, bp_a_1.1, bp_b_2.1, bp_a_2.0, bp_a_2.1, bp_b_1.1
                    );
                    Some(hgvs)
                } else {
                    None
                }
            }
        }
    }
}

pub fn mapping_to_string(mapping: &Mapping) -> String {
    let uk = "UNKNOWN".to_string();
    format!(
        "{}:{}-{}({}:{}-{})",
        mapping.target_name.clone().unwrap_or(uk.clone()),
        mapping.target_start,
        mapping.target_end,
        mapping.query_name.clone().unwrap_or(uk),
        mapping.query_start,
        mapping.query_end
    )
}

fn mappings_to_string(mappings: &Vec<Mapping>) -> String {
    let v = mappings
        .iter()
        .map(mapping_to_string)
        .collect::<Vec<String>>();
    v.join("//")
}

pub fn get_ref_pos(mappings: Vec<Mapping>) -> Result<ContigRef> {
    let mut mappings = mappings;
    mappings.sort_by(|a, b| a.query_start.cmp(&b.query_start));

    if mappings.len() == 1 {
        return Ok(ContigRef::Unique(mappings.get(0).unwrap().clone()));
    } else {
        let mut grouped: Vec<Vec<Mapping>> = group_mappings(&mut mappings)?;
        // let mut grouped: VecDeque<Vec<Mapping>> = group_mappings(&mut mappings)?.into();

        if grouped.len() == 1 {
            let r = grouped.into_iter().flat_map(|e| e).collect();
            return Ok(ContigRef::Ambigous(r));
        } else if grouped.len() >= 2 {
            // let first = grouped.pop_back().unwrap();
            // let last = grouped.pop_front().unwrap();
            let first = grouped.first().unwrap().to_vec();
            let last = grouped.last().unwrap().to_vec();
            grouped.remove(0);
            grouped.remove(grouped.len() - 1);
            assert!(first[0].query_start < last[0].query_start);

            if grouped.len() == 0 {
                if first.len() == 1 && last.len() == 1 {
                    return Ok(ContigRef::Chimeric((
                        first.get(0).unwrap().clone(),
                        last.get(0).unwrap().clone(),
                    )));
                } else if first.len() == 1 {
                    return Ok(ContigRef::RightAmbiguity((
                        first.get(0).unwrap().clone(),
                        last.clone(),
                    )));
                } else if last.len() == 1 {
                    return Ok(ContigRef::LeftAmbiguity((
                        first.clone(),
                        last.get(0).unwrap().clone(),
                    )));
                } else {
                    let all: Vec<Mapping> = vec![first, last].into_iter().flat_map(|e| e).collect();
                    return Ok(ContigRef::Ambigous(all));
                }
            }
            if first.len() == 1 && last.len() == 1 {
                if grouped.len() == 1 {
                    return Ok(ContigRef::ChimericTriple((
                        first.get(0).unwrap().clone(),
                        grouped.get(0).unwrap().get(0).unwrap().clone(),
                        last.get(0).unwrap().clone(),
                    )));
                } else {
                    return Ok(ContigRef::ChimericMultiple((
                        first.get(0).unwrap().clone(),
                        grouped.into_iter().flat_map(|e| e).collect(),
                        last.get(0).unwrap().clone(),
                    )));
                }
            } else if first.len() == 1 {
                let right: Vec<Mapping> = vec![grouped.into_iter().flat_map(|e| e).collect(), last]
                    .into_iter()
                    .flat_map(|e| e)
                    .collect();
                return Ok(ContigRef::RightAmbiguity((
                    first.get(0).unwrap().clone(),
                    right,
                )));
            } else if last.len() == 1 {
                let left: Vec<Mapping> = vec![first, grouped.into_iter().flat_map(|e| e).collect()]
                    .into_iter()
                    .flat_map(|e| e)
                    .collect();
                return Ok(ContigRef::LeftAmbiguity((
                    left,
                    last.get(0).unwrap().clone(),
                )));
            } else {
                let all: Vec<Mapping> =
                    vec![first, grouped.into_iter().flat_map(|e| e).collect(), last]
                        .into_iter()
                        .flat_map(|e| e)
                        .collect();
                return Ok(ContigRef::Ambigous(all));
            }
        } else {
            return Ok(ContigRef::Ambigous(
                grouped.into_iter().flat_map(|e| e).collect(),
            ));
        }
    }
}

impl Genome {
    pub fn new() -> Self {
        Genome {
            chromosomes: HashMap::new(),
        }
    }

    pub fn iter(&self) -> std::collections::hash_map::Iter<'_, String, Chromosome> {
        self.chromosomes.iter()
    }

    pub fn contigs(&self) -> impl Iterator<Item = &Contig> {
        self.chromosomes.iter().flat_map(|(_, e)| e.iter())
    }

    pub fn add_contig(
        &mut self,
        id: String,
        mappings: Vec<Mapping>,
        supporting_records: Option<Vec<Record>>,
        sequence: String,
    ) -> Result<()> {
        let mut mappings = mappings;
        mappings.sort_by(|a, b| a.query_start.cmp(&b.query_start));
        let new_contig = Contig {
            id,
            mappings: mappings.clone(),
            supporting_records,
            sequence,
            contig_ref: get_ref_pos(mappings)?,
        };
        // get the category of Mapping
        match new_contig.contig_ref.clone() {
            ContigRef::Unique(contig_mapping) => {
                match self
                    .chromosomes
                    .get_mut(&contig_mapping.target_name.unwrap())
                {
                    Some(chromosome) => {
                        chromosome.contigs.push(new_contig);
                    }
                    None => (),
                }
            }
            ContigRef::Chimeric((a, b)) => {
                let a_target_name = a.target_name.unwrap();
                let b_target_name = b.target_name.unwrap();
                if a_target_name == b_target_name {
                    if let Some(chromosome) = self.chromosomes.get_mut(&a_target_name) {
                        chromosome.contigs.push(new_contig);
                    } else {
                        self.chromosomes.insert(
                            a_target_name,
                            Chromosome {
                                contigs: vec![new_contig],
                            },
                        );
                    }
                } else {
                    let chimeric_name = format!("{}-{}", a_target_name, b_target_name);
                    if let Some(chromosome) = self.chromosomes.get_mut(&chimeric_name) {
                        chromosome.contigs.push(new_contig);
                    } else {
                        self.chromosomes.insert(
                            chimeric_name,
                            Chromosome {
                                contigs: vec![new_contig],
                            },
                        );
                    }
                }
            }
            ContigRef::ChimericMultiple((left, _, right)) => {
                let left_target_name = left.target_name.unwrap();
                let right_target_name = right.target_name.unwrap();
                if left_target_name == right_target_name {
                    if let Some(chromosome) = self.chromosomes.get_mut(&left_target_name) {
                        chromosome.contigs.push(new_contig);
                    } else {
                        self.chromosomes.insert(
                            left_target_name,
                            Chromosome {
                                contigs: vec![new_contig],
                            },
                        );
                    }
                } else {
                    let chimeric_name = format!("{}-{}", left_target_name, right_target_name);
                    if let Some(chromosome) = self.chromosomes.get_mut(&chimeric_name) {
                        chromosome.contigs.push(new_contig);
                    } else {
                        self.chromosomes.insert(
                            chimeric_name,
                            Chromosome {
                                contigs: vec![new_contig],
                            },
                        );
                    }
                }
            }
            _ => {
                if let Some(chromosome) = self.chromosomes.get_mut("Ambigous") {
                    chromosome.contigs.push(new_contig);
                } else {
                    self.chromosomes.insert(
                        "Ambigous".to_string(),
                        Chromosome {
                            contigs: vec![new_contig],
                        },
                    );
                }
            }
        };

        Ok(())
    }

    pub fn add_contig_from_seq(
        &mut self,
        name: String,
        sequence: &[u8],
        aligner: impl Fn(String) -> Result<Vec<Mapping>>,
    ) -> Result<()> {
        let mappings = aligner(String::from_utf8(sequence.to_vec())?)?;
        // println!("{mappings:?}");
        self.add_contig(name, mappings, None, String::from_utf8(sequence.to_vec())?)?;
        Ok(())
    }

    pub fn write_contigs_sequences(&self, dir: &str) {
        for contig in self.contigs() {
            contig.write_fasta(&format!("{dir}/{}.fasta", contig.id))
        }
        // self.iter().for_each(|(_, chr)| {
        //     chr.iter().for_each(|c| c.write_fasta(&format!("{dir}/{}.fasta", c.id)))
        // });
    }

    pub fn from_contigs_sequences(dir: &str) -> Result<Self> {
        let aligner_url = "http://localhost:4444/align";
        let aligner = aligner_client::dist_align(aligner_url.to_string());
        let mut genome = Self::new();
        // let paths = get_ext_paths(dir, "fasta")?;
        let paths = get_contigs_fa_paths(dir)?;
        for path in paths {
            let fa = read_fasta(path.to_str().unwrap())?;
            for (name, sequence) in fa {
                genome.add_contig_from_seq(name, sequence.as_ref(), &aligner)?;
            }
        }
        Ok(genome)
    }
    //
    // pub fn from_dir_bed(dir: &str)  {
    //
    // }

    // pub fn write_records(&self, file: &str) {
    //     let mut records =  Vec::new();
    //     for (name, chromosome) in self.chromosomes.iter() {
    //         for contig in chromosome.iter() {
    //             if let Some(rec) = &contig.supporting_records {
    //                 records.extend(rec.iter());
    //             };
    //         }
    //     }
    //     // header
    //
    //     // writer
    //     rust_htslib::bam::Writer::from_path(path, header, format)
    // }

    pub fn stats(&self) {
        for (k, v) in self.chromosomes.iter() {
            info!("{}:{}", k, v.contigs.len());
        }
    }

    pub fn chromosome(&self, chromosome: &str) -> Option<Vec<Contig>> {
        self.chromosomes
            .get(chromosome)
            .map(|chr| chr.contigs.clone())
    }
}

impl Chromosome {
    pub fn iter(&self) -> std::slice::Iter<'_, Contig> {
        self.contigs.iter()
    }
}

impl Contig {
    // pub fn sort(&mut self) {
    //     // sorting by target order
    //     self.mappings
    //         .sort_by(|a, b| a.target_start.cmp(&b.target_start));
    // }

    pub fn to_igv(&self, dir_path: &str) -> Result<()> {
        let supporting_records = self.supporting_records.clone().ok_or(anyhow!("no reads"))?;
        let contig_name = if let Some(hgvs) = self.hgvs() {
            hgvs
        } else {
            self.id.clone()
        };
        let contig_dir = format!("{dir_path}/{contig_name}");
        fs::create_dir_all(contig_dir.clone())?;

        let fasta_path = format!("{contig_dir}/contig.fa");
        write_fasta(&fasta_path, &vec![(self.id.clone(), self.sequence.clone())]);
        write_fai(&fasta_path);

        let reads_path = format!("{contig_dir}/reads.fa");
        let n_reads = supporting_records
            .clone()
            .into_iter()
            .map(|r| {
                (
                    String::from_utf8(r.qname().to_vec()).unwrap(),
                    String::from_utf8(r.seq().as_bytes()).unwrap(),
                )
            })
            .collect();
        write_fasta(&reads_path, &n_reads);

        let bam_path = format!("{contig_dir}/{}.bam", self.id);
        write_bam(&fasta_path, &reads_path, &bam_path)?;

        let bed_path = format!("{contig_dir}/contig.bed");
        match &self.contig_ref {
            ContigRef::Chimeric((a, b)) => {
                let d = vec![
                    (
                        self.id.clone(),
                        a.query_start,
                        a.query_end,
                        format!(
                            "{}:{}-{}",
                            a.target_name.clone().unwrap(),
                            a.target_start,
                            a.target_end
                        ),
                    ),
                    (
                        self.id.clone(),
                        b.query_start,
                        b.query_end,
                        format!(
                            "{}:{}-{}",
                            b.target_name.clone().unwrap(),
                            b.target_start,
                            b.target_end
                        ),
                    ),
                ];
                write_bed(&bed_path, &d)?;
            }
            ContigRef::ChimericTriple((a, b, c)) => {
                let d: Vec<(String, i32, i32, String)> = [a, b, c]
                    .iter()
                    .map(|r| {
                        (
                            self.id.clone(),
                            r.query_start,
                            r.query_end,
                            format!(
                                "{}:{}-{}",
                                r.target_name.clone().unwrap(),
                                r.target_start,
                                r.target_end
                            ),
                        )
                    })
                    .collect();
                write_bed(&bed_path, &d)?;
            }
            _ => (),
        }

        Ok(())
    }

    // bug cigar len != seq len
    // pub fn write_bam(&self, path: &str) -> Result<()> {
    //     let aligner = Aligner::builder()
    //         .asm5()
    //         .with_threads(8)
    //         .with_cigar()
    //         .with_seq(self.sequence.as_bytes())
    //         .expect("Unable to build index");
    //
    //     let mut mappings = Vec::new();
    //     let supporting_records = self
    //         .supporting_records
    //         .clone()
    //         .ok_or(anyhow!("no supporting records"))?;
    //     for record in supporting_records.iter() {
    //         let seq = record.seq().as_bytes();
    //         let alignment = aligner
    //             .map(&seq, false, false, None, None)
    //             .expect("Unable to align");
    //         mappings.push(alignment);
    //     }
    //     let mut mappings: Vec<_> = mappings.into_iter().flatten().collect();
    //     mappings.sort_by(|a, b| a.target_start.cmp(&b.target_start));
    //
    //     let mut header = bam::Header::new();
    //     let mut sq_record = rust_htslib::bam::header::HeaderRecord::new(b"SQ");
    //     sq_record.push_tag(b"SN", self.id.clone());
    //     sq_record.push_tag(b"LN", self.sequence.len());
    //     header.push_record(&sq_record);
    //
    //     let mut out = bam::Writer::from_path(path, &header, bam::Format::Bam).unwrap();
    //
    //     // copy reverse reads to new BAM file
    //     for mapping in mappings.iter() {
    //         let record = minimap2::htslib::mapping_to_record(
    //             Some(mapping),
    //             self.sequence.as_bytes(),
    //             header.clone(),
    //             None,
    //             Some(Uuid::new_v4().as_bytes()),
    //         );
    //         let _ = out.write(&record);
    //     }
    //     rust_htslib::bam::index::build(path, None, rust_htslib::bam::index::Type::Bai, 1)?;
    //     Ok(())
    // }

    pub fn hgvs(&self) -> Option<String> {
        self.contig_ref.hgvs()
    }

    pub fn desc(&self) -> Option<String> {
        self.contig_ref.desc()
    }

    pub fn breakpoints_repr(&self) -> Option<Vec<String>> {
        self.contig_ref.breakpoints_repr()
    }

    pub fn write_fasta(&self, fasta_path: &str) {
        write_fasta(fasta_path, &vec![(self.id.clone(), self.sequence.clone())]);
    }
}

fn group_mappings(mappings: &mut [Mapping]) -> Result<Vec<Vec<Mapping>>> {
    // sort alignments by query_start
    mappings.sort_by(|a, b| a.query_start.cmp(&b.query_start));

    let mut alignments: Vec<Vec<Mapping>> = vec![];
    // group by overlapps > 30
    for aln in mappings.iter() {
        let last = alignments.last_mut();
        if let Some(l) = last {
            if l.iter()
                .filter(|a| a.query_end - aln.query_start > 30)
                .count()
                > 0
            {
                l.push(aln.clone());
            } else {
                alignments.push(vec![aln.clone()]);
            }
        } else {
            alignments.push(vec![aln.clone()]);
        }
    }

    Ok(alignments)
}

pub fn write_bed(path: &str, d: &[(String, i32, i32, String)]) -> Result<()> {
    let file = File::create(path).unwrap();
    let mut writer = BufWriter::new(file);
    for (chr, start, end, value) in d.iter() {
        let row = format!(
            "{}\n",
            [
                chr.to_string(),
                start.to_string(),
                end.to_string(),
                value.to_string()
            ]
            .join("\t")
        );
        writer.write_all(row.as_bytes())?;
    }
    Ok(())
}

// unique
pub fn write_fastq(fastq_path: &str, d: &Vec<Record>) -> Result<()> {
    let file = File::create(fastq_path)?;
    let mut writer = BufWriter::new(file);
    for record in d {
        let name = String::from_utf8(record.qname().to_vec()).unwrap();
        writer.write_all(format!("@{name}\n").as_bytes())?;
        let seq = record.seq().as_bytes();
        writer.write_all(&seq)?;
        writer.write_all(b"\n+\n")?;
        let qual = record.qual();
        writer.write_all(qual)?;
    }
    Ok(())
}

pub fn write_fasta(fasta_path: &str, d: &Vec<(String, String)>) {
    let file = File::create(fasta_path).unwrap();
    let mut writer = fasta::writer::Builder::default().build_with_writer(file);
    let mut passed = Vec::new();
    for (name, sequence) in d {
        let name = name.to_string();
        if sequence.is_empty() {
            continue;
        }
        if passed.contains(&name) {
            continue;
        }
        passed.push(name.clone());
        let record = fasta::Record::new(
            fasta::record::Definition::new(name.as_str(), None),
            fasta::record::Sequence::from(sequence.as_bytes().to_vec()),
        );
        writer.write_record(&record).unwrap();
    }
}

pub fn write_fai(path: &str) {
    let mut faidx = Command::new("samtools")
        .arg("faidx")
        .arg(path)
        .spawn()
        .expect("Samtools faidx failed to start");
    faidx.wait().unwrap();
}

pub fn write_bam(ref_path: &str, reads_path: &str, bam_path: &str) -> Result<()> {
    let rg_id = uuid::Uuid::new_v4();

    let mm2 = Command::new("minimap2")
        .arg("-t")
        .arg("128")
        .arg("-ax")
        .arg("map-ont")
        .arg("-R")
        .arg(format!(
            "@RG\\tPL:ONTASM_PROM\\tID:ONTASM_{rg_id}\\tSM:{rg_id}\\tLB:ONTASM_NB_PROM"
        ))
        .arg(ref_path)
        .arg(reads_path)
        .stdout(Stdio::piped())
        .stderr(Stdio::piped())
        .spawn()
        .expect("Minimap2 failed to start");

    let view = Command::new("sambamba")
        .arg("view")
        .arg("-h")
        .arg("-S")
        .arg("-t")
        .arg("20")
        .arg("--format=bam")
        .arg("/dev/stdin")
        .stdin(Stdio::from(mm2.stdout.unwrap()))
        .stdout(Stdio::piped())
        .stderr(Stdio::piped())
        .spawn()
        .expect("Sambamba view failed to start");

    let mut sort = Command::new("sambamba")
        .arg("sort")
        .arg("-t")
        .arg("20")
        .arg("/dev/stdin")
        .arg("-o")
        .arg(bam_path)
        .stderr(Stdio::piped())
        .stdin(Stdio::from(view.stdout.unwrap()))
        .spawn()
        .expect("Sambamba sort failed to start");

    sort.wait().unwrap();
    Ok(())
}

pub fn read_fasta(path: &str) -> Result<Vec<(String, Sequence)>> {
    let mut reader = File::open(path)
        .map(BufReader::new)
        .map(fasta::Reader::new)?;

    let mut res = Vec::new();
    for result in reader.records() {
        let record = result?;
        let u = String::from_utf8(record.name().to_vec())?;
        let s = record.sequence().to_owned();
        res.push((u, s));
    }

    Ok(res)
}

// fn get_ext_paths(dir: &str, ext: &str) -> Result<Vec<PathBuf>> {
//     let paths = std::fs::read_dir(dir)?
//         // Filter out all those directory entries which couldn't be read
//         .filter_map(|res| res.ok())
//         // Map the directory entries to paths
//         .map(|dir_entry| dir_entry.path())
//         // Filter out all paths with extensions other than `csv`
//         .filter_map(|path| {
//             if path.extension().map_or(false, |xt| xt == ext) {
//                 Some(path)
//             } else {
//                 None
//             }
//         })
//         .collect::<Vec<_>>();
//     Ok(paths)
// }

fn get_contigs_fa_paths(dir: &str) -> Result<Vec<PathBuf>> {
    let pattern = format!("{}/**/*_flye.fa", dir);
    let fa_paths: Vec<PathBuf> = glob::glob(&pattern)
        .expect("Failed to read glob pattern")
        .filter_map(Result::ok)
        .collect();

    Ok(fa_paths)
}

pub fn dot_graph(
    graph: &StableGraph<String, String>,
    way: &[NodeIndex],
    value: &str,
    color: &str,
    erase: bool,
) -> String {
    let mut g = graph.clone();

    // erase labels
    if erase {
        g.edge_weights_mut().for_each(|e| *e = "".to_string());
    }

    let mut labels = Vec::new();
    for window in way.windows(2) {
        let edge_id = g.find_edge(window[0], window[1]).unwrap();
        let edge = g.edge_weight_mut(edge_id).unwrap();
        let v = if !edge.is_empty() {
            let mut v: Vec<&str> = edge.split(",").filter(|e| !e.is_empty()).collect();
            v.push(value);
            v.join(", ")
        } else {
            value.to_string()
        };
        labels.push(v.clone());
        *edge = v;
    }

    if erase {
        g.retain_edges(|g, i| g.edge_weight(i).unwrap().to_string() != "".to_string());
        g.retain_nodes(|g, n| g.neighbors_undirected(n).count() > 0);
    }

    let mut dot = Dot::new(&g).to_string().replace("\\\"", "");

    labels.sort_by_key(|b| std::cmp::Reverse(b.len()));
    // labels.sort_by(|a, b| b.len().cmp(&a.len()));
    for label in labels {
        let label_str = format!("label = \"{label}\"");
        dot = dot.replace(
            &format!("{label_str} ]\n"),
            &format!("{label_str}, color = \"{color}\" ]\n"),
        );
    }
    dot
}

/// for non overlapping ways
pub fn dot_graph_biall(
    graph: &StableGraph<String, String>,
    ways: &[Vec<NodeIndex>],
    values: Vec<&str>,
    colors: Vec<&str>,
    erase: bool,
) -> String {
    let mut g = graph.clone();
    // erase labels
    if erase {
        g.edge_weights_mut().for_each(|e| *e = "".to_string());
    }

    let mut labels = Vec::new();
    for ((way, value), color) in ways.iter().zip(values.into_iter()).zip(colors.into_iter()) {
        for window in way.windows(2) {
            let edge_id = g.find_edge(window[0], window[1]).unwrap();
            let edge = g.edge_weight_mut(edge_id).unwrap();
            let v = if !edge.is_empty() {
                let mut v: Vec<&str> = edge.split(",").filter(|e| !e.is_empty()).collect();
                v.push(value);
                v.join(", ")
            } else {
                value.to_string()
            };
            labels.push((v.clone(), color));
            *edge = v;
        }
    }

    // g.retain_edges(|g, i| g.edge_weight(i).unwrap().to_string() != "".to_string());
    g.retain_edges(|g, i| !g.edge_weight(i).unwrap().is_empty());
    g.retain_nodes(|g, n| g.neighbors_undirected(n).count() > 0);

    let mut dot = Dot::new(&g).to_string().replace("\\\"", "");

    labels.sort_by(|a, b| b.0.len().cmp(&a.0.len()));
    for (label, color) in labels {
        let label_str = format!("label = \"{label}\"");
        dot = dot.replace(
            &format!("{label_str} ]\n"),
            &format!("{label_str}, color = \"{color}\" ]\n"),
        );
    }
    dot
}

#[cfg(test)]
mod tests {
    use env_logger::Env;

    use super::*;
    use crate::{
        genomic_graph::GenomicGraph,
        // phase::{variants_phasing, write_phases_bed},
    };

    fn init() {
        let _ = env_logger::Builder::from_env(Env::default().default_filter_or("info"))
            .is_test(true)
            .try_init();
    }

    #[test]
    fn it_works() -> Result<()> {
        let _ = env_logger::builder().is_test(true).try_init();
        let contig_fa = "./data_test/contig_2.fa";
        let aligner_url = "http://localhost:4444/align";

        let mut genome = Genome::new();
        let aligner = aligner_client::dist_align(aligner_url.to_string());

        let sequences = read_fasta(contig_fa)?;
        for (name, seq) in sequences {
            genome.add_contig_from_seq(name.clone(), &seq.as_ref().to_vec(), &aligner)?;
            let mut seqc: Vec<u8> = seq.complement().map(|e| e.unwrap()).collect();
            seqc.reverse();
            genome.add_contig_from_seq(format!("{name}_rev"), &seqc, &aligner)?;
            println!("Sending");
        }
        genome.iter().for_each(|(_, c)| {
            c.iter().for_each(|cont| {
                println!("{}", cont.contig_ref.desc().unwrap());
            });
        });

        Ok(())
    }

    #[test]
    fn test_graph() -> Result<()> {
        init();

        let case = "SALICETTO";
        let chrom = vec!["chr10"];
        info!("This record will be captured by `cargo test`");

        let dir = format!("/data/longreads_basic_pipe/{case}/diag/asm_contis");

        // Load from fasta in dir.
        let genome = Genome::from_contigs_sequences(&dir)?;
        genome.stats();
        let mut genomic_graph = GenomicGraph::from_genome(&genome);

        let sens = vec![true, false];
        let pos = vec![0, i32::MAX];
        let mut all_ways = Vec::new();
        if chrom.len() > 1 {
            (0..4).into_iter().for_each(|i| {
                let start_pos = if i < 2 { 0 } else { i32::MAX };
                let end_pos = pos[i % 2];

                (0..4).into_iter().for_each(|i| {
                    let start_sens = if i < 2 { true } else { false };
                    let end_sens = sens[i % 2];
                    (0..4).into_iter().for_each(|i| {
                        let start_chr = if i < 2 { chrom[0] } else { chrom[1] };
                        let end_chr = chrom[i % 2];
                        let start = (start_sens, start_chr, start_pos);
                        let end = (end_sens, end_chr, end_pos);

                        let (oriented_graph, _integrated_graph, ways) =
                            genomic_graph.ways(start, end);

                        let dot = oriented_graph.dot_graph();
                        println!("dot\n{dot}");

                        for (_i, way) in ways.iter().enumerate() {
                            let s = way
                                .iter()
                                .map(|(_, _, _, s)| s.to_string())
                                .collect::<Vec<String>>()
                                .join("");
                            all_ways.push(s);
                        }
                    });
                });
            });
        } else {
            let start_chr = chrom[0];
            let end_chr = chrom[0];
            let start = (true, start_chr, 0);
            let end = (true, end_chr, i32::MAX);

            let (oriented_graph, _integrated_graph, ways) = genomic_graph.ways(start, end);

            let dot = oriented_graph.dot_graph();
            println!("dot\n{dot}");

            for (_i, way) in ways.iter().enumerate() {
                let s = way
                    .iter()
                    .map(|(_, _, _, s)| s.to_string())
                    .collect::<Vec<String>>()
                    .join("");
                all_ways.push(s);
            }
        }

        all_ways.dedup();
        all_ways
            .iter()
            .enumerate()
            .for_each(|(i, s)| println!("{i}.\t{s}"));

        // let s = Dot::new(&integrated_graph).to_string().replace("\\\"", "");
        // let x11_colors: Vec<String> = vec![
        //     String::from("Red"),
        //     String::from("Green"),
        //     String::from("Blue"),
        //     String::from("Cyan"),
        //     String::from("Magenta"),
        //     String::from("Yellow"),
        //     String::from("DarkRed"),
        //     String::from("DarkGreen"),
        //     String::from("DarkBlue"),
        //     String::from("DarkCyan"),
        //     String::from("DarkMagenta"),
        //     String::from("DarkYellow"),
        //     String::from("LightRed"),
        //     String::from("LightGreen"),
        //     String::from("LightBlue"),
        //     String::from("LightCyan"),
        //     String::from("LightMagenta"),
        //     String::from("LightYellow"),
        //     String::from("Orange"),
        //     String::from("Brown"),
        //     String::from("Beige"),
        // ];
        // let mut s = s.clone();
        // ways.iter().enumerate().for_each(|(i, _)| {
        //     s = s.replace(
        //         &format!("[ label = \"{}\" ]", i + 1),
        //         &format!(
        //             "[ label = \"{}\" color = \"{}\" ]",
        //             i + 1,
        //             x11_colors[i].to_string()
        //         ),
        //     );
        // });
        // println!("{s}");
        //
        // for (i, way) in ways.iter().enumerate() {
        //     let s = way
        //         .iter()
        //         .map(|(_, _, _, s)| s.to_string())
        //         .collect::<Vec<String>>()
        //         .join("");
        //     println!("{}.\t{s}", i + 1);
        // }
        Ok(())
    }

    #[test]
    fn dir() {
        init();
        let id = "ROBIN";
        info!("This record will be captured by `cargo test`");

        let dir = format!("/data/longreads_basic_pipe/{id}/diag/scan/reads",);

        // Load from fasta in dir.
        let genome = Genome::from_contigs_sequences(&dir).unwrap();
        genome.stats();
        let mut res: Vec<String> = genome
            .iter()
            .flat_map(|(_s, chrom)| {
                chrom.iter().filter_map(|c| c.hgvs())
                // .map(|c| println!("{c}"))
            })
            .collect();
        res.sort();
        res.dedup();
        res.iter().for_each(|s| println!("{s}"));

        // println!("{genome:#?}");
    }
}