#!/usr/bin/env bash

# exit when any command fails
set -eo pipefail

input_fasta=false
output_dir=""
input_paf=false
resume=false
map_pct_id=95
n_mappings=false
segment_length=10000
block_length=false
mash_kmer=false
min_match_length=47
sparse_factor=0
transclose_batch=10000000
n_haps=false
block_ratio_min=0
pad_max_depth=100
max_path_jump=0
max_edge_jump=0
target_poa_length=4001,4507
poa_params=false
poa_padding=0.03
do_viz=true
do_layout=true
threads=1
poa_threads=0
mapper=wfmash
no_merge_segments=false
do_stats=false
exclude_delim=false
write_maf=false
consensus_spec=false # Disabled for the moment due to PGGB's #133 and #182 issues
consensus_prefix=Consensus_
no_splits=false
multiqc=false
keep_intermediate_files=false
compress=false
normalize=true
vcf_spec=false

if [ $# -eq 0 ];
then
    show_help=true
fi

# read the options
cmd=$0" "$@
TEMP=`getopt -o i:o:a:p:n:s:l:K:k:f:B:H:j:P:O:Me:t:T:vhASY:G:Q:d:I:R:NrmZV: --long input-fasta:,output-dir:,input-paf:,map-pct-id:,n-mappings:,segment-length:,block-length-min:,mash-kmer:,min-match-length:,sparse-factor:,transclose-batch:,n-haps:,path-jump-max:,subpath-min:,edge-jump-max:,threads:,poa-threads:,skip-viz,do-layout,help,no-merge-segments,do-stats,exclude-delim:,poa-length-target:,poa-params:,poa-padding:,write-maf,consensus-spec:,consensus-prefix:,pad-max-depth:,block-id-min:,block-ratio-min:,no-splits,resume,keep-temp-files,multiqc,compress,vcf-spec: -n 'pggb' -- "$@"`
eval set -- "$TEMP"

# extract options and their arguments into variables.
while true ; do
    case "$1" in
        -i|--input-fasta) input_fasta=$2 ; shift 2 ;;
        -o|--output-dir) output_dir=$2 ; shift 2 ;;
        -a|--input-paf) input_paf=$2 ; shift 2 ;;
        -p|--map-pct-id) map_pct_id=$2 ; shift 2 ;;
        -n|--n-mappings) n_mappings=$2 ; shift 2 ;;
        -s|--segment-length) segment_length=$2 ; shift 2 ;;
        -l|--block-length) block_length=$2 ; shift 2 ;;
        -N|--no-splits) no_splits=true ; shift ;;
        -K|--mash-kmer) mash_kmer=$2 ; shift 2 ;;
        -Y|--exclude-delim) exclude_delim=$2 ; shift 2 ;;
        -k|--min-match-length) min_match_length=$2 ; shift 2 ;;
        -f|--sparse-factor) sparse_factor=$2 ; shift 2 ;;
        -B|--transclose-batch) transclose_batch=$2 ; shift 2 ;;
        -H|--n-haps) n_haps=$2 ; shift 2 ;;
        -d|--pad-max-depth) pad_max_depth=$2 ; shift 2 ;;
        -j|--path-jump-max) max_path_jump=$2 ; shift 2 ;;
        -e|--edge-jump-max) max_edge_jump=$2 ; shift 2 ;;
        -G|--poa-length-target) target_poa_length=$2 ; shift 2 ;;
        -P|--poa-params) poa_params=$2 ; shift 2 ;;
        -O|--poa-padding) poa_padding=$2 ; shift 2 ;;
        -M|--write-maf) write_maf=true ; shift ;;
        #-C|--consensus-spec) consensus_spec=$2 ; shift 2 ;;
        -Q|--consensus-prefix) consensus_prefix=$2 ; shift 2 ;;
        -t|--threads) threads=$2 ; shift 2 ;;
        -T|--poa-threads) poa_threads=$2 ; shift 2 ;;
        -v|--skip-viz) do_viz=false ; do_layout=false; shift ;;
        -S|--do-stats) do_stats=true ; shift ;;
        -m|--multiqc) multiqc=true ; shift ;;
        -r|--resume) resume=true ; shift ;;
        -A|--keep-temp-files) keep_intermediate_files=true ; shift ;;
        -Z|--compress) compress=true ; shift ;;
        -V|--vcf-spec) vcf_spec=$2 ; shift 2 ;;
        -h|--help) show_help=true ; shift ;;
        --) shift ; break ;;
        *) echo "$2" "Internal error!" ; exit 1 ;;
    esac
