#!/usr/bin/env perl
#ABSTRACT: A program to run DADA2 from the CLI
use 5.012;
use warnings;
my $VERSION  = '1.3.1';

BEGIN { 
		
		our $EXIT_QUIET = 0;				# Set to true when some error is found (see END {})
		our $STARTED = 0;					# Set to true when parameters are parsed and pipeline initialized (see END {})
		our $START_TIME = time(); 			# Track starting time
		my ($sec,$min,$hour,$mday,$mon,$year,$wday,$yday,$isdst)=localtime(time);
		our $EXECUTION_TIMESTAMP = sprintf ( "%04d%02d%02d %02d:%02d:%02d",
                                   $year+1900,$mon+1,$mday,$hour,$min,$sec);
}

=pod

=head1 NAME

B<dadaist2> - a shell wrapper for DADA2, to detect representative sequences and
generate a feature table starting from Illumina Paired End reads. 
This is the main program of the I<dadaist2 toolkit> that includes several
wrappers and utilities to streamline the analysis
of metabarcoding reads from the Linux shell to R.

=head1 AUTHOR

Andrea Telatin <andrea.telatin@quadram.ac.uk>

=head1 SYNOPSIS

  dadaist2 [options] -i INPUT_DIR -o OUTPUT_DIR

=cut

# Minichangelog
# 1.0.9  - Seqfu qual, improved cutadapt
# 1.0.0  - Accepted release, documentation updates, GitHub actions 
# 0.9.0  - Improved support tools (addTaxToFasta..), bug fixes
# 0.8.0  - `--dada-pool`
# 0.7.7  - MultiQC, Rhea support, bugfixes...
# 0.7.5  - Crosstalk
# 0.7.3  - Bugfixes
# 0.7.2  - No join
# 0.7.0  - Decipher support
# 0.6.0  - Phyloseq
# 0.5.0  - support for MicrobiomeAnalyst (enabled by default)
# 0.4.0  - Improved verbose mode, added quality parameter, added support to save RDS, added check for dada2 installed, expanded test set, added dadaist2-metadata
# 0.3.00 - Added a reference database downloader (dadaist2-getdb)
# 0.2.00 - Manual, save Rds, further analysis, taxonomy labels
# 0.1.09 - Added autocalculated trimming [qualified_wnd/qualified_min_qual]
# 0.1.07 - Added quality plot

use Getopt::Long;
use Data::Dumper;
use FindBin qw($RealBin);
use File::Basename;
use File::Spec;
use File::Spec::Functions;
use File::Copy;
use Term::ANSIColor qw(color);
use FASTX::Reader;
use FASTX::ScriptHelper;
use JSON::PP;
use Digest::MD5 qw(md5_hex);
use File::Temp;
use Pod::Usage;

my $PROGRAM           = basename($0);

# Log templating
my  $unfilled         = '<p>Not available.</p>';
our $TEMPLATE         = loadTemplate();
our $LOG_TEMPLATE     = getSections($TEMPLATE);

addToTemplate($LOG_TEMPLATE, 'VERSION', $VERSION);
addToTemplate($LOG_TEMPLATE, 'CPU_CORES', getCpu());
addToTemplate($LOG_TEMPLATE, 'RAM_GB', getRamGb());
addToTemplate($LOG_TEMPLATE, 'UNAME', getSysinfo());
 
# Where to place the temporary directory
my $opt_temp_dir     = File::Spec->tmpdir();

my $seqfu_minqual	 = 20;
my $seqfu_wndqual	 = 28.5;
my $opt_verbose      = 1;
my $opt_quiet        = 0;
my $opt_debug        = 0;
my $opt_metadata     = undef;
my $opt_id_separator = '_';	  # [SampleName]_S01_L001_R1_001...
my $opt_for_tag      = '_R1'; # Tag for R1
my $opt_rev_tag      = '_R2'; # Tag for R2
my $opt_seq_prefix   = 'ASV'; # Name for FASTA sequences (MD5 to use seq hash)
my $opt_threads      = num_cpu();			# Default threads = max threads else user input (-t)
my $opt_tax_db       = "skip";
my $decipher_db      = undef;
my $opt_save_rds     = 0;
my $save_rds         = "no";
my $opt_metadata_separator = "\t";
my $opt_notify       = 0;
my $opt_crosstalk    = 0;
my $opt_seqfu        = undef;
my $opt_skip_tree    = 0;

# FASTP Preprocess
# To be autocalculated:
my $opt_notrim;
my $opt_fastp_qc;
my $fastp_trim_front1 = 0;
my $fastp_trim_front2 = 0;
my $fastp_trim_tail1  = 0;
my $fastp_trim_tail2  = 0;
my $fastp_minlen      = 110;
my $opt_skip_qc       = undef;

# AmpliCheck
my $opt_primers       = undef;
my $primer_for        = undef;
my $primer_rev        = undef;

# DADA2
my $do_plots         = 'do_plots';
my $trunc_len_1      = 0;
my $trunc_len_2      = 0;
my $opt_trunc_len_1;
my $opt_trunc_len_2;
my $trunc_start_1    = 0;
my $trunc_start_2    = 0;
my $opt_maxee_1      = 1.0;
my $opt_maxee_2      = 1.5;
my $total_len_1      = 0;
my $total_len_2      = 0;
my $qualified_wnd    = 4;
my $qualified_min_qual = 28;
my $opt_no_trunc;
my $opt_justconcat   = 0;
my $opt_process_pool = 0;
my $opt_maxloss      = 0.2;
my @taxonomy = ();
my $r1_len = 10_000;
my $r2_len = 10_000;
my $tot_reads = 0;
my $opt_trunc_q = 10;
my ($opt_help,  $opt_force,$opt_log_filename, $opt_input_directory,
 $opt_output_directory, $opt_version, $opt_skip_plots, $opt_skip_check_deps);

my $chooser = {
	'filtering' => 'skip',
	'database'  => undef,
};

# Check dependencies
my %dependencies = (
    'fastp'   => 'fastp --version 2>&1',
    'fu-primers' => 'fu-primers --version 2>&1',
    'RScript' => 'Rscript  --version 2>&1',
    'dada2 (lib)'  => 'echo "library(dada2)" | R --no-save > /dev/null',

    'exporter' => File::Spec->catfile($RealBin, 'dadaist2-exporter') . ' --version',
);


 
my $parsed = GetOptions(
# Main parameters
 'i|input-directory=s'             => \$opt_input_directory, #
 'o|output-directory=s'            => \$opt_output_directory, #
 'r|d|database=s'                  => \$opt_tax_db, #
 'm|metadata=s'                    => \$opt_metadata, #
 'l|log-file=s'                    => \$opt_log_filename, #
 't|threads=i'                     => \$opt_threads, #
 'tmp-dir=s'                       => \$opt_temp_dir, #
 'force'                           => \$opt_force, #
# Input reads
 '1|for-tag=s'                     => \$opt_for_tag, #
 '2|rev-tag=s'                     => \$opt_rev_tag, #
 's|id-separator=s'                => \$opt_id_separator, #
# ASV generation
 'p|prefix=s'                      => \$opt_seq_prefix, #
 'q|min-quality=f'                 => \$qualified_min_qual, #
 'trunc-qual=f'                    => \$opt_trunc_q, #
 'save-rds'                        => \$opt_save_rds, #
 's1|trim-primer-for=i'            => \$fastp_trim_front1, #
 's2|trim-primer-rev=i'            => \$fastp_trim_front2, #
 'no-trim'                         => \$opt_notrim, #
 'trunc-len-1=i'                   => \$opt_trunc_len_1,
 'trunc-len-2=i'                   => \$opt_trunc_len_2,
 'j|its|just-concat'               => \$opt_justconcat,
 'dada-pool'                       => \$opt_process_pool,
 'max-loss=f'                      => \$opt_maxloss, # 
 'crosstalk'                       => \$opt_crosstalk, #
 'primers=s'                       => \$opt_primers, #
 'maxee1=f'                        => \$opt_maxee_1, #
 'maxee2=f'                        => \$opt_maxee_2, #
 'no-trunc'                        => \$opt_no_trunc, #
 'skip-plots'                      => \$opt_skip_plots,
 'skip-qc'                         => \$opt_skip_qc,
 'fastp'                           => \$opt_fastp_qc,
 'skip-tree'                       => \$opt_skip_tree,
# Misc
 'no-check'                        => \$opt_skip_check_deps, #
 'quiet'                           => \$opt_quiet,
 'verbose'                         => \$opt_verbose,
 'debug'                           => \$opt_debug,
 'popup'                           => \$opt_notify,
 'version'                         => \$opt_version,
 'h|help'                          => \$opt_help,
) || die " Parameters error, type --help for more info.\n";

die " ERROR: Unexpected parameters: ", join(", ", @ARGV), "\n" if ($ARGV[0]);

if (not $opt_skip_tree) {
	$dependencies{'fasttree'} = 'fasttree 2>&1 | head -n 1 | cut -c 11-';
    $dependencies{'clustalo'} = 'clustalo --version 2>&1';
}  

=pod

=head1 PARAMETERS

=head2 Main Parameters

=over 4

=item I<-i>, I<--input-directory> DIRECTORY

Directory containing the paired end files in FASTQ format, gzipped or not.

=item I<-o>, I<--output-directory> DIRECTORY

Output directory (will be created). 
It is recommended recomment using a new directory for each run.

=item I<-d>, I<--database> DATABASE

Reference database in gzipped FASTA format. 
Optional (default: skip) but highly recommended.

=item I<-m>, I<--metadata> FILE

Metadata file in TSV format, first column must match sample IDs. 
If not supplied a template will be autogenerated using C<dadaist2-metadata>.

=item I<-t>, I<--threads> INT

Number of threads (default: 2)

=item I<--primers> FOR:REV

Strip primers with I<cutadapt>, supply both sequences separated by a colon.

=item I<-j>, I<--just-concat>

Do not try merging paired end reads, just concatenate. 

=item I<--fastp>

Perform the legacy "fastp" QC

=item I<--no-trim>

Do not trim primers (using fastp). 
Equivalent to C<--trim-primer-for 0 --trim-primer-rev 0>.

=item I<--force>

Will overwrite output folder if it already exists, and will
attempt to produce MicrobiomeAnalyst and Rhea folders even when
DADA2 filters too many reads.

=item I<--dada-pool>

Pool samples in DADA2 analysis (experimental)

=back

=head2 Input reads