done

# if a PAF is given as input, we set the segment_length by to an arbitrary number, since we will skip the wfmash alignment process
if [[ $input_paf != false ]];
then
    if [[ $segment_length == false ]];
    then
        segment_length=42
    fi
    mapper=EXTERNAL
fi

if [[
       "$input_fasta" == false
    || $n_mappings == false
   ]];
then
    show_help=true
    >&2 echo "ERROR: mandatory arguments -i and -n"
fi

if (( "$n_mappings" < 2 ));
then
    show_help=true
    >&2 echo "ERROR: -n must be greater than or equal to 2"
fi

if [[ $poa_threads == 0 ]];
then
    poa_threads=$threads
fi

block_length_cmd=""
if [[ $block_length != false ]];
then
    block_length_cmd="-l $block_length"
fi
# wfmash auto-sets our block length based on its parameters

mash_kmer_cmd=""
if [[ $mash_kmer != false ]];
then
    mash_kmer_cmd="-k $mash_kmer"
fi
# wfmash auto-sets our kmer length based on its parameters

# our partial order alignment parameters
# poa param suggestions from minimap2
# - asm5, --poa-params 1,19,39,3,81,1, ~0.1 divergence
# - asm10, --poa-params 1,9,16,2,41,1, ~1 divergence
# - asm20, --poa-params 1,4,6,2,26,1, ~5% divergence
# between asm10 and asm20 ~ 1,7,11,2,33,1
poa_params_cmd=""
if [[ $poa_params == false ]];
then
    poa_params_cmd="-P 1,19,39,3,81,1"
else
    if [[ $poa_params == "asm5" ]]; then
        poa_params_cmd="-P 1,19,39,3,81,1"
    elif [[ $poa_params == "asm10" ]]; then
        poa_params_cmd="-P 1,9,16,2,41,1"
    elif [[ $poa_params == "asm15" ]]; then
        poa_params_cmd="-P 1,7,11,2,33,1"
    elif [[ $poa_params == "asm20" ]]; then
        poa_params_cmd="-P 1,4,6,2,26,1"
    else
        poa_params_cmd="-P $poa_params"
    fi
fi

if [[ $n_haps == false ]];
then
    n_haps=$n_mappings
fi

if [ $show_help ];
then
    padding=`printf %${#0}s` # prints as many spaces as the length of $0
    echo "usage: $0 -i <input-fasta> -n <n-mappings> [options]"
    echo "options:"
    echo "   [wfmash]"
    echo "    -i, --input-fasta FILE      input FASTA/FASTQ file"
    echo "    -s, --segment-length N      segment length for mapping [default: 10k]"
    echo "    -l, --block-length N        minimum block length filter for mapping [default: 5*segment-length]"
    echo "    -N, --no-split              disable splitting of input sequences during mapping [enabled by default]"
    echo "    -p, --map-pct-id PCT        percent identity for mapping/alignment [default: 95]"
    echo "    -n, --n-mappings N          number of mappings to retain for each segment"
    echo "    -K, --mash-kmer N           kmer size for mapping [default: 19]"
    echo "    -Y, --exclude-delim C       skip mappings between sequences with the same name prefix before"
    echo "                                the given delimiter character [default: all-vs-all and !self]"
    echo "   [seqwish]"
    echo "    -k, --min-match-len N       filter exact matches below this length [default: 47]"
    echo "    -f, --sparse-factor N       keep this randomly selected fraction of input matches [default: no sparsification]"
    echo "    -B, --transclose-batch      number of bp to use for transitive closure batch [default: 10000000]"
    echo "   [smoothxg]"
    echo "    -H, --n-haps N              number of haplotypes, if different than that set with -n [default: n-mappings]"
    echo "    -j, --path-jump-max         maximum path jump to include in block [default: 0]"
    echo "    -e, --edge-jump-max N       maximum edge jump before breaking [default: 0 / off]"
    echo "    -G, --poa-length-target N,M target sequence length for POA, first pass = N, second pass = M [default: 4001,4507]"
    echo "    -P, --poa-params PARAMS     score parameters for POA in the form of match,mismatch,gap1,ext1,gap2,ext2"
    echo "                                may also be given as presets: asm5, asm10, asm15, asm20"
    echo "                                [default: 1,19,39,3,81,1 = asm5]"
    echo "    -O, --poa-padding N         pad each end of each sequence in POA with N*(longest_poa_seq) bp [default: 0.03]"
    echo "    -d, --pad-max-depth N       depth/haplotype at which we don't pad the POA problem [default: 100]"
    echo "    -M, --write-maf             write MAF output representing merged POA blocks [default: off]"
    echo "    -Q, --consensus-prefix P    use this prefix for consensus path names [default: Consensus_]"
    #echo "    -C, --consensus-spec SPEC   consensus graph specification: write consensus graphs to"
    #echo "                                BASENAME.cons_[spec].gfa; where each spec contains at least a min_len parameter"
    #echo "                                (which defines the length of divergences from consensus paths to preserve in the"
    #echo "                                output), optionally a file containing reference paths to preserve in the output,"
    #echo "                                a flag (y/n) indicating whether we should also use the POA consensus paths, a"
    #echo "                                minimum coverage of consensus paths to retain (min_cov), and a maximum allele"
    #echo "                                length (max_len, defaults to 1e6); implies -a; example:"
    #echo "                                cons,100,1000:refs1.txt:n,1000:refs2.txt:y:2.3:1000000,10000"
    #echo "                                [default: off]"
    echo "   [odgi]"
    echo "    -v, --skip-viz              don't render visualizations of the graph in 1D and 2D [default: make them]"
    echo "    -S, --stats                 generate statistics of the seqwish and smoothxg graph [default: off]"
    echo "   [vg]"
    echo "    -V, --vcf-spec SPEC         specify a set of VCFs to produce with SPEC = REF:DELIM[,REF:DELIM]*"
    echo "                                the paths matching ^REF are used as a reference, while the sample haplotypes"
    echo "                                are derived from path names, e.g. when DELIM=# and with '-V chm13:#',"
    echo "                                a path named HG002#1#ctg would be assigned to sample HG002 phase 1 [default: off]"
    echo "   [multiqc]"
    echo "    -m, --multiqc               generate MultiQC report of graphs' statistics and visualizations,"
    echo "                                automatically runs odgi stats [default: off]"
    echo "   [general]"
    echo "    -o, --output-dir PATH       output directory"
    echo "    -a, --input-paf FILE        input PAF file; the wfmash alignment step is skipped"
    echo "    -r, --resume PATH           do not overwrite existing output from wfmash, seqwish, smoothxg in given directory"
    echo "                                [default: start pipeline from scratch in a new directory]"
    echo "    -t, --threads N             number of compute threads to use in parallel steps"
    echo "    -T, --poa-threads N         number of compute threads to use during POA (set lower if you OOM during smoothing)"
    echo "    -A, --keep-temp-files       keep intermediate graphs"
    echo "    -Z, --compress              compress alignment (.paf), graph (.gfa, .og), and MSA (.maf) outputs with pigz,"
    echo "                                and variant (.vcf) outputs with bgzip"
    echo "    -h, --help                  this text"
    echo
    echo "Use wfmash, seqwish, smoothxg, and odgi to build and display a pangenome graph."
    exit
fi

# Alignment
mapper_letter='W'
n_mappings_minus_1=$( echo "$n_mappings - 1" | bc )
paf_spec=$mapper_letter-s$segment_length-l$block_length-p$map_pct_id-n$n_mappings_minus_1-K$mash_kmer

if [[ $no_merge_segments == true ]];
then
    merge_cmd=-M
    paf_spec="$paf_spec"-M
fi

if [[ $no_splits == true ]];
then
    split_cmd=-N
    paf_spec="$paf_spec"-N
fi

prefix_paf="$input_fasta".$(echo $paf_spec | sha256sum | head -c 7)
if [[ $input_paf != false ]];
then
    prefix_paf=$input_paf