We recommend to prepare a polished directory of input reads having
filenames like Samplename_R1.fastq.gz and Samplename_R2.fastq.gz.

Filename starting by numbers are not accepted.

=over 4

=item I<-1>, I<--for-tag> TAG

Tag to recognize forward reads (default: _R1)

=item I<-2>, I<--rev-tag> TAG

Tag to recognize reverse reads (default: _R2)

=item I<-s>, I<--id-separator> CHAR

Character separating the sample name from the rest of the filename (default: _)

=back 

=head2 Metabarcoding processing

=over 4

=item I<--trunc-len-1> and I<--trunc-len-2> INT

Position at which truncate reads (forward and reverse, respectively).

=item I<-q>, I<--min-qual> FLOAT

Minimum average quality for DADA2 truncation (default: 28)

=item I<--no-trunc>

Do not truncate reads at the end
 (required for non-overlapping amplicons, like ITS)


=item I<--maxee1>, and I<--maxee2> FLOAT

Maximum Expected Errors in R1 and R2, respectively (default: 1.0 and 1.5)

=item I<--trunc-qual> FLOAT

DADA2 truncate quality (default: 10)

=item I<-s1>, I<--trim-primer-for> INT

Trim primer from R1 read specifying the number of bases. Similarly
use C<-s2> (C<--trim-primer-rev>) to remove the front bases from the
reverse pair (R2). Default: 20 bases each side.

=item I<--save-rds>

Save a copy of the RDS file (default: off)

=item I<--max-loss> FLOAT

After DADA2 run, check the amount of reads globally remaining from input
to non-chimeric, abort if the ratio is below threshold (default: 0.2)

=back

=head2 Other parameters

=over 4

=item I<--crosstalk> 

Remove crosstalk using the UNCROSS2 algorithm 
as described here L<https://doi.org/10.1101/400762>.

=item I<-p>, I<--prefix> STRING

Prefix for the output FASTA file, if "MD5" is specified, the sequence MD5 hash
will be used instead. Default is "ASV".

=item I<-l>, I<--log-file> FILE

Filename for the program log.

=item I<--tmp-dir> DIR

Where to place the temporary directory (default are system temp dir
or C<$TMPDIR>).


=item I<--skip-tree>

Do not generate tree. I<Experimental|Not recommended>.

=item I<--skip-plots>

Do not generate quality plots.

=item I<--popup>

Display popup notifications (tested on MacOS and Ubuntu)

=item I<--quiet>

Reduce verbosity

=item I<--verbose> and I<--debug>

Increase reported information

=back


=cut
$opt_verbose = 0 if ($opt_quiet);
$opt_version && version();
$opt_help    && pod2usage({-exitval => 0, -verbose => 2});
$do_plots    = 'skip' if ($opt_skip_plots);

if ($opt_notrim) {
	$fastp_trim_front1 = 0;
	$fastp_trim_front2 = 0;
}
$opt_verbose = 1 if ($opt_debug);



# Check essential parameters
my $fastp_threads = $opt_threads > 8 ? 8 : $opt_threads;

usage("Missing parameter(s): input directory (-i DIR).")  if (not defined $opt_input_directory);
usage("Missing parameter(s): output directory (-o DIR).") if (not defined $opt_output_directory);

# [*] Splash screen - here starts the user communication
my $ascii_art_dadaist = <<'EOF';
    ____            __      _      __ ___ 
   / __ \____ _____/ /___ _(_)____/ /|__ \
  / / / / __ `/ __  / __ `/ / ___/ __/_/ /
 / /_/ / /_/ / /_/ / /_/ / (__  ) /_/ __/ 
/_____/\__,_/\__,_/\__,_/_/____/\__/____/ 
                                          
EOF
say STDERR color('bold magenta'), $ascii_art_dadaist, color('reset'), $VERSION
 if not $opt_quiet;



# [-] Check taxonomy database
if ( $opt_tax_db ne 'skip' and ! -e "$opt_tax_db") {
  die " FATAL ERROR: Taxonomy database not found: $opt_tax_db\n"
}

# Prepare temporary directories
my $temp_dir = File::Temp->newdir( 'dadaist2_XXXXXX',
		CLEANUP => 0,
		DIR => $opt_temp_dir);
my $dada2_temp = File::Spec->catdir("$temp_dir", 'dada2');

# Prepare output directory
$opt_output_directory = File::Spec->rel2abs($opt_output_directory) . "/";
$opt_output_directory .= '/' if ($opt_output_directory !~/\/$/);

if ( ! -d "$opt_output_directory") {
	mkdir "$opt_output_directory"    || die " FATAL ERROR:\n Unable to create directory $opt_output_directory/\n $!\n";
} else {
	if (not directoryIsEmpty($opt_output_directory)) {
		print STDERR "\n[WARNING] Output directory found.\n"; 
		print STDERR " This is a warning but in future releases this might require to specify --force to proceed.\n";
		
	}
}
mkdir File::Spec->catdir("$opt_output_directory", 'qc')    || die " FATAL ERROR:\n Unable to create directory $opt_output_directory/qc\n $!\n";
die "ERROR: Unable to make ", File::Spec->catdir("$opt_output_directory", 'qc'), " directory.\n" unless (-d File::Spec->catdir("$opt_output_directory", 'qc') );



# Log file
my $log_filename =  $opt_log_filename // File::Spec->catfile($opt_output_directory, 'dadaist.log');
my $S = FASTX::ScriptHelper->new({
	verbose    => $opt_verbose,
	logfile    => $log_filename,
});
our $STARTED = 1;
$S->verbose("$PROGRAM $VERSION");
notify_message("Starting dadaist2");
# [-] Check input directory
if ( ! -d "$opt_input_directory") {
  $S->verbose(" FATAL ERROR: Input directory not found: $opt_input_directory");
}
# [***] Remove primers if provided 					[ >>>>>    <<<<< ]
if ($opt_primers) {
	($primer_for, $primer_rev) = split /:/, $opt_primers;
	die " ERROR: Primers should be in the format --primers FORWARD:REVERSE. Semicolon not found.\n" 
	  unless defined $primer_rev;
	if ($opt_seqfu) {
  		$chooser->{'filtering'}  = 'seqfu';
		$S->verbose("Filtering with SeqFu: ON");
	} else {
		$chooser->{'filtering'}  = 'cutadapt';
		$dependencies{'cutadapt'} = 'cutadapt --version';
		$S->verbose("Filtering with Cutadapt: ON");
	}
} elsif ($opt_fastp_qc) {
	$chooser->{'filtering'} = "fastp";
	$S->verbose("Filtering with Fastp: ON");
} elsif ($opt_skip_qc) {
	# TODO - Check primer management here /* warning is temporary */
	$S->verbose(" [WARNING] Primer sequence provided, but skipping QC implies no trimming!");
	$chooser->{'filtering'}  = 'skip';
	$S->verbose("Preparing to skip QC: " . $chooser->{'filtering'}) if ($opt_debug);


}	

if ($opt_save_rds) {
  $save_rds = 'save';	
  $S->verbose("Parameter: --save-rds OK");
}

# If metadata is not supplied, check Metadata maker
# TODO - Switch to "seqfu metadata"
if (! $opt_metadata) {
  $dependencies{'dadaist2-metadata'} = File::Spec->catfile($RealBin, 'dadaist2-metadata') . ' --version';
}

if ($opt_crosstalk ) {
	$dependencies{'crosstalk'} = File::Spec->catfile($RealBin, 'dadaist2-crosstalk') . ' --version';
}

# Check taxonomy database
if ($opt_tax_db ne 'skip') {
		$dependencies{'assign-taxonomy'} = "$RealBin/dadaist2-assigntax --version";
		if ($opt_tax_db =~/RData$/i ) {
			
			# This is a Decipher database
			$dependencies{'DECIPHER'} = 'echo "library(DECIPHER)" | R --no-save > /dev/null';
			
			if (-e $opt_tax_db) {
				$decipher_db = $opt_tax_db;
				$opt_tax_db = 'skip'; # TODO REFACTOR
				$S->verbose("DECIPHER Taxonomy database found: $decipher_db");
			} else {
				$S->verbose("Decipher database not found: $opt_tax_db");
				die " FATAL ERROR: Decipher database not found: $opt_tax_db\n";
			}
			
		} else {
			# This is a DADA2 Naive classifier database
			if (valid_gzipped_fasta($opt_tax_db)) {
				$chooser->{'database'} = 'dada2';
				$dependencies{'Taxonomy'} = File::Spec->catfile($RealBin, 'dadaist2-assigntax') . ' --version';
				$S->verbose("Taxonomy database found: $opt_tax_db");
			} else {
				$S->verbose("Taxonomy database not valid: $opt_tax_db. Taxonomy will be skipped.");
				$opt_tax_db = 'skip';
			}
		}
		$S->verbose("Parameter: taxonomy-type: " . $chooser->{'database'});	
		$S->verbose("Parameter: taxonomy-db: " . $opt_tax_db);
} else {
	$S->verbose("[WARNING] Reference-free denoising is not fully supported.")
}

## Parameters
addToTemplate($LOG_TEMPLATE, 'INPUT_PATH', File::Spec->rel2abs( $opt_input_directory));
addToTemplate($LOG_TEMPLATE, 'OUTPUT_PATH', File::Spec->rel2abs( $opt_output_directory));
addToTemplate($LOG_TEMPLATE, 'PARAMETERS', paramToHtml('Database', File::Spec->rel2abs( $opt_tax_db)), 'append');
addToTemplate($LOG_TEMPLATE, 'PARAMETERS', paramToHtml('Metadata', File::Spec->rel2abs( $opt_metadata)), 'append');
addToTemplate($LOG_TEMPLATE, 'PARAMETERS', paramToHtml('Threads', File::Spec->rel2abs( $opt_threads)), 'append');

if ($opt_verbose) {
	say STDERR " * ", color('reset'), "Input directory: ", color('blue'), $opt_input_directory, color('reset');
	say STDERR " * ", color('reset'), "Output directory: ", color('blue'),  $opt_output_directory, color('reset');	
	say STDERR " * ", color('reset'), "Metadata: ", color('blue'),  $opt_metadata // '<not provided>', color('reset');
	say STDERR " * ", color('reset'), "Reference database: ", color('blue'),  $opt_tax_db, color('reset');
	say STDERR " * ", color('reset'), "Threads: ", color('blue'),  $opt_threads, color('reset');
	say STDERR " * ", color('reset'), "Temporary directory: ", color('blue'),  $temp_dir, color('reset');
	say STDERR " * ", color('reset'), "QC strategy: ", color('blue'),  $chooser->{'filtering'}, color('reset');
}