fi

if [[ $exclude_delim != false ]];
then
    exclude_cmd="-Y "$exclude_delim
else
    exclude_cmd=-X
fi

# Graph induction
prefix_seqwish="$prefix_paf".$(echo k$min_match_length-f$sparse_factor-B$transclose_batch | sha256sum | head -c 7)

# Normalization
block_id_min=$(echo "scale=4; $map_pct_id / 100.0" | bc)
prefix_smoothed="$prefix_seqwish".$(echo h$n_haps-G$target_poa_length-j$max_path_jump-e$max_edge_jump-d$pad_max_depth-I$block_id_min-R$block_ratio_min-p$poa_params-O$poa_padding | sha256sum | head -c 7)


fmt="%C\n%Us user %Ss system %P cpu %es total %MKb max memory"
timer=$(which time)

if [[ "$output_dir" != "" ]]; then
	if [ ! -e "$output_dir" ]; then
		mkdir "$output_dir"
	fi
	prefix_paf="$output_dir"/$(basename "$prefix_paf")
	prefix_seqwish="$output_dir"/$(basename "$prefix_seqwish")
	prefix_smoothed="$output_dir"/$(basename "$prefix_smoothed")
fi

date=`date "+%m-%d-%Y_%H:%M:%S"`
log_file="$prefix_smoothed".$date.log
param_file="$prefix_smoothed".$date.params.yml

# write parameters to log_file:
echo -e "Starting pggb on `date`\n" > "$log_file"
echo -e "Command: $cmd\n" >> "$log_file"
echo -e "PARAMETERS\n" >> "$log_file"
cat <<EOT | tee -a "$log_file" "$param_file" >/dev/null
general:
  input-fasta:        $input_fasta
  output-dir:         $output_dir
  resume:             $resume
  compress:           $compress
  threads:            $threads
wfmash:
  mapping-tool:       $mapper
  no-splits:          $no_splits
  segment-length:     $segment_length
  block-length:       $block_length
  no-merge-segments:  $no_merge_segments
  map-pct-id:         $map_pct_id
  n-mappings:         $n_mappings
  mash-kmer:          $mash_kmer
  exclude-delim:      $exclude_delim
seqwish:
  min-match-len:      $min_match_length
  sparse-factor:      $sparse_factor
  transclose-batch:   $transclose_batch
smoothxg:
  n-haps:             $n_haps
  block_id_min:       $block_id_min
  path-jump-max:      $max_path_jump
  edge-jump-max:      $max_edge_jump
  poa-length-target:  $target_poa_length
  poa-params:         $poa_params_cmd
  write-maf:          $write_maf
  consensus-prefix:   $consensus_prefix
  consensus-spec:     $consensus_spec
  pad-max-depth:      $pad_max_depth
  block-id-min:       $block_id_min
  block-ratio-min:    $block_ratio_min
  poa_threads:        $poa_threads
  poa_padding:        $poa_padding
odgi:
  viz:                $do_viz
  layout:             $do_layout
  stats:              $do_stats
gfaffix:
  normalize:          $normalize
vg:
  deconstruct:        $vcf_spec  
reporting:
  multiqc:            $multiqc
EOT

echo -e "\nRunning pggb\n" >> "$log_file"
if [[ ! -s $prefix_paf.$mapper.paf || $resume == false ]]; then
  if [[ "$mapper" == "wfmash" ]];
  then
          ($timer -f "$fmt" wfmash \
              $exclude_cmd \
              -s $segment_length \
              $block_length_cmd \
              $merge_cmd \
              $split_cmd \
              $mash_kmer_cmd \
              -p $map_pct_id \
              -n $n_mappings_minus_1 \
              -t $threads \
              "$input_fasta" "$input_fasta" \
              > "$prefix_paf".$mapper.paf) 2> >(tee -a "$log_file")
  fi
fi

# correct prefix_paf and ajdust input here accordingly
seqwish_paf="$prefix_paf".$mapper.paf
if [[ $input_paf != false ]];
then
    seqwish_paf=$input_paf
fi
if [[ ! -s $prefix_seqwish.seqwish.gfa || $resume == false ]]; then
    $timer -f "$fmt" seqwish \
        -t $threads \
        -s "$input_fasta" \
        -p "$seqwish_paf" \
        -k $min_match_length \
        -f $sparse_factor \
        -g "$prefix_seqwish".seqwish.gfa \
        -B $transclose_batch \
        -P \
        2> >(tee -a "$log_file")
fi

if [[ $consensus_spec != false ]]; then
    # for merging consensus (currently problematic) we should add "-M -J 1 -G 150" here
    consensus_params="-C "$prefix_smoothed".cons,"$consensus_spec
else
    consensus_params="-V"
fi

if [[ $write_maf != false ]]; then
    maf_params="-m "$prefix_smoothed".smooth.maf"
fi

# how many times will we smooth?
smooth_iterations=$(echo $target_poa_length | tr ',' '\n' | wc -l)

keep_temp=""
if [[ $keep_intermediate_files == true ]]; then
    keep_temp="-K"
fi

for i in $(seq 1 $smooth_iterations);
do
    input_gfa="$prefix_seqwish".seqwish.gfa
    if [[ $i != 1 ]]; then
        input_gfa="$prefix_smoothed".smooth.$(echo $i - 1 | bc).gfa
    fi
    if [[ $i != $smooth_iterations ]]; then
        if [[ ! -s $prefix_smoothed.smooth.$i.gfa || $resume == false ]]; then
            poa_length=$(echo $target_poa_length | cut -f $i -d, )
            $timer -f "$fmt" smoothxg \
                   -t $threads \
                   -T $poa_threads \
                   -g "$input_gfa" \
                   -w $(echo "$poa_length * $n_haps" | bc) \
                   $keep_temp \
                   -X 100 \
                   -I $block_id_min \
                   -R $block_ratio_min \
                   -j $max_path_jump \
                   -e $max_edge_jump \
                   -l $poa_length \
                   $poa_params_cmd \
                   -O $poa_padding \
                   -Y $(echo "$pad_max_depth * $n_haps" | bc) \
                   -d 0 -D 0 \
                   -V \
                   -o "$prefix_smoothed".smooth.$i.gfa \
                   2> >(tee -a "$log_file")
        fi
    else
        if [[ ! -s $prefix_smoothed.smooth.gfa || $resume == false ]]; then
            poa_length=$(echo $target_poa_length | cut -f $i -d, )
            $timer -f "$fmt" smoothxg \
                   -t $threads \
                   -T $poa_threads \
                   -g "$input_gfa" \
                   -w $(echo "$poa_length * $n_haps" | bc) \
                   $keep_temp \
                   -X 100 \
                   -I $block_id_min \
                   -R $block_ratio_min \
                   -j $max_path_jump \
                   -e $max_edge_jump \
                   -l $poa_length \
                   $poa_params_cmd \
                   -O $poa_padding \
                   -Y $(echo "$pad_max_depth * $n_haps" | bc) \
                   -d 0 -D 0 \
                   $maf_params \
                   -Q $consensus_prefix \
                   $consensus_params \
                   -o "$prefix_smoothed".smooth.gfa \
                   2> >(tee -a "$log_file")
        fi
    fi
done

prefix_final_graph="$prefix_smoothed".smooth.final
if [[ $normalize == true ]];
then
    # Remove redundancy and sort
    ( $timer -f "$fmt" gfaffix "$prefix_smoothed".smooth.gfa -o "$prefix_final_graph".gfa | $timer -f "$fmt" pigz >"$prefix_final_graph".affixes.tsv.gz ) 2> >(tee -a "$log_file")
    ( $timer -f "$fmt" odgi build -t $threads -P -g "$prefix_final_graph".gfa -o - -O \
        | $timer -f "$fmt" odgi sort -P -p Ygs -t $threads -i - -o "$prefix_final_graph".og ) 2> >(tee -a "$log_file")
    ( $timer -f "$fmt" odgi view -i "$prefix_final_graph".og -g >"$prefix_final_graph".gfa ) 2> >(tee -a "$log_file")
else
    mv "$prefix_smoothed".smooth.gfa "$prefix_final_graph".gfa
    $timer -f "$fmt" odgi build -t $threads -P -g "$prefix_final_graph".gfa -o "$prefix_final_graph".og 2> >(tee -a "$log_file")