# Check Seqfu Qual

if ( 
		(! defined $opt_trunc_len_1 or ! defined $opt_trunc_len_2) 
 		and
		( $chooser->{'filtering'} eq 'skip' or $chooser->{'filtering'} eq 'cutadapt' ) 
	) {
	# Seqfu Qual
	my $max_len = 500; #TODO: change with reads length?
	my $cmd_for = qq(seqfu qual --wnd 4 --min-qual $seqfu_minqual --wnd-qual $seqfu_wndqual --max 5000 --maxlen $max_len --skip 3 "$opt_input_directory"/*${opt_for_tag}*);
	my $cmd_rev = qq(seqfu qual --wnd 4 --min-qual $seqfu_minqual --wnd-qual $seqfu_wndqual --max 5000 --maxlen $max_len --skip 3 "$opt_input_directory"/*${opt_rev_tag}*);
	$S->verbose("QC: Checking quality profile with SeqFu");
	my $QC_for = $S->run($cmd_for);
	my $QC_rev = $S->run($cmd_rev); 
	my ($qualified_for, $low_samples_for) = get_qualified_from_seqfu($QC_for->{'stdout'});
	my ($qualified_rev, $low_samples_rev) = get_qualified_from_seqfu($QC_rev->{'stdout'});

	if ($$low_samples_for[0]) {
		$S->verbose("[WARNING] Low quality R1 samples: ", $low_samples_for)
	}
	if ($$low_samples_rev[0]) {
		$S->verbose("[WARNING] Low quality R2 samples: ", $low_samples_rev)
	}
	
	$S->verbose("SeqFu quality truncation at (trunc-len-1 and trunc-len-2): $qualified_for - $qualified_rev");

	$opt_trunc_len_1 = $qualified_for;
	$opt_trunc_len_2 = $qualified_rev; 
}


# Dependency check
$S->verbose("Checking dependencies");
check_required_binaries(\%dependencies) unless ($opt_skip_check_deps);


# Metadata?
my $meta_samples;
if ($opt_metadata) {
	if (not -e "$opt_metadata") {
		$S->verbose("Metadata file (-m) not found: $opt_metadata");
		die "ERROR: Metadata file not found: $opt_metadata\n";
	}
  $meta_samples = check_metadata($opt_metadata);
  copy($opt_metadata,
		 File::Spec->catfile($opt_output_directory,'metadata.tsv')) || die " ERROR:\n Unable to copy metadù file to $opt_output_directory.\n";

} else {
  $opt_metadata = File::Spec->catfile($opt_output_directory, 'metadata.tsv');
  my $meta_cmd = qq($RealBin/dadaist2-metadata -i "$opt_input_directory" -o "$opt_metadata");
  $S->run($meta_cmd, {candie => 1});
  $meta_samples = check_metadata($opt_metadata);
}
$temp_dir = prepare_temporary_directories($temp_dir);

$S->verbose("Threads: $opt_threads");
$S->verbose("Output directory: $opt_output_directory");

my ($files_hash, $is_paired, $samples_count) = get_file_reads($opt_input_directory);

# Add check on metadata-file consistency with input files


for my $file_sample (sort keys %$files_hash) {
   if (defined $meta_samples->{$file_sample}) {
   		$S->verbose("Checked metadata for $file_sample");
   } else {
   		$S->verbose("[WARNING] Sample $file_sample not found in metadata file but found in input directory");
        $S->verbose("Metadata contains: ", join(', ', keys %$meta_samples));
		die;
   }
 
}

my $library_type = $is_paired > 0 ? 'paired-end' : 'single-end';
$S->verbose("Input directory \"$opt_input_directory\": $samples_count found ($library_type) ");

# R output dir
my $opt_R =  File::Spec->catdir($opt_output_directory, 'R');
mkdir "$opt_R";


# Input sanitation with "fastp"
my @failed_samples = ();
my $counter_samples = 0;
my $tot = scalar keys %{ $files_hash };
for my $sample (sort keys %{ $files_hash } ) {
	$counter_samples++;
	$S->verbose("($counter_samples/$tot) Processing $sample: " . $chooser->{'filtering'});
	my $json_file = File::Spec->catfile($temp_dir, "${sample}.json");
	my $filt_cmd;
	my $out;

	if ( $chooser->{'filtering'} eq 'fastp' and not defined $primer_for) {
		# [1.0] Removed max ns 1

		if ($$files_hash{$sample}{rev}) {
			# Paired reads
			$filt_cmd = qq(fastp -w $fastp_threads -i "$$files_hash{$sample}{for}" -I "$$files_hash{$sample}{rev}" -o "$temp_dir"/for/${sample}_R1.fastq.gz -O  "$temp_dir"/rev/${sample}_R2.fastq.gz  ) .
						qq( --trim_front1 $fastp_trim_front1 --trim_front2  $fastp_trim_front2 ).
						qq( --detect_adapter_for_pe --disable_quality_filtering  ).
						qq( --json "$json_file" --html /dev/null ).
						qq( --trim_tail1 $fastp_trim_tail1 --trim_tail2 $fastp_trim_tail2 --length_required $fastp_minlen );
		} else {
			# Single end
			$filt_cmd = qq(fastp -w $fastp_threads -i "$$files_hash{$sample}{for}" -o "$temp_dir"/for/${sample}_R1.fastq.gz   ) .
					qq( --trim_front1 $fastp_trim_front1  ).
					qq( --disable_quality_filtering  ).
					qq( --json "$json_file" --html /dev/null ).
					qq( --trim_tail1 $fastp_trim_tail1 - --length_required $fastp_minlen );
					die " FATAL ERROR: Single end processing is not supported at the moment.\n";
		}

		# Run FASTP
		$out = $S->run( $filt_cmd, { candie => 1 } );
		if ($out->{exit} != 0) {
			die " FATAL ERROR: 'fastp' failed processing sample $sample.\nCommand: $out->{cmd}\nOut: $out->{stdout}\nErr: $out->{stderr}";
		} else {
			addToTemplate($LOG_TEMPLATE, 'TRIMMING_LOG', runToHtml($out), 1);
		}
	# Use Amplicheck (cutadapt)
	} elsif ($chooser->{'filtering'} eq 'cutadapt' and defined $primer_for) {
		my $rev_rc = iupac_rc($primer_rev);
		my $for_rc = iupac_rc($primer_for);

		$filt_cmd = qq(cutadapt -a "$primer_for...$rev_rc" -A "$primer_rev...$for_rc" ).
		qq( --minimum-length $fastp_minlen --discard-untrimmed ) .
		qq( -o "$temp_dir"/for/${sample}_R1.fastq.gz -p "$temp_dir"/rev/${sample}_R2.fastq.gz "$$files_hash{$sample}{for}" "$$files_hash{$sample}{rev}");

		my $cutadapt = $S->run($filt_cmd, { candie => 1});
		if ($cutadapt->{exit} != 0) {
			$S->verbose("Cutadapt error: ", $cutadapt->{stdout}, "\n", $cutadapt->{stderr});
			die;
		}
	} elsif ($chooser->{'filtering'} eq 'seqfu' and defined $primer_for) {
		if ($$files_hash{$sample}{rev}) {
			$filt_cmd = qq(fu-primers -1 "$$files_hash{$sample}{for}" -2 "$$files_hash{$sample}{rev}" ).
									qq( --primer-for "$primer_for" --primer-rev "$primer_rev" | seqfu deinterleave --for-ext "_R1.fastq" --rev-ext "_R2.fastq" -o "$temp_dir"/for/${sample} -);
			$S->verbose("\$ " . $filt_cmd) if ($opt_debug);
			$out = $S->run( $filt_cmd, { candie => 1 } );
			addToTemplate($LOG_TEMPLATE, 'TRIMMING_LOG', runToHtml($out), 1);
			$S->run(qq(gzip "$temp_dir"/for/${sample}*));
			move(
				File::Spec->catfile(
					File::Spec->catdir($temp_dir, 'for'),
					"${sample}_R2.fastq.gz"),
				File::Spec->catdir($temp_dir, 'rev')
			);
		} else {
			die " Single end processing not supported at the moment.\n";
		}
	} elsif ( $chooser->{'filtering'} eq 'skip' ) {
		my $for_dir = File::Spec->catdir($temp_dir, 'for');
		my $rev_dir = File::Spec->catdir($temp_dir, 'rev');
		my $for_name = basename($$files_hash{$sample}{for});
		my $rev_name = basename($$files_hash{$sample}{rev});
		
		if (not defined $opt_trunc_len_1 or not defined $opt_trunc_len_2) {
			$S->verbose("ERROR: Missing boundaries --trunc-len-1 and --trunc-len-2 to operate in QC-less mode.");
			exit;
		}
		$trunc_len_1 = $opt_trunc_len_1;
		$trunc_len_2 = $opt_trunc_len_2;
		$for_name =~s/\.(fq|FQ|FASTQ)/.fastq/;
		$rev_name =~s/\.(fq|FQ|FASTQ)/.fastq/;
		#my $dest_for = File::Spec->catfile($for_dir, $for_name);
		#my $dest_rev = File::Spec->catfile($rev_dir, $rev_name);

		my $dest_for = File::Spec->catfile($for_dir, "${sample}_R1.fastq.gz");
		my $dest_rev = File::Spec->catfile($rev_dir, "${sample}_R2.fastq.gz");

		if ($fastp_trim_front1 > 0 or $fastp_trim_front2 > 0) {
			$S->verbose("Fixed primer trimming");
			my $cmd1 = qq(seqfu cat --trim-front $fastp_trim_front1 "$$files_hash{$sample}{for}" | gzip -c > "$dest_for" );
			my $cmd2 = qq(seqfu cat --trim-front $fastp_trim_front2 "$$files_hash{$sample}{rev}" | gzip -c > "$dest_rev" );
			$S->run($cmd1, { candie => 0 });
			$S->run($cmd2, { candie => 0 });
		} else {
			$S->verbose("Copying input reads for DADA2");
			
			
			copy(
					$$files_hash{$sample}{for},
					$dest_for
				);
			copy(
					$$files_hash{$sample}{rev},
					$dest_rev
				);

		}
		
		if ($for_name !~ /gz$/) {
			$S->run( qq(gzip "$dest_for") );
		}
		if ($rev_name !~ /gz$/) {
			$S->run( qq(gzip "$dest_rev") );
		}
		
	} else {
		$S->verbose("Unexpected setting: filtering should be performed either via trimming (fastp) or primer removal (SeqFu): <". $chooser->{'filtering'} .">\n");
		die;
	}

	

	
  if ( $chooser->{'filtering'} eq 'fastp' ) {
		copy($json_file, 
			File::Spec->catdir("$opt_output_directory", 'qc/')) or 
			die " ERROR:\n Unable to copy JSON stats <$json_file> to <$opt_output_directory>.\n";
		my $summary = load_json_from_file( $json_file );
		my $raw_reads  = $summary->{summary}->{before_filtering}->{total_reads};
		my $pass_reads = $summary->{summary}->{after_filtering}->{total_reads};
		my $ins_size   = $summary->{insert_size}->{peak};
		my $q30        = $summary->{summary}->{after_filtering}->{q30_rate};
		my $passed     = 0;
		my @qual_R1    = @{ $summary->{read1_after_filtering}->{quality_curves}->{mean} };
		my $fastp_qc_pass = 1;
		($trunc_start_1, $trunc_len_1, $total_len_1)  = get_qualified_positions(@qual_R1);
		$fastp_qc_pass = 0 unless (defined $trunc_start_1);

		my @qual_R2;
		if ($$files_hash{$sample}{rev}) {
			@qual_R2    = @{ $summary->{read2_after_filtering}->{quality_curves}->{mean} };
			($trunc_start_2, $trunc_len_2, $total_len_2)  = get_qualified_positions(@qual_R2);
			$fastp_qc_pass = 0 unless (defined $trunc_start_2);
		}
 
		if (not $fastp_qc_pass) {
			$S->verbose("Sample failed QC: " . $sample);
			push(@failed_samples, $sample);
			#TODO
		}
	 

		$tot_reads += int( $pass_reads / 2);
		$r1_len     = $summary->{read1_after_filtering}->{total_cycles}
			if ($summary->{read1_after_filtering}->{total_cycles} < $r1_len);
		$r2_len     = $summary->{read2_after_filtering}->{total_cycles}
			if ($summary->{read2_after_filtering}->{total_cycles} < $r2_len);
		$passed += sprintf("%.2f", 100*$pass_reads/$raw_reads) if ($raw_reads);
		$S->verbose(qq($pass_reads/$raw_reads (${passed}%) reads kept.\n * Average insert size: $ins_size bp\n * Q30: $q30\n * Qualified region: [$trunc_start_1 - $trunc_len_1]/$total_len_1, [$trunc_start_2 - $trunc_len_2]/$total_len_2));

	} elsif ($chooser->{'filtering'} eq 'seqfu' ) {
		$trunc_start_1 = 0;
		$trunc_len_1 = 0;
		$trunc_start_2 = 0;
		$trunc_len_2 = 0;
		$S->verbose($out->{stderr}) if ($out->{stderr});
	}
}

# Failed samples
$S->verbose(scalar @failed_samples . " failed samples: " . join(",", @failed_samples))
  if (scalar @failed_samples);

notify_message("Read trimming finished, starting DADA2.");

my $dada2_min_reads = $tot_reads < 10000 ? $tot_reads : 10000;
my @dada_paired_parameter_labels = (
'forward_reads', 'reverse_reads', 'feature_table_output', 'stats_output',
'filt_forward', 'filt_reverse', 'truncLenF', 'truncLenR',
'trimLeftF', 'trimLeftR', 'maxEEF', 'maxEER',
'truncQ', 'chimeraMethod', 'minFold','threads',
'nreads_learn','baseDir', 'doPlots', 'taxonomyDb', 'saveRDS', 'noMerge', 'processPool'
);


if  ($opt_no_trunc) {
	$S->verbose("Skipping truncation (trunc-len-1 and trunc-len-2 are set to 0)");
	$trunc_len_1 = 0;
	$trunc_len_2 = 0;
}

if ($opt_debug){
	$S->writelog("Counting input reads for DADA2");
	my $forcount = $S->run("seqfu stats " . File::Spec->catfile($temp_dir, 'for', '*.gz'));
	my $revcount = $S->run("seqfu stats " . File::Spec->catfile($temp_dir, 'rev', '*.gz'));
	$S->writelog("Forward: " . $forcount->{stdout});
	$S->writelog("Reverse: " . $revcount->{stdout});

}
my @dada_paired_args = (
# 1) File path to directory with the FORWARD .fastq.gz files to be processed.
#    Ex: path/to/dir/with/FWD_fastqgzs
	"${temp_dir}/for",

# 2) File path to directory with the REVERSE .fastq.gz files to be processed.
#    Ex: path/to/dir/with/REV_fastqgzs
	"${temp_dir}/rev",

# 3) File path to output tsv file. If already exists, will be overwritten.
#    Ex: path/to/output_file.tsv
	"${temp_dir}/dada2/dada2.tsv",

# 4) File path to tracking tsv file. If already exists, will be overwritte.
#    Ex: path/to/tracking_stats.tsv
	"${temp_dir}/dada2/stats.tsv",

# 5) File path to directory to write the filtered FORWARD .fastq.gz files. These files are intermediate
#               for the full workflow. Currently they remain after the script finishes. Directory must
#               already exist.
#    Ex: path/to/dir/with/FWD_fastqgzs/filtered
	"${temp_dir}/for/filtered",

# 6) File path to directory to write the filtered REVERSE .fastq.gz files. These files are intermediate
#               for the full workflow. Currently they remain after the script finishes. Directory must
#               already exist.
#    Ex: path/to/dir/with/REV_fastqgzs/filtered
	"${temp_dir}/rev/filtered",


### FILTERING ARGUMENTS ###
#
# 7) truncLenF - The position at which to truncate forward reads. Forward reads shorter
#               than truncLenF will be discarded.
#               Special values: 0 - no truncation or length filtering.
#    Ex: 240
    $trunc_len_1,

# 8) truncLenR - The position at which to truncate reverse reads. Reverse reads shorter
#               than truncLenR will be discarded.
#               Special values: 0 - no truncation or length filtering.
#    Ex: 160
    $trunc_len_2,

# 9) trimLeftF - The number of nucleotides to remove from the start of
#               each forward read. Should be less than truncLenF.
#    Ex: 0
    $trunc_start_1,

# 10) trimLeftR - The number of nucleotides to remove from the start of
#               each reverse read. Should be less than truncLenR.
#    Ex: 0
    $trunc_start_2,

# 11) maxEEF - Forward reads with expected errors higher than maxEEF are discarded.
#               Both forward and reverse reads are independently tested.
#    Ex: 2.0
    $opt_maxee_1,

# 12) maxEER - Reverse reads with expected errors higher than maxEER are discarded.
#               Both forward and reverse reads are independently tested.
#    Ex: 2.0
    $opt_maxee_2,

# 13) truncQ - Reads are truncated at the first instance of quality score truncQ.
#                If the read is then shorter than truncLen, it is discarded.
#    Ex: 2
    $opt_trunc_q,
### CHIMERA ARGUMENTS ###
#
# 14) chimeraMethod - The method used to remove chimeras. Valid options are:
#               none: No chimera removal is performed.
#               pooled: All reads are pooled prior to chimera detection.
#               consensus: Chimeras are detect in samples individually, and a consensus decision
#                           is made for each sequence variant.
#    Ex: consensus
    'consensus',
# 15) minParentFold - The minimum abundance of potential "parents" of a sequence being
#               tested as chimeric, expressed as a fold-change versus the abundance of the sequence being
#               tested. Values should be greater than or equal to 1 (i.e. parents should be more
#               abundant than the sequence being tested).
#    Ex: 1.0
    1.0,
### SPEED ARGUMENTS ###
#
# 16) nthreads - The number of threads to use.
#                 Special values: 0 - detect available and use all.
#    Ex: 1
     $opt_threads,
# 17) nreads_learn - The minimum number of reads to learn the error model from.
#                 Special values: 0 - Use all input reads.
#    Ex: 1000000
     $dada2_min_reads,

# 18)  dir (general use),
	${temp_dir},

# 19) make plots (if "do_plots" will make pdfs)
    ${do_plots},

# 20) Taxonomy database || 'skip'
	${opt_tax_db},

# 21) Save rds. 'save' => yes
	${save_rds},

# 22) Just concat 0=false
    ${opt_justconcat},

# 23) Pool
    ${opt_process_pool}
);


## RUN DADA2
$S->verbose("Running DADA2...");
my $cmd_dada_paired = "Rscript --vanilla $RealBin/D2-dada.R " . join(" ", @dada_paired_args);

$S->verbose($cmd_dada_paired) if ($opt_debug);
$S->verbose("Dada2 script parameters:");
my $i_index = 0;
for my $arg (@dada_paired_args) {
  my $label = $dada_paired_parameter_labels[$i_index];
  $i_index++;
	if ($opt_verbose) {
	  say STDERR  " * [$i_index]", color('blue'), " $label: " , color('reset'), "$arg";
	}
	$S->writelog("Dada2 parameter\t$label\t$arg");
}

# RUN_DADA2
my $dada2_execution = $S->run($cmd_dada_paired, { candie => 1 });
addToTemplate($LOG_TEMPLATE, 'DADA2_LOG', runToHtml($dada2_execution));
notify_message("DADA2 finished. Finalizing...");

if ($opt_debug){
	$S->writelog("Counting filtered reads from DADA2");
	my $forcount = $S->run("seqfu stats " . File::Spec->catfile("${temp_dir}", 'for', 'filtered', '*.gz'));
	my $revcount = $S->run("seqfu stats " . File::Spec->catfile("${temp_dir}", 'rev', 'filtered', '*.gz'));
	$S->writelog("Forward: " . $forcount->{stdout});
	$S->writelog("Reverse: " . $revcount->{stdout});

}
# Check exit status of DADA2 wrapper
if ( $dada2_execution->{'exit'} != 0 ) {
    say STDERR color('red'), " ~~ DADA2 ERROR ~~ ";
	$S->verbose("DADA2 Failed:\nERR>\n$dada2_execution->{stderr}\nOUT>\n$dada2_execution->{stdout}");
	$S->verbose("Input files: " . "${temp_dir}/{for,rev}");
	our $EXIT_QUIET = 0;
	exit 1;
} else {
	$S->verbose("DADA2 output:\n" . color('yellow') . $dada2_execution->{stdout} . color('reset')) if ($opt_debug);
	if ($opt_debug) {
	   $S->verbose(" * Temporary directory kept:  $temp_dir/{for,rev}");
	} else {
	   $S->run(qq(rm -rf "$temp_dir"/for "$temp_dir"/rev));
	}
}


# Copy output files

copy(File::Spec->catfile($dada2_temp, 'stats.tsv'),
	 File::Spec->catfile($opt_output_directory, 'dada2_stats.tsv')) || die " ERROR:\n Unable to copy stats.tsv file from $dada2_temp to $opt_output_directory.\n";


if (! $opt_skip_plots ) {
	copy(File::Spec->catfile("$temp_dir", 'quality_R1.pdf'),
		 File::Spec->catfile(File::Spec->catdir("$opt_output_directory", 'qc'), 'quality_R1.pdf')) || 
		 	die " ERROR:\n Unable to copy quality_R1.pdf file from <$temp_dir> to <$opt_output_directory>.\n";
	copy(File::Spec->catfile("$temp_dir", 'quality_R2.pdf'),
		 File::Spec->catfile(File::Spec->catdir("$opt_output_directory", 'qc'), 'quality_R2.pdf')) || 
		 	die " ERROR:\n Unable to copy quality_R2.pdf file from <$temp_dir> to <$opt_output_directory>.\n";

}


## PROCESS DADA2 MAIN OUTPUT
my $dada2_tsv   = File::Spec->catfile($dada2_temp, 'dada2.tsv');
my $dada2_file  = File::Spec->catfile($opt_output_directory, 'feature-table.tsv');
my $repseq_file = File::Spec->catfile($opt_output_directory, 'rep-seqs.fasta');

# Where to store DADA2 tsv output
my $final_dada2_tsv = File::Spec->catfile($opt_output_directory, 'dada2_raw.tsv');
 
if ($opt_justconcat) {
	my $join_cmd = qq(dadaist2-mergeseqs -i "$dada2_tsv" -o "$final_dada2_tsv");
	$S->run($join_cmd, { candie => 0});
} else {
	copy($dada2_tsv,
		$final_dada2_tsv) || 
		die " ERROR:\n Unable to copy feature table file from $dada2_temp to $opt_output_directory.\n";
}

my @header = ();
open my $DADA_TSV,   '<', "$dada2_tsv"   || die " ERROR:\n Unable to open <$dada2_tsv>\n";
open my $dada2_FO,   '>', "$dada2_file"  || die " ERROR:\n Unable to write to <$dada2_file>\n";
open my $repseqs_FO, '>', "$repseq_file" || die " ERROR:\n Unable to write to <$repseq_file>\n";
my $feature_counter = 0;
$S->verbose("DADA2 Finished.");


# Dada2 taxonomy
if (defined $chooser->{'database'} and $chooser->{'database'} eq 'dada2') {
	my $taxonomy_file = File::Spec->catfile("$temp_dir", 'taxonomy.tsv');
	$S->verbose("Converting " . $chooser->{'database'}.  " taxonomy output: " . $taxonomy_file);
	die " Taxonomy not found: $taxonomy_file\n" unless (-e $taxonomy_file);
  @taxonomy = loadTaxonomyArray($taxonomy_file);
	die " Error loading taxonomy: $taxonomy_file\n"if (scalar @taxonomy == 0);
}

# Create OTU Table and FASTA file from TSV
while ( my $line = readline($DADA_TSV) ) {
	if ($line=~/^#/) {
		$line=~s/_R1.fastq.gz//g;
		print {$dada2_FO} $line;
		chomp($line);
		@header = split /\t/, $line;
	} else {
		chomp($line);
		$feature_counter++;
		my ($sequence, @values) = split /\t/, $line;
		my $name;
    my $comment = '';
		if ($opt_seq_prefix eq 'MD5') {
			$name = md5_hex($sequence);
		} else {
			$name = $opt_seq_prefix . $feature_counter;
      #$comment = "\t" . $taxonomy[$feature_counter - 1] if ($chooser->{'database'} eq 'dada2'); ##TODO
		}

		say {$repseqs_FO} '>' , $name, $comment, "\n", $sequence;
		say {$dada2_FO} $name, "\t", join("\t", @values);
	}

}

$S->verbose("$feature_counter representative sequences found.");

## Clean cross-talk
if ($opt_crosstalk) {
	$S->verbose("Removing crosstalk");
	my $temp_table = File::Spec->catfile($opt_output_directory, 'feature-table.bak');
	move($dada2_file, $temp_table) || die "Unable to copy intermediate file $dada2_file to $temp_table.\n";
	my $uncross = $S->run(qq($RealBin/dadaist2-crosstalk -i "$temp_table" -o "$dada2_file" 2>&1));
	if ($opt_debug) {
		say STDERR " * Raw table: ", $temp_table;
		say STDERR $uncross->{stdout};
	} else {
		unlink "$temp_table";
	}
}

## CHECK DADA2 DATA LOSS
my $dada_failed = 0;
my $data_loss = stats_check(File::Spec->catfile($dada2_temp, 'stats.tsv'), $opt_maxloss);
if ($data_loss->{pass} == 0) {
	say STDERR color('red'), "DADA2 ERROR:", color('reset');
	$S->verbose("DADA2 filtered too many reads: " . $data_loss->{percentage} . "% from total " . $data_loss->{"reads_input"} . " to " . $data_loss->{"reads_non-chimeric"});
	$dada_failed = 1;
} else {
	$S->verbose("DADA2 filtered " . $data_loss->{percentage} . "% from total " . $data_loss->{"reads_input"} . " to " . $data_loss->{"reads_non-chimeric"});
}
if ( defined $decipher_db ) {
	# Generate 'taxonomy.tsv' using the DECIPHER database
	
	$S->verbose("Assigning taxonomy using DECIPHER: " . $decipher_db);
	die  "$repseq_file?\n" if not -e "$repseq_file";
	my $decipher_cmd = qq($RealBin/dadaist2-assigntax --mode -i "$repseq_file" -o "$temp_dir" -r "$decipher_db" -t "$opt_threads");
	

	my $taxRun = $S->run($decipher_cmd, { candie => 1});
	addToTemplate($LOG_TEMPLATE, 'TAXONOMY_LOG', runToHtml($taxRun));
}

if ($opt_tax_db ne 'skip' and $chooser->{'database'} eq 'decipher') {
	my $f = File::Spec->catfile("$temp_dir", 'taxonomy.tsv');
	$S->verbose("Converting " . $chooser->{'database'}.  " taxonomy output: " . $f);
	die " Taxonomy not found: $f\n" unless (-e $f);
  	@taxonomy = loadTaxonomyArray($f);
	die " Error loading taxonomy: $f\n" if (scalar @taxonomy == 0);
}


my %output_files = (
 'feature-table' => "$dada2_file",
 'rep-seqs'      => "$repseq_file",
);

## Generate TAX OTUS
if ($opt_tax_db ne 'skip' and defined $chooser->{'database'}) {
	my $repseq_tax_file = File::Spec->catfile($opt_output_directory, 'rep-seqs-tax.fasta');
	open (my $I, '<', File::Spec->catfile($opt_output_directory, 'dada2_raw.tsv'));
	open (my $O, '>', "$repseq_tax_file") || die " ERROR:\n Unable to write to <$repseq_tax_file>\n";
	if (scalar @taxonomy == 0) {
		die " Error loading taxonomy [db:$opt_tax_db and db_type:" . $chooser->{'database'} . "]\n";
	}
	my $seq_counter = 0;
	while (my $line = readline( $I  )) {
		chomp($line);
		next if ($line =~/^#/);
		$seq_counter++;
		my @fields = split /\t/, $line;
		my $sequence = $fields[0];
		my $comment = '';
		if (defined $taxonomy[$seq_counter - 1]) {
			 $comment = "\t" . $taxonomy[$seq_counter - 1];
		} else {
			$S->verbose("No taxonomy found for sequence $seq_counter [$taxonomy[0]...]");
		}
		my $name;
		if ($opt_seq_prefix eq 'MD5') {
			$name = md5_hex($sequence);
		} else {
			$name = $opt_seq_prefix . $seq_counter;
		}
		say {$O} '>', $name, $comment, "\n", $sequence;
		
	}
} else {
	
}

# Copy optional files
if ($opt_tax_db ne 'skip' and not defined $decipher_db) {
	# Taxonomy generated via DADA2
  $output_files{'dada-taxonomy-table'} = File::Spec->catfile("$opt_output_directory", 'taxonomy.txt');
	$S->verbose("Copying DADA2 taxonomy: " . File::Spec->catfile("$temp_dir", 'taxonomy.tsv')) if $opt_debug;
	copy(File::Spec->catfile("$temp_dir", 'taxonomy.tsv'),
	     File::Spec->catfile("$opt_output_directory", 'taxonomy.txt')) || 
			 die " ERROR:\n Unable to copy taxonomy.txt file from " . File::Spec->catfile("$temp_dir", 'taxonomy.txt') . " to $opt_output_directory.\n";
} elsif (defined $decipher_db) {
	# Taxonomy was generated via DECIPHER
	$output_files{'decipher-taxonomy-table'} = File::Spec->catfile("$opt_output_directory", 'taxonomy.txt');
	copy(File::Spec->catfile("$temp_dir", 'taxonomy.tsv'),
	     File::Spec->catfile("$opt_output_directory", 'taxonomy.txt')) || 
			 die " ERROR:\n Unable to copy taxonomy.txt file from " . File::Spec->catfile("$temp_dir", 'taxonomy.txt') . " to $opt_output_directory.\n";
}

if ($save_rds eq 'save') {
  $output_files{'feature-table-RDS'} = File::Spec->catfile("$opt_output_directory", 'feature-table.rds');
  copy(File::Spec->catfile("$dada2_temp", 'dada2.rds'),
	     File::Spec->catfile("$opt_output_directory", 'feature-table.rds')) || die " ERROR:\n Unable to copy RDS file from " . File::Spec->catfile("$temp_dir",  'dada2.rds') . " to $opt_output_directory.\n";

}

# Make Tree
my $msa_output  = File::Spec->catfile("$opt_output_directory", 'rep-seqs.msa');
my $tree_output = File::Spec->catfile("$opt_output_directory", 'rep-seqs.tree');
my $msa_cmd     = qq(clustalo -i "$repseq_file" -o "$msa_output" --outfmt=fasta --force);
my $tree_cmd    = qq(fasttree -nt -gtr -no2nd -spr 4 -quiet "$msa_output" > "$tree_output");

if (not $opt_skip_tree) {
	# TREE
	$S->verbose("Multiple sequence alignment and tree generation");
	my $msa_exec    = $S->run($msa_cmd, { candie => 1 });
	my $tree_exec   = $S->run($tree_cmd , { candie => 1 });
	if ( $msa_exec->{'exit'} != 0 or $tree_exec->{'exit'} != 0) {
	say STDERR color('red'), " *** ERROR *** " if ($opt_verbose);
		$S->verbose("Clustalo execution failed.\n cmd> $msa_cmd\noutput> $msa_exec->{stderr}\nstderr> $msa_exec->{stdout}");
		exit 1;
	} elsif ( $tree_exec->{'exit'} != 0) {
	say STDERR color('red'), " *** ERROR *** " if ($opt_verbose);
	$S->verbose("FastTree execution failed.\n cmd> $tree_cmd\noutput> $tree_exec->{stderr}\nstderr> $tree_exec->{stdout}");
		exit 1;
	} else {
		$S->verbose("Feature tree generated");
	$output_files{'multiple-alignment'} = $msa_output;
	$output_files{'features-tree'} = $tree_output;
	}
}


if (not $dada_failed or $opt_force) {
		# Export MicrobiomeAnalyst + Rhea

		my $export_cmd = qq($RealBin/dadaist2-exporter -i "$opt_output_directory");
		my $exported = $S->run($export_cmd, { candie => 1});

		if ($exported->{'exit'} == 0) {
			$S->verbose("Exporting MicrobiomeAnalyst");
			$output_files{'mba-files'} = File::Spec->catdir($opt_output_directory, 'MicrobiomeAnalyst');
		} else {
			$S->verbose("MicrobiomeAnalyst not ready.");
		}

		# Generate phyloseq
		my $ps_out = File::Spec->catfile("$opt_R", 'phyloseq.rds');
		my $ps_cmd = qq($RealBin/dadaist2-phyloseqMake -i "$opt_output_directory" -o "$ps_out");
		my $ps_exec = $S->run($ps_cmd, {candie => 1});

		if ($ps_exec->{exit} != 0) {
			$S->verbose("PhyloSeq file not generated: " . $ps_out);
			#if ($opt_debug) {
				$S->verbose("Diagnostics: " . $ps_exec->{stdout} . "\n" . $ps_exec->{stderr});
			#}
		} else {
			$S->verbose("Generating PhyloSeq object");
			$output_files{'phyloseq'} = $ps_out;
		}


		# Experimental: Rhea
		my $Rhea_dir   = File::Spec->catdir($opt_output_directory, 'Rhea');
		my $Rhea_table = File::Spec->catfile($Rhea_dir, 'OTUs-Table.tab');
		my $Rhea_norm  = File::Spec->catfile($Rhea_dir, 'OTUs_Table-norm.tab');
		my $rhea_norm  = qq($RealBin/dadaist2-normalize -i "$Rhea_table" -o "$Rhea_dir");
		my $rhea_alpha = qq($RealBin/dadaist2-alpha -i "$Rhea_norm" -o "$Rhea_dir");
		my $RheaNormRun= $S->run($rhea_norm,  { candie => 1});
		my $RheaAlphaRun=$S->run($rhea_alpha, { candie => 1});
		if ($RheaNormRun->{exit} == 0 and $RheaAlphaRun->{exit} == 0) {
			$S->verbose("Rhea normalization/alpha finished.");
			$output_files{'rhea'} = $Rhea_dir;
		} else {
			$S->verbose("Rhea input files prepared, but normalization and alpha div. skipped.");
		}
}


# Print relevant output files produced
my $output_message = "Dadaist finished, output files saved:\n";
for my $f (sort keys %output_files) {
  $output_message .= " * $f: " . $output_files{$f} . "\n";
}
$S->verbose($output_message);

if (notify_message("Dadaist2 finished.")) {
	$S->run(qq(open "$opt_output_directory"), { candie => 1});
}

# Delete, or not, temp dir
our $EXIT_QUIET = 1;
if ($opt_debug) {
  $S->verbose("Temporary directory _not_ deleted: $temp_dir");
} else {
  $S->verbose("Cleaning up");
  $S->run(qq(rm -rf "$temp_dir"));

}


sub version {
	say "$PROGRAM v$VERSION";
	exit 0;
}

sub usage {
	say STDERR<<END;
USAGE:
 dadaist2 [options] [-d DB] -i INPUT_DIR -o OUTPUT_DIR

 dadaist2 --help for full manual
END

 if ($_[0]) {
	our $EXIT_QUIET = 1;
 	say "ERROR: $_[0]";
 	exit 1;
 }
}

sub check_metadata {
  my $file = shift @_;
  my $I;
  if (not open ($I, '<', $file)) {
    $S->verbose("Unable to load metadata from $file");
    exit 1;
  }
  my $count_lines = 0;
  my %samples = ();
  while (my $line = readline($I)) {
    $count_lines++;
    next if ($line =~/^#/);
    my @fields = split /$opt_metadata_separator/, $line;
    $samples{ $fields[0] }++;
  }
  return \%samples;
}

sub loadTaxonomyArray {
  # load taxonomy from files
  my $file = $_[0];
  my @results = ();
  return () if (! -e "$file");
  open(my $I, '<', "$file") || return ();
  my $c = 0;
  while (my $line = readline($I)) {
    chomp($line);
    $c++;
    next if ($c == 1);
    if ($line =~/^\d+/) {
			my @local = ();
      my @fields = split / /, $line;
      shift(@fields);
      
			my $newstr = join(" ", @fields);
			while ($newstr =~/(".*?"|\S+)/g) {
				push(@local, $1)
			}
			push(@results, join(';', @local));
    }
  }
	
  return @results;
}

sub get_qualified_from_seqfu { 

	my ($output_raw) = @_;
	my @lines = split /\n/, $output_raw;
	my %sample_qual = (); 
	my $tot = 0;
	my $sum = 0;
	my @discarded_samples = ();
	for my $line (@lines) {
		my ($sample, undef, undef, undef, undef, $qual) = split /\t/, $line;
		$sample_qual{$sample}  = $qual;
		$sum += $qual;
		$tot += 1;
	}
  	if ($tot == 0) {
		die "ERROR: Qualified quality from SeqFu: 0 samples received.\n";
	}

	my $avg = sprintf("%.2f", $sum/$tot);
	my $ths = $avg - ($avg * 0.25);
	
	$tot = 0;
	$sum = 0;
	for my $sample_ID (sort { $sample_qual{$a} <=> $sample_qual{$b}} keys %sample_qual) {
		if ( $sample_qual{$sample_ID} < $ths ) {
			push(@discarded_samples, $sample_ID)
		} else {
			$tot+=1;
			$sum+=$sample_qual{$sample_ID};

		}
	}

	my $qualified_avg = sprintf("%.2f", $sum/$tot);
	return (int($qualified_avg), \@discarded_samples);
	

}

sub get_qualified_positions {
	# Scan a list of qualityes and return the boundaries of qualified quality
	my $len = @_;
	my $start = undef;
	my $end = undef;

	for (my $pos = 0; $pos < ($len - $qualified_wnd); $pos++) {
		my $avg = avg_array( @_[$pos..$pos+$qualified_wnd] );

		if ($avg > $qualified_min_qual) {
			$start = $pos unless defined $start;
		} else {
			next if not defined $start;
			$end = $pos unless defined $end;
		}


	}
	$end = $len unless defined $end;

	return($start, $end, $len);
}
sub avg_array {
	my $agg = 0;
	$agg += $_ for @_;
	return $agg/@_;
}
sub sum_array {
	my $agg = 0;
	$agg += $_ for @_;
	return $agg
}

sub load_json_from_file {
	my $file = shift @_;
	my $json_read = $S->run(qq(cat "$file"));


	my $data;
	eval {
		$data = decode_json $json_read->{stdout};
	};
	if ($@) {
		die " FATAL ERROR: Unable to decode JSON from <$file>:\n$@\n";
	}
	return $data;

}
sub get_file_reads {
	my ($dir) = @_;
	my $is_paired = undef;
	$dir = File::Spec->rel2abs($dir);
	my @files = <"$dir"/*.*>;
	my %samples;
	my $counter_for = 0;
	my $counter_rev = 0;
	for my $file (sort @files) {
		next if (substr($file, 0, 1) eq '.');
		next if ($file !~/\.(fastq|fq)/i);
		my ($id) = split /$opt_id_separator/, basename($file);

		if ($id=~/^\d/) {
			$S->verbose("[Warning] Sample id ($id) starting with a digit can introduce discrepancies in R analyses.");
			$S->verbose("[Warning] This will introduce an error in the next release.");
		
		}

		if ( $samples{$id}{'for'} and $samples{$id}{'rev'} ) {
			$S->verbose("Error: Sample ID <$id> was already found. The delimiter ($opt_id_separator) is probably contained in the filename ID.");
			die "FATAL ERROR: Duplicateed SampleID <$id> in file $file.\n";
		}
		if ($file =~/$opt_for_tag/) {
			$counter_for++;
			$samples{$id}{'for'} = $file;
		} elsif ($file =~/$opt_rev_tag/) {
			$counter_rev++;
			$samples{$id}{'rev'} = $file;
		} else {
			$S->verbose("Skipping file <$file>: missing $opt_for_tag/$opt_rev_tag");
		}
	}
	if ($counter_for == 0 and $counter_rev == 0) {
		die "FATAL ERROR: No samples found in <$dir>.\n";
	} elsif ($counter_for == $counter_rev) {
		$is_paired = 1;
	} elsif ($counter_rev == 0 ) {
		$is_paired = 0;
	} else {
		die "FATAL ERROR: $counter_rev paired samples found, but $counter_for forward pairs found.\n";
	}
	return (\%samples, $is_paired, $counter_rev);
}

sub iupac_rc  {
	my $seq = $_[0];
	$seq =~ tr/acgtrymkbdhvACGTRYMKBDHV/tgcayrkmvhdbTGCAYRKMVHDB/;
	$seq = reverse($seq);
	return $seq; 
}

sub valid_gzipped_fasta {
	my $file = shift @_;
	my $reader;
	if ($file =~/\.f.*\.gz$/) {
		return 1;
	} else {
		return 0;
	}

	
	# my $check = eval {
	# 	$reader = FASTX::Reader->new({ filename => "$file"});
	# };
	# if ($@) {
	# 	$S->verbose("Reference $file reading error: unable to open:\n$@");
	# 	return 0;
	# } else {
	# 	if ($reader->{compressed} == 1) {
	# 		return 1;
	# 	} else {
	# 		$S->verbose("Reference $file should be compressed with GZip!");
	# 	}

	# }

	exit;
}


sub stats_check {
  my ($file, $maxloss) = @_;
  my $data;
  my @labels = ('input'  ,'filtered'  ,'denoised'  ,'merged'  ,'non-chimeric');
  my $I;
  if (not open ( $I, '<', $file) ) {
    $data->{pass} = 0;
    $data->{label} = 'fileNotFound';
  }
  my $c = 0;
  my %counts = ();
  my %loss   = (); 
  
  $data->{pass} = 1;
  $data->{step} = '';

  # Parse DADA2 stats
  while (my $line = readline($I)) {
    chomp($line);
    $c++;
    my @fields = split /\t/, $line;
    if ($c == 1) {  
      #x,input,filtered,denoised,merged,non-chimeric
      return 0 if ($fields[1] ne 'input' and $fields[4] ne 'merged');
    } else {
      my $sample = shift @fields;

      for my $label (@labels) {
        $counts{$label} += shift @fields;
      }
    }
  }

  # Check at each step the data loss
  for my $label (@labels) {
    $loss{$label} += $counts{$label} / $counts{'input'} if ($counts{'input'});
    if ($loss{$label} < $maxloss) {
      $data->{pass} = 0;
      $data->{step} = $label;
      $data->{percentage} = sprintf("%.4f", 100 * $counts{$label} / $counts{'input'});
    } 
    $data->{"reads_$label"} = $counts{$label};
  }

  # Set data loss if pass
  if ($data->{pass} == 1) {
    $data->{step} = 'non-chimeric';
    $data->{percentage} = sprintf("%.4f", 100 * $counts{$data->{step}} / $counts{'input'}) if ($counts{'input'});
  }

  
  return $data;
}

sub check_required_binaries {

	for my $bin (sort keys  %{ $_[0]} ) {
		my $run = $S->run("${ $_[0]}{$bin}", { candie => 1});
		addToTemplate($LOG_TEMPLATE, 'INIT_LOG', runToHtml($run), 'append');
		if ($run->{exit} == 0 and $run->{stdout} !~/not found/) {
			my $pass = $run->{stdout} =~/\S/ ? $run->{stdout} : '<pass>';
			if ($opt_debug) {
				$S->verbose("Checking $bin: $pass") 
			} else {
				say STDERR " * $bin: ", color('yellow'), $pass, color('reset');
			}

		} else {
			if ($run->{exit} != 0 ) {
				die "FATAL ERROR:\nDependency '$bin' returned non zero exit value.\nShell: ${ $_[0]}{$bin}\n";
			} else {
				die "FATAL ERROR:\nDependency '$bin' not found or not.\nShell: ${ $_[0]}{$bin}\n";

			}
		}

	}
}

sub prepare_temporary_directories {
 my $temp_dir = shift @_;
 $temp_dir = File::Spec->rel2abs($temp_dir);
	mkdir "$temp_dir/for"            || die " FATAL ERROR:\n Unable to create directory $temp_dir/for\n $!\n";
	mkdir "$temp_dir/rev"            || die " FATAL ERROR:\n Unable to create directory $temp_dir/rev\n $!\n";
	mkdir "$temp_dir/for/filtered"   || die " FATAL ERROR:\n Unable to create directory $temp_dir/for/filtered\n $!\n";
	mkdir "$temp_dir/rev/filtered"   || die " FATAL ERROR:\n Unable to create directory $temp_dir/rev/filtered\n $!\n";
	mkdir "$dada2_temp"              || die " FATAL ERROR:\n Unable to create directory $temp_dir/dada2\n $!\n";
  $S->verbose("Temporary directory: $temp_dir");
  return $temp_dir;
}

sub notify_message {
	my ($message) = @_;
	if ($opt_notify) {
		my $timestamp = getLoggingTime();
		my $cmd = "osascript -e \"display notification \\\"$message\\\" with title \\\"Dadaist2\\\" subtitle \\\"\\\" \" ";
		$cmd .= "|| zenity --notification --text \"DADAIST2: $message\" ";
		$cmd .= "|| notify-send \"DADAIST2: $message\"";
		my $i = $S->run($cmd,
		{candie => 1});

		if ($i->{'exit'} == 0) {
			return 1;
		} else {
			return 0;
		}
	} 
}

sub getLoggingTime {

    my ($sec,$min,$hour,$mday,$mon,$year,$wday,$yday,$isdst)=localtime(time);
    my $nice_timestamp = sprintf ( "%04d%02d%02d %02d:%02d:%02d",
                                   $year+1900,$mon+1,$mday,$hour,$min,$sec);
    return $nice_timestamp;
}
sub getSections {
	my $template = $_[0];
	my $log_structure = {};
	while ($template =~/\{\{([\w_]+)\}\}/g) {
		my $key = $1;
			$log_structure->{$key} = $unfilled;
	}
	return $log_structure;
}
sub num_cpu {
  if ( $^O =~ m/linux/i ) {
    my($num) = qx(grep -c ^processor /proc/cpuinfo);
    return $1 if $num =~ m/^(\d+)/;
  }
  elsif ( $^O =~ m/darwin/i ) {
    my($num) = qx(system_profiler SPHardwareDataType | grep Cores);
    return $1 if $num =~ /.*Cores: (\d+)/;
  }
  return 1;
}

sub addToTemplate {
	my ($template, $var, $value, $append) = @_;
	if (defined $template->{"$var"}) {
		if ($template->{"$var"} ne $unfilled and not defined $append) {
			# Trying to overwrite variable without appending
			say STDERR "WARNING: Replacing '$var' section in LOG.\n";
		}
		if (defined $append) {
			if ($template->{"$var"} eq $unfilled) {
					# Initialize
					$template->{"$var"} = '<br>' . $value;
			} else {
					# Append
					$template->{"$var"} .= '<br>' . $value;
			}
			
		} else {
			$template->{"$var"} = $value;
		}
		
		return $template;
	} else {
		say STDERR Dumper $template;
		die "Attempting to add to template '$var': not found.\n";
	}
}
sub getRamGb {
	if ($^O eq 'darwin') {
		my $mac = `system_profiler SPHardwareDataType 2>/dev/null | grep  "Memory" `;
		if (length($mac) > 0 and $? == 0) {
			$mac =~/(\d+)/;
			return $1;
		} else {
			return 'Unknown';
		}
	} elsif ($^O eq 'linux')  {
		my $linux = `free -m  2>/dev/null | grep ^Mem`;
		if (length($linux) > 0 and $? == 0) {
			$linux =~/(\d+)/;
			return $1 / 1024;
		} else {
			return 'Unknown';
		}
	} else {
	return 'Unsupported OS';
	}
}

sub getSysinfo {
	if ($^O eq 'darwin') {
		my @mac = `system_profiler SPHardwareDataType 2>/dev/null | grep  ": " `;
		my $data = '<table>';
		if ($#mac > 0 and $? == 0) {
			for my $line (@mac) {
				chomp($line);
				my ($key, $val) =  $line =~/\s*([^:]+)\s*:\s*([^:]+)/;
				$data .= qq(<tr><td>$key</td><td>$val</td></tr>\n);
			}
			return "$data</table>\n";
		}
	} elsif ($^O eq 'linux')  {
		my $linux = `uname -a 2>/dev/null`;
		if (length($linux) > 0 and $? == 0) {
			return $linux;
		}
	} else {
		return 'Unsupported OS';
	}
	return '';
}

sub getCpu {
	if ($^O eq 'darwin') {
		my $mac = `system_profiler SPHardwareDataType 2>/dev/null | grep  "Total Number of Cores" `;
		if (length($mac) > 0 and $? == 0) {
			$mac =~/(\d+)/;
			return $1;
		}
	} elsif ($^O eq 'linux')  {
		my $linux = `cat /proc/cpuinfo 2>/dev/null | grep -c processor`;
		if (length($linux) > 0 and $? == 0) {
			$linux =~/(\d+)/;
			return $1;
		}
	} else {
	return 'Unsupported OS';
	}
}

END {
    # cleanup
		my $END_TIME = time();
		my $run_time = formatsec($END_TIME - our $START_TIME);

		addToTemplate($LOG_TEMPLATE, 'EXECUTION_TIME', $run_time);
		addToTemplate($LOG_TEMPLATE, 'EXECUTION_TIMESTAMP', our $EXECUTION_TIMESTAMP);
		my $html_log = makeTemplate($TEMPLATE, $LOG_TEMPLATE);

		if ($STARTED and defined $opt_output_directory and defined $html_log) {
			open (my $O, '>', "$opt_output_directory/dadaist2.html");
			print {$O} $html_log;
			close $O;
		}
		if (our $EXIT_QUIET == 0 and our $STARTED > 0) {
			say STDERR "Dadaist2 execution finished (", formatsec($run_time),")";
		}
}

sub directoryIsEmpty {
	my $dir = $_[0];
	my $empty = 1;
	if (opendir(my $D, "$dir")) {
		while (my $file = readdir($D)) {
			if ($file ne '.' and $file ne '..') {
				$empty = 0;
				last;
			}
		}
		closedir($D);
	}
	return $empty;
}

sub formatsec { 
  my $time = shift; 
	if ($time =~/\D/) {
		return $time;
	}
  my $days = int($time / 86400); 
   $time -= ($days * 86400); 
  my $hours = int($time / 3600); 
   $time -= ($hours * 3600); 
  my $minutes = int($time / 60); 
  my $seconds = $minutes < 1 ? sprintf("%.2f", $time - ($minutes*60) ) : $time % 60;
 
  $days = $days < 1 ? '' : $days .'d '; 
  $hours = $hours < 1 ? '' : $hours .'h '; 
  $minutes = $minutes < 1 ? '' : $minutes . 'm '; 
 
  $time = $days . $hours . $minutes . $seconds . 's'; 
  return $time; 
}

sub runToHtml {
	my $run = $_[0];
	my $answer = '<p>';
	my $cmd = $run->{cmd} // '<em>not_defined</em>';
	my @params = split /\s+/, $cmd;
	my $bin = $params[0];
	while ($cmd =~/(\s)(\-\w)(\s)/g) {
		$cmd =~s/(\s)(\-\w)(\s)/ <br><b>$2<\/b> /;
	}
	while ($cmd =~/(\s)(\--[A-Za-z0-9_-]+)(\s)/g) {
		$cmd =~s/(\s)(\--[A-Za-z0-9_-]+)(\s)/ <br><b>$2<\/b> /;
	}
	my $exitstatus = $run->{exit} // '<em>not_defined</em>';
	my $span = formatsec($run->{duration}) // '<em>unknown</em>';
	my $time = $run->{time} // '<em>unknown</em>';
	$answer .= '<h3>Shell command <span style="color: #b44">' . $bin . '</span></h3>';
	$answer .= qq(<p>Executed at <em>$time</em> for <em>$span</em> and  returned $exitstatus </p>');
	$answer .= "<p><strong>Command</strong><pre>"  . $cmd . "</pre></p>\n";;
	$answer .= "<p><strong>Standard output</strong><pre>" . $run->{stdout} // '<em>None.</em>' . "</pre></p>\n";
	$answer .= "<p><strong>Standard error</strong><pre>" . $run->{stderr}  // '<em>None.</em>' . "</pre></p>\n";
	$answer .= "</p>\n";
}
sub paramToHtml {
	my ($param, $value) = @_;
	return qq(<tr>
    	<td>$param</td>
    	<td>$value</td>
  	</tr>);
 
}
sub makeTemplate {
	my ($template, $keys)  = @_;
	if (length($template) < 100) {
		return '';
	}
	for my $placeholder (keys %{ $keys }) {
		my $value = $keys->{$placeholder};
		
		$template =~s/\{\{$placeholder\}\}/$value/g;
	}
	return $template;
}
sub loadTemplate {
	return q(<!DOCTYPE html>
<html>
<head>
<meta name="viewport" content="width=device-width, initial-scale=1">
 
<style>
h1, h2, h3, h4, h5, h6, table, td, p { font-family: Helvetica, Verdana, Arial; }
pre { white-space: pre-wrap; 
overflow-x: auto;
background-color: lightyellow;
font-family: Consolas, "Andale Mono WT", "Andale Mono", "Lucida Console", "Lucida Sans Typewriter", "DejaVu Sans Mono", "Bitstream Vera Sans Mono", "Liberation Mono", "Nimbus Mono L", Monaco, "Courier New", Courier, monospace; }
h1{ color: #b44;}
.accordion {
  background-color: #eee;
  color: #444;
  cursor: pointer;
  padding: 18px;
  width: 100%;
  border: none;
  text-align: left;
  outline: none;
  font-size: 15px;
  transition: 0.4s;
}

.active, .accordion:hover {
  background-color: #ccc;
}

.accordion:after {
  content: '\002B';
  color: #777;
  font-weight: bold;
  float: right;
  margin-left: 5px;
}

.active:after {
  content: "\2212";
}

.panel {
  padding: 0 18px;
  background-color: white;
  max-height: 0;
  overflow: hidden;
  transition: max-height 0.2s ease-out;
}

/* Table */

table {
  font-family: arial, sans-serif;
  border-collapse: collapse;
  width: 100%;
}

td, th {
  border: 1px solid #eee;
  text-align: left;
  padding: 8px;
}
th {
  background-color: #ddd;
}
tr:nth-child(odd) {
  background-color: #eee;
}

.bar .bar-item{padding:8px 16px;float:left;width:auto;border:none;display:block;outline:0}
.bar .button{white-space:normal}
.bar .dropdown-hover,.bar .dropdown-click{position:static;float:left}
.bar-block .bar-item{width:100%;display:block;padding:8px 16px;text-align:left;border:none;white-space:normal;float:none;outline:0}
.bar-block .dropdown-hover .button,.bar-block .dropdown-click .button{width:100%;text-align:left;padding:8px 16px}
.bar-block .dropdown-hover .dropdown-content,.bar-block .dropdown-click .dropdown-content{min-width:100%}
.bar-block .dropdown-hover,.bar-block .dropdown-click{width:100%}
.bar-block.center .bar-item{text-align:center}.block{display:block;width:100%}
.bar-block.center .bar-item{text-align:center}.block{display:block;width:100%}
.bar{width:100%;overflow:hidden}.center .bar{display:inline-block;width:auto}
.bar{width:100%;overflow:hidden}.center .bar{display:inline-block;width:auto}
.black,.hover-black:hover{color:#fff!important;background-color:#000!important}
.btn,.button{-webkit-touch-callout:none;-webkit-user-select:none;-khtml-user-select:none;-moz-user-select:none;-ms-user-select:none;user-select:none}   
.btn,.button{border:none;display:inline-block;padding:8px 16px;vertical-align:middle;overflow:hidden;text-decoration:none;color:inherit;background-color:inherit;text-align:center;cursor:pointer;white-space:nowrap}
.button:hover{color:#000!important;background-color:#ccc!important}
.cell-row:before,.cell-row:after,.clear:after,.clear:before,.bar:before,.bar:after{content:"";display:table;clear:both}
.cell-row:before,.cell-row:after,.clear:after,.clear:before,.bar:before,.bar:after{content:"";display:table;clear:both}
.cell-row:before,.cell-row:after,.clear:after,.clear:before,.bar:before,.bar:after{content:"";display:table;clear:both} 
.container:after,.container:before,.panel:after,.panel:before,.row:after,.row:before,.row-padding:after,.row-padding:before,
.disabled,.btn:disabled,.button:disabled{cursor:not-allowed;opacity:0.3}.disabled *,:disabled *{pointer-events:none}
.display-topleft{position:absolute;left:0;top:0}.display-topright{position:absolute;right:0;top:0}
.dropdown-hover.mobile,.dropdown-hover.mobile .btn,.dropdown-hover.mobile .button,.dropdown-click.mobile,.dropdown-click.mobile .btn,.dropdown-click.mobile .button{width:100%}}
.dropdown-hover:hover > .button:first-child,.dropdown-click:hover > .button:first-child{background-color:#ccc;color:#000}
.hide-small{display:none!important}.mobile{display:block;width:100%!important}.bar-item.mobile,.dropdown-hover.mobile,.dropdown-click.mobile{text-align:center}
.tooltip,.display-container{position:relative}.tooltip .text{display:none}.tooltip:hover .text{display:inline-block}

.responsive {
  width: 100%;
  height: auto;
}
</style>
</head>
<body>
<h1>Dadaist2 {{VERSION}} execution log</h1>
<p>
This is not a report with the results of the pipeline, but the log of the events occurred to generate it.
</p>
<button class="accordion">Summary</button>
<div class="panel">

<!--summary -->
 

<div class="bar black">
  <button class="bar-item button" onclick="openCity('Summary')">Summary</button>
  <button class="bar-item button" onclick="openCity('Parameters')">Parameters</button>
  <button class="bar-item button" onclick="openCity('Host')">Machine</button>
</div>

<div id="Summary" class="container display-container city">
  <span onclick="this.parentElement.style.display='none'"
  class="button large display-topright">&times;</span>

  <h2>Summary</h2>
  <p>This reports contain technical details on the execution of 
  <a href="https://quadram-institute-bioscience.github.io/dadaist2">Dadaist2</a>.</p>

<table>
  <tr>
    <th>Parameter</th>
    <th>Value</th>
  </tr>
  <tr>
    <td>Input</td>
    <td>{{INPUT_PATH}}</td>
  </tr>
  <tr>
    <td>Output</td>
    <td>{{OUTPUT_PATH}}</td>
  </tr>
  <tr>
    <td>Timestamp</td>
    <td>{{EXECUTION_TIMESTAMP}}</td>
  </tr>
  <tr>
    <td>Execution time</td>
    <td>{{EXECUTION_TIME}}</td>
  </tr>

</table>


</div>

<div id="Parameters" class="container display-container city" style="display:none">
  <span onclick="this.parentElement.style.display='none'"
  class="button large display-topright">&times;</span>
  <h2>Parameters</h2>  
	<table>
  <tr>
    <th>Parameter</th>
    <th>Value</th>
  </tr>
	{{PARAMETERS}}
	</table>
</div>

<div id="Host" class="container display-container city" style="display:none">
  <span onclick="this.parentElement.style.display='none'"
  class="button large display-topright">&times;</span>
  <h2>Machine details</h2>
<table>
  <tr>
    <th>Timestamp</th>
    <th>Event</th>
  </tr>
  <tr>
    <td>Processors</td>
    <td>{{CPU_CORES}}</td>
  </tr>
  <tr>
    <td>RAM</td>
    <td>{{RAM_GB}} Gb</td>
  </tr>
  <tr>
    <td>OS</td>
    <td>{{UNAME}}</td>
  </tr>
</table>
</div>

 <br>
</div>
<!--/summary-->





<button class="accordion">Read trimming (process)</button>
<div class="panel">
  {{TRIMMING_LOG}}
</div>

<button class="accordion">DADA2 Denoising (process)</button>
<div class="panel">
  {{DADA2_LOG}}
</div>

<button class="accordion">Taxonomy (process)</button>
<div class="panel">
  {{TAXONOMY_LOG}}
</div>


<button class="accordion">Initial checks</button>
<div class="panel">
	{{INIT_LOG}}
</div>


<script>
var acc = document.getElementsByClassName("accordion");
var i;

for (i = 0; i < acc.length; i++) {
  acc[i].addEventListener("click", function() {
    this.classList.toggle("active");
    var panel = this.nextElementSibling;
    if (panel.style.maxHeight) {
      panel.style.maxHeight = null;
    } else {
      panel.style.maxHeight = panel.scrollHeight + "px";
    } 
  });
}
</script>


<script>
function openCity(cityName) {
  var i;
  var x = document.getElementsByClassName("city");
  for (i = 0; i < x.length; i++) {
    x[i].style.display = "none";  
  }
  document.getElementById(cityName).style.display = "block";  
}
</script>
</body>



</html>

)
}

__END__

=pod

=head1 SOURCE CODE AND DOCUMENTATION

The program is freely available at L<https://quadram-institute-bioscience.github.io/dadaist2>
released under the MIT licence. The website contains the full I<DOCUMENTATION> and we recommend 
checking for updates.

The paper describing Dadaist2 was published in:

Ansorge, R.; Birolo, G.; James, S.A.; Telatin, A. I<Dadaist2: A Toolkit to Automate and Simplify Statistical Analysis and Plotting of Metabarcoding Experiments>.
Int. J. Mol. Sci. 2021, 22, 5309. L<https://doi.org/10.3390/ijms22105309>

=cut