fi

if [[ $multiqc == true ]];
then
    do_stats=true
fi

if [[ $do_stats == true ]];
then
    $timer -f "$fmt" odgi build -t $threads -P -g "$prefix_seqwish".seqwish.gfa -o "$prefix_seqwish".seqwish.og 2> >(tee -a "$log_file")
    odgi stats -i "$prefix_seqwish".seqwish.og -m > "$prefix_seqwish".seqwish.og.stats.yaml 2>&1 | tee -a "$log_file"
    odgi stats -i "$prefix_final_graph".og -m  > "$prefix_final_graph".og.stats.yaml 2>&1 | tee -a "$log_file"
    if [[ $consensus_spec != false ]]; then
        for consensus_graph in "$prefix_smoothed"*.cons*.gfa; do
            odgi build -t $threads -P -g "$consensus_graph" -o "$consensus_graph".og 2> >(tee -a "$log_file")
            odgi stats -i "$consensus_graph".og -m >"$consensus_graph".og.stats.yaml 2>&1 | tee -a "$log_file"
        done
    fi
fi

if [[ $do_viz == true ]];
then
    # big problem: this assumes that there is no "Consensus_" in the input sequences
    $timer -f "$fmt" odgi viz -i "$prefix_final_graph".og \
                    -o "$prefix_final_graph".og.viz_multiqc.png \
                    -x 1500 -y 500 -a 10 -I $consensus_prefix \
                    2> >(tee -a "$log_file")
    $timer -f "$fmt" odgi viz -i "$prefix_final_graph".og \
                    -o "$prefix_final_graph".og.viz_pos_multiqc.png \
                    -x 1500 -y 500 -a 10 -u -d -I $consensus_prefix \
                    2> >(tee -a "$log_file")
    $timer -f "$fmt" odgi viz -i "$prefix_final_graph".og \
                    -o "$prefix_final_graph".og.viz_depth_multiqc.png \
                    -x 1500 -y 500 -a 10 -m -I $consensus_prefix \
                    2> >(tee -a "$log_file")
    $timer -f "$fmt" odgi viz -i "$prefix_final_graph".og \
                    -o "$prefix_final_graph".og.viz_inv_multiqc.png \
                    -x 1500 -y 500 -a 10 -z -I $consensus_prefix \
                    2> >(tee -a "$log_file")
fi

if [[ $do_layout == true ]];
then
    # the 2D layout is "smoother" when we chop the nodes of the graph to a fixed maximum length
    #$timer -f "$fmt" odgi chop -i "$prefix_final_graph".og -c 100 -o ""$prefix_final_graph".chop.og \
    #    2> >(tee -a "$log_file")

    # adding `-N g` to this call can help when rendering large, complex graphs that aren't globally linear
    $timer -f "$fmt" odgi layout -i "$prefix_final_graph".og \
                       -o "$prefix_final_graph".og.lay \
                       -T "$prefix_final_graph".og.lay.tsv \
                       -t $threads -P \
                       2> >(tee -a "$log_file")

    # this can be configured to draw the graph in different ways, based on the same layout
    # here we draw in default mode
    $timer -f "$fmt" odgi draw -i "$prefix_final_graph".og \
                     -c "$prefix_final_graph".og.lay \
                     -p "$prefix_final_graph".og.lay.draw.png \
                     -H 1000 \
                     2> >(tee -a "$log_file")
    # this attempts to add paths
    $timer -f "$fmt" odgi draw -i "$prefix_final_graph".og \
                     -c "$prefix_final_graph".og.lay \
                     -p "$prefix_final_graph".og.lay.draw_multiqc.png \
                     -C -w 20 \
                     -H 1000 \
                     2> >(tee -a "$log_file")
fi

if [[ $vcf_spec != false ]];
then
    for s in $( echo "$vcf_spec" | tr ',' ' ' );
    do
        ref=$(echo "$s" | cut -f 1 -d: )
        delim=$(echo "$s" | cut -f 2 -d: )
        echo "[vg::deconstruct] making VCF with reference=$ref and delim=$delim"
        vcf="$prefix_final_graph".$(echo $ref | tr '/|' '_').vcf
        ( TEMPDIR=$(pwd) $timer -f "$fmt" vg deconstruct -P $ref \
                 -H $delim -e -a -t $threads "$prefix_final_graph".gfa >"$vcf" ) 2> >(tee -a "$log_file")
        bcftools stats "$vcf" > "$vcf".stats
    done
fi

multiqc_out_dir=$(dirname "$input_fasta")

multiqc_config="# Report section config for nice titles and descriptions
custom_data:
  odgi_viz:
    section_name: ODGI 1D visualization
    description: This image shows a 1D rendering of the built pangenome graph. The graph nodes are arranged from left to right, forming the pangenome sequence. The colored bars represent the paths versus the pangenome sequence in a binary matrix. The path names are placed on the left. The black lines under the paths are the links, which represent the graph topology.
  odgi_viz_pos:
    section_name: ODGI 1D visualization by path position
    description: This shows a 1D rendering of the built pangenome graph where the paths are colored according to their nucleotide position. Light grey means a low path position, black is the highest path position.
  odgi_viz_inv:
    section_name: ODGI 1D visualization by path orientation
    description: This image shows a 1D rendering of the built pangenome graph where the paths are colored by orientation. Forward is black, reverse is red.
  odgi_viz_depth:
    section_name: ODGI 1D visualization by node depth
    description: This shows a 1D rendering of the built pangenome graph where the paths are colored according to path depth. Using the Spectra color palette with 4 levels of path depths, white indicates no depth, while grey, red, and yellow indicate depth 1, 2, and greater than or equal to 3, respectively.
  odgi_draw:
    section_name: ODGI 2D drawing
    description: This image shows a 2D rendering of the built pangenome graph.    

# Custom search patterns to find the image outputs
sp:
  odgi_draw:
    fn: \"*draw_multiqc.png\"
  odgi_viz:
    fn: \"*viz_multiqc.png\"
  odgi_viz_pos:
    fn: \"*viz_pos_multiqc.png\"
  odgi_viz_inv:
    fn: \"*viz_inv_multiqc.png\"
  odgi_viz_depth:
    fn: \"*viz_depth_multiqc.png\"
  testing_name:
    fn: \"*draw.png\"    
ignore_images: false

# Make the custom content stuff come after the ODGI module output
module_order:
  - odgi
  - custom_content

# Set the order that the custom content plots should come in
custom_content:
  order:
    - odgi_viz
    - odgi_viz_pos
    - odgi_viz_inv
    - odgi_viz_depth
    - odgi_draw
fn_clean_exts:
  - \".gfa\""

if [[ $multiqc == true ]];
then
    echo "$multiqc_config" > "$output_dir"/multiqc_config.yaml
    if [[ $output_dir == "" ]];
    then
        $timer -f "$fmt" multiqc "$multiqc_out_dir" \
        -s \
        -o "$multiqc_out_dir" \
        -c "$output_dir"/multiqc_config.yaml \
        2> >(tee -a "$log_file")
    else
        $timer -f "$fmt" multiqc "$output_dir" \
        -s \
        -o "$output_dir" \
        -c "$output_dir"/multiqc_config.yaml \
        2> >(tee -a "$log_file")
    fi
fi

if [[ $keep_intermediate_files != true ]];
then
  rm -f "$prefix_seqwish".seqwish.{gfa,og}
  rm -f "$prefix_smoothed".smooth.[0-9]*.{gfa,og}

  if [[ $normalize == true ]];
  then
    rm -f "$prefix_smoothed".smooth.{gfa,og}
  fi
fi

if [[ $compress == true ]];
then
    if [[ $input_paf == false ]];
    then
        pigz -f -q -p $threads "$prefix_paf"*.paf -v
    fi
    pigz -f -q -p $threads "$prefix_seqwish"*.{gfa,og} -v
    if [[ $write_maf != false ]]; then
      pigz -f -q -p $threads "$prefix_seqwish"*.maf -v
    fi
    if [[ $vcf_spec != false ]]; then
      ls "$prefix_seqwish"*.vcf | while read f; do bgzip $f -f -@ $threads; done
    fi
    if [[ $do_layout == true ]]; then
        pigz -f -q -p $threads "$prefix_final_graph".og.lay.tsv -v
    fi
fi
