changeset 0:e5953c54cfb5 draft

planemo upload for repository https://github.com/gpovysil/VariantAnalyzerGalaxy/tree/master/tools/variant_analyzer commit ee4a8e6cf290e6c8a4d55f9cd2839d60ab3b11c8
author mheinzl
date Sun, 04 Oct 2020 17:19:39 +0000
parents
children 2a505d46f682
files mut2read.py mut2read.xml mut2sscs.py mut2sscs.xml read2mut.py read2mut.xml test-data/Aligned_Families_test_data_VA.tabular test-data/DCS_Mutations_test_data_VA.tabular test-data/DCS_test_data_VA.bam test-data/DCS_test_data_VA.bam.bai test-data/Interesting_Reads_test_data_VA.fastq test-data/Interesting_Reads_test_data_VA.trim.bam test-data/Interesting_Reads_test_data_VA.trim.bam.bai test-data/SSCS_counts_test_data_VA.json test-data/SSCS_test_data_VA.bam test-data/SSCS_test_data_VA.bam.bai test-data/mutant_reads_summary_short_trim_test_data_VA.xlsx test-data/tag_count_dict_test_data_VA.json va_macros.xml
diffstat 19 files changed, 1641 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mut2read.py	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,155 @@
+#!/usr/bin/env python
+
+"""mut2read.py
+
+Author -- Gundula Povysil
+Contact -- povysil@bioinf.jku.at
+
+Takes a tabular file with mutations and a BAM file as input and prints
+all tags of reads that carry the mutation to a user specified output file.
+Creates fastq file of reads of tags with mutation.
+
+=======  ==========  =================  ================================
+Version  Date        Author             Description
+0.2.1    2019-10-27  Gundula Povysil    -
+=======  ==========  =================  ================================
+
+USAGE: python mut2read.py DCS_Mutations.tabular DCS.bam Aligned_Families.tabular Interesting_Reads.fastq
+                          tag_count_dict.json
+"""
+
+import argparse
+import json
+import os
+import sys
+
+import numpy as np
+import pysam
+
+
+def make_argparser():
+    parser = argparse.ArgumentParser(description='Takes a tabular file with mutations and a BAM file as input and prints all tags of reads that carry the mutation to a user specified output file and creates a fastq file of reads of tags with mutation.')
+    parser.add_argument('--mutFile',
+                        help='TABULAR file with DCS mutations.')
+    parser.add_argument('--bamFile',
+                        help='BAM file with aligned DCS reads.')
+    parser.add_argument('--familiesFile',
+                        help='TABULAR file with aligned families.')
+    parser.add_argument('--outputFastq',
+                        help='Output FASTQ file of reads with mutations.')
+    parser.add_argument('--outputJson',
+                        help='Output JSON file to store collected data.')
+    return parser
+
+
+def mut2read(argv):
+    parser = make_argparser()
+    args = parser.parse_args(argv[1:])
+
+    file1 = args.mutFile
+    file2 = args.bamFile
+    file3 = args.familiesFile
+    outfile = args.outputFastq
+    json_file = args.outputJson
+
+    if os.path.isfile(file1) is False:
+        sys.exit("Error: Could not find '{}'".format(file1))
+
+    if os.path.isfile(file2) is False:
+        sys.exit("Error: Could not find '{}'".format(file2))
+
+    if os.path.isfile(file3) is False:
+        sys.exit("Error: Could not find '{}'".format(file3))
+
+    # read mut file
+    with open(file1, 'r') as mut:
+        mut_array = np.genfromtxt(mut, skip_header=1, delimiter='\t', comments='#', dtype=str)
+
+    # read dcs bam file
+    # pysam.index(file2)
+    bam = pysam.AlignmentFile(file2, "rb")
+
+    # get tags
+    tag_dict = {}
+    cvrg_dict = {}
+
+    if mut_array.shape == (1,13):
+        mut_array = mut_array.reshape((1, len(mut_array)))
+
+    for m in range(len(mut_array[:, 0])):
+        print(str(m + 1) + " of " + str(len(mut_array[:, 0])))
+        chrom = mut_array[m, 1]
+        stop_pos = mut_array[m, 2].astype(int)
+        chrom_stop_pos = str(chrom) + "#" + str(stop_pos)
+        ref = mut_array[m, 9]
+        alt = mut_array[m, 10]
+
+        dcs_len = []
+
+        for pileupcolumn in bam.pileup(chrom.tostring(), stop_pos - 2, stop_pos, max_depth=100000000):
+
+            if pileupcolumn.reference_pos == stop_pos - 1:
+                count_alt = 0
+                count_ref = 0
+                count_indel = 0
+                count_n = 0
+                count_other = 0
+                count_lowq = 0
+                print("unfiltered reads=", pileupcolumn.n, "filtered reads=", len(pileupcolumn.pileups),
+                      "difference= ", len(pileupcolumn.pileups) - pileupcolumn.n)
+                for pileupread in pileupcolumn.pileups:
+                    if not pileupread.is_del and not pileupread.is_refskip:
+                        # query position is None if is_del or is_refskip is set.
+                        nuc = pileupread.alignment.query_sequence[pileupread.query_position]
+                        dcs_len.append(len(pileupread.alignment.query_sequence))
+                        if nuc == alt:
+                            count_alt += 1
+                            tag = pileupread.alignment.query_name
+                            if tag in tag_dict:
+                                tag_dict[tag][chrom_stop_pos] = alt
+                            else:
+                                tag_dict[tag] = {}
+                                tag_dict[tag][chrom_stop_pos] = alt
+                        elif nuc == ref:
+                            count_ref += 1
+                        elif nuc == "N":
+                            count_n += 1
+                        elif nuc == "lowQ":
+                            count_lowq += 1
+                        else:
+                            count_other += 1
+                    else:
+                        count_indel += 1
+                dcs_median = np.median(np.array(dcs_len))
+                cvrg_dict[chrom_stop_pos] = (count_ref, count_alt, dcs_median)
+
+                print("coverage at pos %s = %s, ref = %s, alt = %s, other bases = %s, N = %s, indel = %s, low quality = %s, median length of DCS = %s\n" %
+                      (pileupcolumn.pos, count_ref + count_alt, count_ref, count_alt, count_other, count_n,
+                       count_indel, count_lowq, dcs_median))
+    bam.close()
+
+    with open(json_file, "w") as f:
+        json.dump((tag_dict, cvrg_dict), f)
+
+    # create fastq from aligned reads
+    with open(outfile, 'w') as out:
+        with open(file3, 'r') as families:
+            for line in families:
+                line = line.rstrip('\n')
+                splits = line.split('\t')
+                tag = splits[0]
+
+                if tag in tag_dict:
+                    str1 = splits[4]
+                    curr_seq = str1.replace("-", "")
+                    str2 = splits[5]
+                    curr_qual = str2.replace(" ", "")
+
+                    out.write("@" + splits[0] + "." + splits[1] + "." + splits[2] + "\n")
+                    out.write(curr_seq + "\n")
+                    out.write("+" + "\n")
+                    out.write(curr_qual + "\n")
+
+
+if __name__ == '__main__':
+    sys.exit(mut2read(sys.argv))
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mut2read.xml	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,65 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<tool id="mut2read" name="DCS mutations to tags/reads:" version="1.0.1" profile="19.01">
+    <description>Extracts all tags that carry a mutation in the duplex consensus sequence (DCS)</description>
+    <macros>
+        <import>va_macros.xml</import>
+    </macros>
+    <expand macro="requirements"/>
+    <command><![CDATA[
+        ln -s '$file2' bam_input.bam &&
+        ln -s '${file2.metadata.bam_index}' bam_input.bam.bai &&
+        python '$__tool_directory__/mut2read.py' 
+        --mutFile '$file1'
+        --bamFile bam_input.bam
+        --familiesFile '$file3'
+        --outputFastq '$output_fastq' 
+        --outputJson '$output_json'
+    ]]>
+    </command>
+    <inputs>
+        <param name="file1" type="data" format="tabular" label="DCS Mutation File" optional="false" help="TABULAR file with DCS mutations. See Help section below for a detailed explanation."/>
+        <param name="file2" type="data" format="bam" label="DCS BAM File" optional="false" help="BAM file with aligned DCS reads."/>
+        <param name="file3" type="data" format="tabular" label="Aligned Families File" optional="false" help="TABULAR file with aligned families."/>
+    </inputs>
+    <outputs>
+        <data name="output_fastq" format="fastq" label="${tool.name} on ${on_string}: FASTQ"/>
+        <data name="output_json" format="json" label="${tool.name} on ${on_string}: JSON"/>
+    </outputs>
+    <tests>
+        <test>
+            <param name="file1" value="DCS_Mutations_test_data_VA.tabular"/>
+            <param name="file2" value="DCS_test_data_VA.bam"/>
+            <param name="file3" value="Aligned_Families_test_data_VA.tabular"/>
+            <output name="output_fastq" file="Interesting_Reads_test_data_VA.fastq" lines_diff="136"/>
+            <output name="output_json" file="tag_count_dict_test_data_VA.json" lines_diff="2"/>
+        </test>
+    </tests>
+    <help> <![CDATA[
+**What it does**
+
+Takes a tabular file with mutations, a BAM file of aligned DCS reads, and a 
+tabular file with aligned families as input and prints all tags of reads that 
+carry a mutation to a user specified output file and creates a fastq file of 
+reads of tags with a mutation.
+
+**Input** 
+
+**Dataset 1:** Tabular file with duplex consesus sequence (DCS) mutations as 
+generated by the **Variant Annotator** tool.
+
+**Dataset 2:** BAM file of aligned DCS reads. This file can be obtained by the 
+tool `Map with BWA-MEM <https://arxiv.org/abs/1303.3997>`_.
+
+**Dataset 3:** Tabular file with reads as produced by the 
+**Du Novo: Align families** tool of the `Du Novo Analysis Pipeline 
+<https://doi.org/10.1186/s13059-016-1039-4>`_
+
+**Output**
+
+The output is a json file containing dictonaries of the tags of reads containing mutations 
+in the DCS and a fastq file of all reads of these tags.
+
+    ]]> 
+    </help>
+    <expand macro="citation" />
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mut2sscs.py	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,133 @@
+#!/usr/bin/env python
+
+"""mut2sscs.py
+
+Author -- Gundula Povysil
+Contact -- povysil@bioinf.jku.at
+
+Takes a tabular file with mutations from DCS and a BAM file of SSCS as input
+and extracts all tags of reads that carry the mutation.
+Calculates statistics about number of ab/ba/duplex per mutation.
+
+=======  ==========  =================  ================================
+Version  Date        Author             Description
+0.2.1    2019-10-27  Gundula Povysil    -
+=======  ==========  =================  ================================
+
+USAGE: python mut2sscs.py DCS_Mutations.tabular SSCS.bam SSCS_counts.json
+
+"""
+
+from __future__ import division
+
+import argparse
+import json
+import os
+import sys
+
+import numpy as np
+import pysam
+
+
+def make_argparser():
+    parser = argparse.ArgumentParser(description='Takes a tabular file with mutations and a BAM file as input and prints all tags of reads that carry the mutation to a user specified output file.')
+    parser.add_argument('--mutFile',
+                        help='TABULAR file with DCS mutations.')
+    parser.add_argument('--bamFile',
+                        help='BAM file with aligned SSCS reads.')
+    parser.add_argument('--outputJson',
+                        help='Output JSON file to store SSCS counts.')
+    return parser
+
+
+def mut2sscs(argv):
+    parser = make_argparser()
+    args = parser.parse_args(argv[1:])
+
+    file1 = args.mutFile
+    file2 = args.bamFile
+    sscs_counts_json = args.outputJson
+
+    if os.path.isfile(file1) is False:
+        sys.exit("Error: Could not find '{}'".format(file1))
+
+    if os.path.isfile(file2) is False:
+        sys.exit("Error: Could not find '{}'".format(file2))
+
+    # 1. read mut file
+    with open(file1, 'r') as mut:
+        mut_array = np.genfromtxt(mut, skip_header=1, delimiter='\t', comments='#', dtype=str)
+
+    # 2 read SSCS bam file
+    # pysam.index(file2)
+    bam = pysam.AlignmentFile(file2, "rb")
+
+    # get tags
+    mut_pos_dict = {}
+    ref_pos_dict = {}
+    if mut_array.shape == (1,13):
+        mut_array = mut_array.reshape((1, len(mut_array)))
+
+    for m in range(0, len(mut_array[:, 0])):
+        print(str(m + 1) + " of " + str(len(mut_array[:, 0])))
+        chrom = mut_array[m, 1]
+        stop_pos = mut_array[m, 2].astype(int)
+        chrom_stop_pos = str(chrom) + "#" + str(stop_pos)
+        ref = mut_array[m, 9]
+        alt = mut_array[m, 10]
+
+        for pileupcolumn in bam.pileup(chrom.tostring(), stop_pos - 2, stop_pos, max_depth=1000000000):
+            if pileupcolumn.reference_pos == stop_pos - 1:
+                count_alt = 0
+                count_ref = 0
+                count_indel = 0
+                print("unfiltered reads=", pileupcolumn.n, "filtered reads=", len(pileupcolumn.pileups),
+                      "difference= ", len(pileupcolumn.pileups) - pileupcolumn.n)
+                for pileupread in pileupcolumn.pileups:
+                    if not pileupread.is_del and not pileupread.is_refskip:
+                        tag = pileupread.alignment.query_name
+                        abba = tag[-2:]
+                        # query position is None if is_del or is_refskip is set.
+                        if pileupread.alignment.query_sequence[pileupread.query_position] == alt:
+                            count_alt += 1
+                            if chrom_stop_pos in mut_pos_dict:
+                                if abba in mut_pos_dict[chrom_stop_pos]:
+                                    mut_pos_dict[chrom_stop_pos][abba] += 1
+                                else:
+                                    mut_pos_dict[chrom_stop_pos][abba] = 1
+                            else:
+                                mut_pos_dict[chrom_stop_pos] = {}
+                                mut_pos_dict[chrom_stop_pos][abba] = 1
+                        elif pileupread.alignment.query_sequence[pileupread.query_position] == ref:
+                            count_ref += 1
+                            if chrom_stop_pos in ref_pos_dict:
+                                if abba in ref_pos_dict[chrom_stop_pos]:
+                                    ref_pos_dict[chrom_stop_pos][abba] += 1
+                                else:
+                                    ref_pos_dict[chrom_stop_pos][abba] = 1
+                            else:
+                                ref_pos_dict[chrom_stop_pos] = {}
+                                ref_pos_dict[chrom_stop_pos][abba] = 1
+                        else:
+                            count_indel += 1
+
+                print("coverage at pos %s = %s, ref = %s, alt = %s, indel = %s,\n" %
+                      (pileupcolumn.pos, count_ref + count_alt, count_ref, count_alt, count_indel))
+
+        # if mutation is in DCS file but not in SSCS, then set counts to NA
+        if chrom_stop_pos not in mut_pos_dict.keys():
+            mut_pos_dict[chrom_stop_pos] = {}
+            mut_pos_dict[chrom_stop_pos]["ab"] = 0
+            mut_pos_dict[chrom_stop_pos]["ba"] = 0
+            ref_pos_dict[chrom_stop_pos] = {}
+            ref_pos_dict[chrom_stop_pos]["ab"] = 0
+            ref_pos_dict[chrom_stop_pos]["ba"] = 0
+    bam.close()
+
+    # save counts
+    with open(sscs_counts_json, "w") as f:
+        json.dump((mut_pos_dict, ref_pos_dict), f)
+
+
+if __name__ == '__main__':
+    sys.exit(mut2sscs(sys.argv))
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/mut2sscs.xml	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,59 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<tool id="mut2sscs" name="DCS mutations to SSCS stats:" version="1.0.1" profile="19.01">
+    <description>Extracts all tags from the single stranded consensus sequence (SSCS) bam file that carry a mutation at the same position a mutation is called in the duplex consensus sequence (DCS) and calculates their frequencies</description>
+    <macros>
+        <import>va_macros.xml</import>
+    </macros>
+    <expand macro="requirements"/>
+    <command><![CDATA[
+        ln -s '$file2' bam_input.bam &&
+        ln -s '${file2.metadata.bam_index}' bam_input.bam.bai &&
+        python '$__tool_directory__/mut2sscs.py' 
+        --mutFile '$file1'
+        --bamFile bam_input.bam
+        --outputJson '$output_json'
+    ]]>
+    </command>
+    <inputs>
+        <param name="file1" type="data" format="tabular" label="DCS Mutation File" optional="false" help="TABULAR file with DCS mutations. See Help section below for a detailed explanation."/>
+        <param name="file2" type="data" format="bam" label="SSCS BAM File" optional="false" help="BAM file with aligned SSCS reads."/>
+    </inputs>
+    <outputs>
+        <data name="output_json" format="json" label="${tool.name} on ${on_string}: JSON"/>
+    </outputs>
+    <tests>
+        <test>
+            <param name="file1" value="DCS_Mutations_test_data_VA.tabular"/>
+            <param name="file2" value="SSCS_test_data_VA.bam"/>
+            <output name="output_json" file="SSCS_counts_test_data_VA.json" lines_diff="2"/>
+        </test>
+    </tests>
+    <help> <![CDATA[
+**What it does**
+
+Takes a tabular file with DCS mutations and a BAM file of aligned SSCS reads 
+as input and writes statistics about tags of reads that carry a mutation in the 
+SSCS at the same position a mutation is called in the DCS to a user specified output file..
+
+**Input** 
+
+**Dataset 1:** Tabular file with duplex consesus sequence (DCS) mutations as 
+generated by the **Variant Annotator** tool.
+
+**Dataset 2:** BAM file of aligned single stranded consensus sequence (SSCS) 
+reads. This file can be obtained by the tool `Map with BWA-MEM 
+<https://arxiv.org/abs/1303.3997>`_.
+
+**Dataset 3:** Tabular file with reads as produced by the 
+**Du Novo: Align families** tool of the `Du Novo Analysis Pipeline 
+<https://doi.org/10.1186/s13059-016-1039-4>`_
+
+**Output**
+
+The output is a json file containing dictonaries with stats of tags that carry a mutation in the SSCS 
+at the same position a mutation is called in the DCS.
+
+    ]]> 
+    </help>
+    <expand macro="citation" />
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/read2mut.py	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,953 @@
+#!/usr/bin/env python
+
+"""read2mut.py
+
+Author -- Gundula Povysil
+Contact -- povysil@bioinf.jku.at
+
+Looks for reads with mutation at known
+positions and calculates frequencies and stats.
+
+=======  ==========  =================  ================================
+Version  Date        Author             Description
+0.2.1    2019-10-27  Gundula Povysil    -
+=======  ==========  =================  ================================
+
+
+USAGE: python read2mut.py --mutFile DCS_Mutations.tabular --bamFile Interesting_Reads.trim.bam
+                          --inputJson tag_count_dict.json --sscsJson SSCS_counts.json
+                          --outputFile mutant_reads_summary_short_trim.xlsx --thresh 10 --phred 20 --trim 10 --chimera_correction
+
+"""
+
+from __future__ import division
+
+import argparse
+import json
+import operator
+import os
+import re
+import sys
+
+import numpy as np
+import pysam
+import xlsxwriter
+
+
+def make_argparser():
+    parser = argparse.ArgumentParser(description='Takes a tabular file with mutations, a BAM file and JSON files as input and prints stats about variants to a user specified output file.')
+    parser.add_argument('--mutFile',
+                        help='TABULAR file with DCS mutations.')
+    parser.add_argument('--bamFile',
+                        help='BAM file with aligned raw reads of selected tags (FASTQ created by mut2read.py - trimming with Trimmomatic - alignment with bwa).')
+    parser.add_argument('--inputJson',
+                        help='JSON file with data collected by mut2read.py.')
+    parser.add_argument('--sscsJson',
+                        help='JSON file with SSCS counts collected by mut2sscs.py.')
+    parser.add_argument('--outputFile',
+                        help='Output xlsx file of mutation details.')
+    parser.add_argument('--thresh', type=int, default=0,
+                        help='Integer threshold for displaying mutations. Only mutations occuring less than thresh times are displayed. Default of 0 displays all.')
+    parser.add_argument('--phred', type=int, default=20,
+                        help='Integer threshold for Phred score. Only reads higher than this threshold are considered. Default 20.')
+    parser.add_argument('--trim', type=int, default=10,
+                        help='Integer threshold for assigning mutations at start and end of reads to lower tier. Default 10.')
+    parser.add_argument('--chimera_correction', action="store_true",
+                        help='Add additional tier for chimeric variants and correct the variant frequencies')
+
+    return parser
+
+
+def safe_div(x, y):
+    if y == 0:
+        return None
+    return x / y
+
+
+def read2mut(argv):
+    parser = make_argparser()
+    args = parser.parse_args(argv[1:])
+    file1 = args.mutFile
+    file2 = args.bamFile
+    json_file = args.inputJson
+    sscs_json = args.sscsJson
+    outfile = args.outputFile
+    thresh = args.thresh
+    phred_score = args.phred
+    trim = args.trim
+    chimera_correction = args.chimera_correction
+
+    if os.path.isfile(file1) is False:
+        sys.exit("Error: Could not find '{}'".format(file1))
+    if os.path.isfile(file2) is False:
+        sys.exit("Error: Could not find '{}'".format(file2))
+    if os.path.isfile(json_file) is False:
+        sys.exit("Error: Could not find '{}'".format(json_file))
+    if thresh < 0:
+        sys.exit("Error: thresh is '{}', but only non-negative integers allowed".format(thresh))
+    if phred_score < 0:
+        sys.exit("Error: phred is '{}', but only non-negative integers allowed".format(phred_score))
+    if trim < 0:
+        sys.exit("Error: trim is '{}', but only non-negative integers allowed".format(thresh))
+
+    # 1. read mut file
+    with open(file1, 'r') as mut:
+        mut_array = np.genfromtxt(mut, skip_header=1, delimiter='\t', comments='#', dtype=str)
+
+    # 2. load dicts
+    with open(json_file, "r") as f:
+        (tag_dict, cvrg_dict) = json.load(f)
+
+    with open(sscs_json, "r") as f:
+        (mut_pos_dict, ref_pos_dict) = json.load(f)
+
+    # 3. read bam file
+    # pysam.index(file2)
+    bam = pysam.AlignmentFile(file2, "rb")
+
+    # 4. create mut_dict
+    mut_dict = {}
+    mut_read_pos_dict = {}
+    mut_read_dict = {}
+    reads_dict = {}
+    if mut_array.shape == (1, 13):
+        mut_array = mut_array.reshape((1, len(mut_array)))
+
+    for m in range(0, len(mut_array[:, 0])):
+        print(str(m + 1) + " of " + str(len(mut_array[:, 0])))
+        #    for m in range(0, 5):
+        chrom = mut_array[m, 1]
+        stop_pos = mut_array[m, 2].astype(int)
+        chrom_stop_pos = str(chrom) + "#" + str(stop_pos)
+        ref = mut_array[m, 9]
+        alt = mut_array[m, 10]
+        mut_dict[chrom_stop_pos] = {}
+        mut_read_pos_dict[chrom_stop_pos] = {}
+        reads_dict[chrom_stop_pos] = {}
+
+        for pileupcolumn in bam.pileup(chrom.tostring(), stop_pos - 2, stop_pos, max_depth=1000000000):
+            if pileupcolumn.reference_pos == stop_pos - 1:
+                count_alt = 0
+                count_ref = 0
+                count_indel = 0
+                count_n = 0
+                count_other = 0
+                count_lowq = 0
+                n = 0
+                print("unfiltered reads=", pileupcolumn.n, "filtered reads=", len(pileupcolumn.pileups),
+                      "difference= ", len(pileupcolumn.pileups) - pileupcolumn.n)
+                for pileupread in pileupcolumn.pileups:
+                    n += 1
+                    if not pileupread.is_del and not pileupread.is_refskip:
+                        tag = pileupread.alignment.query_name
+                        nuc = pileupread.alignment.query_sequence[pileupread.query_position]
+                        phred = ord(pileupread.alignment.qual[pileupread.query_position]) - 33
+                        if phred < phred_score:
+                            nuc = "lowQ"
+                        if tag not in mut_dict[chrom_stop_pos]:
+                            mut_dict[chrom_stop_pos][tag] = {}
+                        if nuc in mut_dict[chrom_stop_pos][tag]:
+                            mut_dict[chrom_stop_pos][tag][nuc] += 1
+                        else:
+                            mut_dict[chrom_stop_pos][tag][nuc] = 1
+                        if tag not in mut_read_pos_dict[chrom_stop_pos]:
+                            mut_read_pos_dict[chrom_stop_pos][tag] = np.array(pileupread.query_position) + 1
+                            reads_dict[chrom_stop_pos][tag] = len(pileupread.alignment.query_sequence)
+                        else:
+                            mut_read_pos_dict[chrom_stop_pos][tag] = np.append(
+                                mut_read_pos_dict[chrom_stop_pos][tag], pileupread.query_position + 1)
+                            reads_dict[chrom_stop_pos][tag] = np.append(
+                                reads_dict[chrom_stop_pos][tag], len(pileupread.alignment.query_sequence))
+
+                        if nuc == alt:
+                            count_alt += 1
+                            if tag not in mut_read_dict:
+                                mut_read_dict[tag] = {}
+                                mut_read_dict[tag][chrom_stop_pos] = alt
+                            else:
+                                mut_read_dict[tag][chrom_stop_pos] = alt
+                        elif nuc == ref:
+                            count_ref += 1
+                        elif nuc == "N":
+                            count_n += 1
+                        elif nuc == "lowQ":
+                            count_lowq += 1
+                        else:
+                            count_other += 1
+                    else:
+                        count_indel += 1
+
+                print("coverage at pos %s = %s, ref = %s, alt = %s, other bases = %s, N = %s, indel = %s, low quality = %s\n" % (pileupcolumn.pos, count_ref + count_alt, count_ref, count_alt, count_other, count_n, count_indel, count_lowq))
+
+    for read in bam.fetch(until_eof=True):
+        if read.is_unmapped:
+            pure_tag = read.query_name[:-5]
+            nuc = "na"
+            for key in tag_dict[pure_tag].keys():
+                if key not in mut_dict:
+                    mut_dict[key] = {}
+                if read.query_name not in mut_dict[key]:
+                    mut_dict[key][read.query_name] = {}
+                if nuc in mut_dict[key][read.query_name]:
+                    mut_dict[key][read.query_name][nuc] += 1
+                else:
+                    mut_dict[key][read.query_name][nuc] = 1
+    bam.close()
+
+    # 5. create pure_tags_dict
+    pure_tags_dict = {}
+    for key1, value1 in sorted(mut_dict.items()):
+        i = np.where(np.array(['#'.join(str(i) for i in z)
+                               for z in zip(mut_array[:, 1], mut_array[:, 2])]) == key1)[0][0]
+        ref = mut_array[i, 9]
+        alt = mut_array[i, 10]
+        pure_tags_dict[key1] = {}
+        for key2, value2 in sorted(value1.items()):
+            for key3, value3 in value2.items():
+                pure_tag = key2[:-5]
+                if key3 == alt:
+                    if pure_tag in pure_tags_dict[key1]:
+                        pure_tags_dict[key1][pure_tag] += 1
+                    else:
+                        pure_tags_dict[key1][pure_tag] = 1
+
+    # 6. create pure_tags_dict_short with thresh
+    if thresh > 0:
+        pure_tags_dict_short = {}
+        for key, value in sorted(pure_tags_dict.items()):
+            if len(value) < thresh:
+                pure_tags_dict_short[key] = value
+    else:
+        pure_tags_dict_short = pure_tags_dict
+
+    #whole_array = []
+    #for k in pure_tags_dict.values():
+    #    whole_array.extend(k.keys())
+
+    # 7. output summary with threshold
+    workbook = xlsxwriter.Workbook(outfile)
+    ws1 = workbook.add_worksheet("Results")
+    ws2 = workbook.add_worksheet("Allele frequencies")
+    ws3 = workbook.add_worksheet("Tiers")
+
+    format1 = workbook.add_format({'bg_color': '#BCF5A9'})  # green
+    format2 = workbook.add_format({'bg_color': '#FFC7CE'})  # red
+    format3 = workbook.add_format({'bg_color': '#FACC2E'})  # yellow
+
+    header_line = ('variant ID', 'tier', 'tag', 'mate', 'read pos.ab', 'read pos.ba', 'read median length.ab',
+                   'read median length.ba', 'DCS median length',
+                   'FS.ab', 'FS.ba', 'FSqc.ab', 'FSqc.ba', 'ref.ab', 'ref.ba', 'alt.ab', 'alt.ba',
+                   'rel. ref.ab', 'rel. ref.ba', 'rel. alt.ab', 'rel. alt.ba',
+                   'na.ab', 'na.ba', 'lowq.ab', 'lowq.ba', 'trim.ab', 'trim.ba',
+                   'SSCS alt.ab', 'SSCS alt.ba', 'SSCS ref.ab', 'SSCS ref.ba',
+                   'other mut', 'chimeric tag')
+    ws1.write_row(0, 0, header_line)
+
+    counter_tier11 = 0
+    counter_tier12 = 0
+    counter_tier21 = 0
+    counter_tier22 = 0
+    counter_tier23 = 0
+    counter_tier24 = 0
+    counter_tier31 = 0
+    counter_tier32 = 0
+    counter_tier41 = 0
+    counter_tier42 = 0
+
+    if chimera_correction:
+        counter_tier43 = 0
+
+    counter_tier5 = 0
+
+    row = 1
+    tier_dict = {}
+    for key1, value1 in sorted(mut_dict.items()):
+        counts_mut = 0
+        if key1 in pure_tags_dict_short.keys():
+            i = np.where(np.array(['#'.join(str(i) for i in z)
+                                   for z in zip(mut_array[:, 1], mut_array[:, 2])]) == key1)[0][0]
+            ref = mut_array[i, 9]
+            alt = mut_array[i, 10]
+            dcs_median = cvrg_dict[key1][2]
+            whole_array = pure_tags_dict_short[key1].keys()
+
+            tier_dict[key1] = {}
+            if chimera_correction:
+                values_tier_dict = [("tier 1.1", 0), ("tier 1.2", 0), ("tier 2.1", 0), ("tier 2.2", 0), ("tier 2.3", 0), ("tier 2.4", 0), ("tier 3.1", 0), ("tier 3.2", 0), ("tier 4.1", 0), ("tier 4.2", 0), ("tier 4.3", 0), ("tier 5", 0)]
+            else:
+                values_tier_dict = [("tier 1.1", 0), ("tier 1.2", 0), ("tier 2.1", 0), ("tier 2.2", 0), ("tier 2.3", 0), ("tier 2.4", 0), ("tier 3.1", 0), ("tier 3.2", 0), ("tier 4.1", 0), ("tier 4.2", 0), ("tier 5", 0)]
+            
+            for k, v in values_tier_dict:
+                tier_dict[key1][k] = v
+
+            used_keys = []
+            if 'ab' in mut_pos_dict[key1].keys():
+                sscs_mut_ab = mut_pos_dict[key1]['ab']
+            else:
+                sscs_mut_ab = 0
+            if 'ba' in mut_pos_dict[key1].keys():
+                sscs_mut_ba = mut_pos_dict[key1]['ba']
+            else:
+                sscs_mut_ba = 0
+            if 'ab' in ref_pos_dict[key1].keys():
+                sscs_ref_ab = ref_pos_dict[key1]['ab']
+            else:
+                sscs_ref_ab = 0
+            if 'ba' in ref_pos_dict[key1].keys():
+                sscs_ref_ba = ref_pos_dict[key1]['ba']
+            else:
+                sscs_ref_ba = 0
+            for key2, value2 in sorted(value1.items()):
+                add_mut14 = ""
+                add_mut23 = ""
+                if (key2[:-5] in pure_tags_dict_short[key1].keys()) and (key2[:-5] not in used_keys) and (key1 in tag_dict[key2[:-5]].keys()):
+                    if key2[:-5] + '.ab.1' in mut_dict[key1].keys():
+                        total1 = sum(mut_dict[key1][key2[:-5] + '.ab.1'].values())
+                        if 'na' in mut_dict[key1][key2[:-5] + '.ab.1'].keys():
+                            na1 = mut_dict[key1][key2[:-5] + '.ab.1']['na']
+                            # na1f = na1/total1
+                        else:
+                            # na1 = na1f = 0
+                            na1 = 0
+                        if 'lowQ' in mut_dict[key1][key2[:-5] + '.ab.1'].keys():
+                            lowq1 = mut_dict[key1][key2[:-5] + '.ab.1']['lowQ']
+                            # lowq1f = lowq1 / total1
+                        else:
+                            # lowq1 = lowq1f = 0
+                            lowq1 = 0
+                        if ref in mut_dict[key1][key2[:-5] + '.ab.1'].keys():
+                            ref1 = mut_dict[key1][key2[:-5] + '.ab.1'][ref]
+                            ref1f = ref1 / (total1 - na1 - lowq1)
+                        else:
+                            ref1 = ref1f = 0
+                        if alt in mut_dict[key1][key2[:-5] + '.ab.1'].keys():
+                            alt1 = mut_dict[key1][key2[:-5] + '.ab.1'][alt]
+                            alt1f = alt1 / (total1 - na1 - lowq1)
+                        else:
+                            alt1 = alt1f = 0
+                        total1new = total1 - na1 - lowq1
+                        if (key2[:-5] + '.ab.1') in mut_read_dict.keys():
+                            k1 = mut_read_dict[(key2[:-5] + '.ab.1')].keys()
+                            add_mut1 = len(k1)
+                            if add_mut1 > 1:
+                                for k, v in mut_read_dict[(key2[:-5] + '.ab.1')].items():
+                                    if k != key1:
+                                        if len(add_mut14) == 0:
+                                            add_mut14 = str(k) + "_" + v
+                                        else:
+                                            add_mut14 = add_mut14 + ", " + str(k) + "_" + v
+                        else:
+                            k1 = []
+                    else:
+                        total1 = total1new = na1 = lowq1 = 0
+                        ref1 = alt1 = ref1f = alt1f = 0
+                        k1 = []
+
+                    if key2[:-5] + '.ab.2' in mut_dict[key1].keys():
+                        total2 = sum(mut_dict[key1][key2[:-5] + '.ab.2'].values())
+                        if 'na' in mut_dict[key1][key2[:-5] + '.ab.2'].keys():
+                            na2 = mut_dict[key1][key2[:-5] + '.ab.2']['na']
+                            # na2f = na2 / total2
+                        else:
+                            # na2 = na2f = 0
+                            na2 = 0
+                        if 'lowQ' in mut_dict[key1][key2[:-5] + '.ab.2'].keys():
+                            lowq2 = mut_dict[key1][key2[:-5] + '.ab.2']['lowQ']
+                            # lowq2f = lowq2 / total2
+                        else:
+                            # lowq2 = lowq2f = 0
+                            lowq2 = 0
+                        if ref in mut_dict[key1][key2[:-5] + '.ab.2'].keys():
+                            ref2 = mut_dict[key1][key2[:-5] + '.ab.2'][ref]
+                            ref2f = ref2 / (total2 - na2 - lowq2)
+                        else:
+                            ref2 = ref2f = 0
+                        if alt in mut_dict[key1][key2[:-5] + '.ab.2'].keys():
+                            alt2 = mut_dict[key1][key2[:-5] + '.ab.2'][alt]
+                            alt2f = alt2 / (total2 - na2 - lowq2)
+                        else:
+                            alt2 = alt2f = 0
+                        total2new = total2 - na2 - lowq2
+                        if (key2[:-5] + '.ab.2') in mut_read_dict.keys():
+                            k2 = mut_read_dict[(key2[:-5] + '.ab.2')].keys()
+                            add_mut2 = len(k2)
+                            if add_mut2 > 1:
+                                for k, v in mut_read_dict[(key2[:-5] + '.ab.2')].items():
+                                    if k != key1:
+                                        if len(add_mut23) == 0:
+                                            add_mut23 = str(k) + "_" + v
+                                        else:
+                                            add_mut23 = add_mut23 + ", " + str(k) + "_" + v
+                        else:
+                            k2 = []
+                    else:
+                        total2 = total2new = na2 = lowq2 = 0
+                        ref2 = alt2 = ref2f = alt2f = 0
+                        k2 = []
+
+                    if key2[:-5] + '.ba.1' in mut_dict[key1].keys():
+                        total3 = sum(mut_dict[key1][key2[:-5] + '.ba.1'].values())
+                        if 'na' in mut_dict[key1][key2[:-5] + '.ba.1'].keys():
+                            na3 = mut_dict[key1][key2[:-5] + '.ba.1']['na']
+                            # na3f = na3 / total3
+                        else:
+                            # na3 = na3f = 0
+                            na3 = 0
+                        if 'lowQ' in mut_dict[key1][key2[:-5] + '.ba.1'].keys():
+                            lowq3 = mut_dict[key1][key2[:-5] + '.ba.1']['lowQ']
+                            # lowq3f = lowq3 / total3
+                        else:
+                            # lowq3 = lowq3f = 0
+                            lowq3 = 0
+                        if ref in mut_dict[key1][key2[:-5] + '.ba.1'].keys():
+                            ref3 = mut_dict[key1][key2[:-5] + '.ba.1'][ref]
+                            ref3f = ref3 / (total3 - na3 - lowq3)
+                        else:
+                            ref3 = ref3f = 0
+                        if alt in mut_dict[key1][key2[:-5] + '.ba.1'].keys():
+                            alt3 = mut_dict[key1][key2[:-5] + '.ba.1'][alt]
+                            alt3f = alt3 / (total3 - na3 - lowq3)
+                        else:
+                            alt3 = alt3f = 0
+                        total3new = total3 - na3 - lowq3
+                        if (key2[:-5] + '.ba.1') in mut_read_dict.keys():
+                            add_mut3 = len(mut_read_dict[(key2[:-5] + '.ba.1')].keys())
+                            if add_mut3 > 1:
+                                for k, v in mut_read_dict[(key2[:-5] + '.ba.1')].items():
+                                    if k != key1 and k not in k2:
+                                        if len(add_mut23) == 0:
+                                            add_mut23 = str(k) + "_" + v
+                                        else:
+                                            add_mut23 = add_mut23 + ", " + str(k) + "_" + v
+                    else:
+                        total3 = total3new = na3 = lowq3 = 0
+                        ref3 = alt3 = ref3f = alt3f = 0
+
+                    if key2[:-5] + '.ba.2' in mut_dict[key1].keys():
+                        total4 = sum(mut_dict[key1][key2[:-5] + '.ba.2'].values())
+                        if 'na' in mut_dict[key1][key2[:-5] + '.ba.2'].keys():
+                            na4 = mut_dict[key1][key2[:-5] + '.ba.2']['na']
+                            # na4f = na4 / total4
+                        else:
+                            # na4 = na4f = 0
+                            na4 = 0
+                        if 'lowQ' in mut_dict[key1][key2[:-5] + '.ba.2'].keys():
+                            lowq4 = mut_dict[key1][key2[:-5] + '.ba.2']['lowQ']
+                            # lowq4f = lowq4 / total4
+                        else:
+                            # lowq4 = lowq4f = 0
+                            lowq4 = 0
+                        if ref in mut_dict[key1][key2[:-5] + '.ba.2'].keys():
+                            ref4 = mut_dict[key1][key2[:-5] + '.ba.2'][ref]
+                            ref4f = ref4 / (total4 - na4 - lowq4)
+                        else:
+                            ref4 = ref4f = 0
+                        if alt in mut_dict[key1][key2[:-5] + '.ba.2'].keys():
+                            alt4 = mut_dict[key1][key2[:-5] + '.ba.2'][alt]
+                            alt4f = alt4 / (total4 - na4 - lowq4)
+                        else:
+                            alt4 = alt4f = 0
+                        total4new = total4 - na4 - lowq4
+                        if (key2[:-5] + '.ba.2') in mut_read_dict.keys():
+                            add_mut4 = len(mut_read_dict[(key2[:-5] + '.ba.2')].keys())
+                            if add_mut4 > 1:
+                                for k, v in mut_read_dict[(key2[:-5] + '.ba.2')].items():
+                                    if k != key1 and k not in k1:
+                                        if len(add_mut14) == 0:
+                                            add_mut14 = str(k) + "_" + v
+                                        else:
+                                            add_mut14 = add_mut14 + ", " + str(k) + "_" + v
+                    else:
+                        total4 = total4new = na4 = lowq4 = 0
+                        ref4 = alt4 = ref4f = alt4f = 0
+
+                    read_pos1 = read_pos2 = read_pos3 = read_pos4 = -1
+                    read_len_median1 = read_len_median2 = read_len_median3 = read_len_median4 = 0
+
+                    if key2[:-5] + '.ab.1' in mut_read_pos_dict[key1].keys():
+                        read_pos1 = np.median(mut_read_pos_dict[key1][key2[:-5] + '.ab.1'])
+                        read_len_median1 = np.median(reads_dict[key1][key2[:-5] + '.ab.1'])
+                    if key2[:-5] + '.ab.2' in mut_read_pos_dict[key1].keys():
+                        read_pos2 = np.median(mut_read_pos_dict[key1][key2[:-5] + '.ab.2'])
+                        read_len_median2 = np.median(reads_dict[key1][key2[:-5] + '.ab.2'])
+                    if key2[:-5] + '.ba.1' in mut_read_pos_dict[key1].keys():
+                        read_pos3 = np.median(mut_read_pos_dict[key1][key2[:-5] + '.ba.1'])
+                        read_len_median3 = np.median(reads_dict[key1][key2[:-5] + '.ba.1'])
+                    if key2[:-5] + '.ba.2' in mut_read_pos_dict[key1].keys():
+                        read_pos4 = np.median(mut_read_pos_dict[key1][key2[:-5] + '.ba.2'])
+                        read_len_median4 = np.median(reads_dict[key1][key2[:-5] + '.ba.2'])
+
+                    used_keys.append(key2[:-5])
+                    counts_mut += 1
+                    if (alt1f + alt2f + alt3f + alt4f) > 0.5:
+                        if total1new == 0:
+                            ref1f = alt1f = None
+                            alt1ff = -1
+                        else:
+                            alt1ff = alt1f
+                        if total2new == 0:
+                            ref2f = alt2f = None
+                            alt2ff = -1
+                        else:
+                            alt2ff = alt2f
+                        if total3new == 0:
+                            ref3f = alt3f = None
+                            alt3ff = -1
+                        else:
+                            alt3ff = alt3f
+                        if total4new == 0:
+                            ref4f = alt4f = None
+                            alt4ff = -1
+                        else:
+                            alt4ff = alt4f
+
+                        beg1 = beg4 = beg2 = beg3 = 0
+
+                        details1 = (total1, total4, total1new, total4new, ref1, ref4, alt1, alt4, ref1f, ref4f, alt1f, alt4f, na1, na4, lowq1, lowq4, beg1, beg4)
+                        details2 = (total2, total3, total2new, total3new, ref2, ref3, alt2, alt3, ref2f, ref3f, alt2f, alt3f, na2, na3, lowq2, lowq3, beg2, beg3)
+                        
+                        chrom, pos = re.split(r'\#', key1)
+                        var_id = '-'.join([chrom, pos, ref, alt])
+                        sample_tag = key2[:-5]
+                        array2 = np.unique(whole_array)  # remove duplicate sequences to decrease running time
+                        # exclude identical tag from array2, to prevent comparison to itself
+                        same_tag = np.where(array2 == sample_tag)
+                        index_array2 = np.arange(0, len(array2), 1)
+                        index_withoutSame = np.delete(index_array2, same_tag)  # delete identical tag from the data
+                        array2 = array2[index_withoutSame]
+                        if len(array2) != 0:  # only perform chimera analysis if there is more than 1 variant
+                            array1_half = sample_tag[0:int(len(sample_tag) / 2)]  # mate1 part1
+                            array1_half2 = sample_tag[int(len(sample_tag) / 2):int(len(sample_tag))]  # mate1 part 2
+                            array2_half = np.array([ii[0:int(len(ii) / 2)] for ii in array2])  # mate2 part1
+                            array2_half2 = np.array([ii[int(len(ii) / 2):int(len(ii))] for ii in array2])  # mate2 part2
+
+                            min_tags_list_zeros = []
+                            chimera_tags = []
+                            for mate_b in [False, True]:
+                                i = 0  # counter, only used to see how many HDs of tags were already calculated
+                                if mate_b is False:  # HD calculation for all a's
+                                    half1_mate1 = array1_half
+                                    half2_mate1 = array1_half2
+                                    half1_mate2 = array2_half
+                                    half2_mate2 = array2_half2
+                                elif mate_b is True:  # HD calculation for all b's
+                                    half1_mate1 = array1_half2
+                                    half2_mate1 = array1_half
+                                    half1_mate2 = array2_half2
+                                    half2_mate2 = array2_half
+                                # calculate HD of "a" in the tag to all "a's" or "b" in the tag to all "b's"
+                                dist = np.array([sum(map(operator.ne, half1_mate1, c)) for c in half1_mate2])
+                                min_index = np.where(dist == dist.min())  # get index of min HD
+                                # get all "b's" of the tag or all "a's" of the tag with minimum HD
+                                min_tag_half2 = half2_mate2[min_index]
+                                min_tag_array2 = array2[min_index]  # get whole tag with min HD
+                                min_value = dist.min()
+                                # calculate HD of "b" to all "b's" or "a" to all "a's"
+                                dist_second_half = np.array([sum(map(operator.ne, half2_mate1, e))
+                                                             for e in min_tag_half2])
+
+                                dist2 = dist_second_half.max()
+                                max_index = np.where(dist_second_half == dist_second_half.max())[0]  # get index of max HD
+                                max_tag = min_tag_array2[max_index]
+
+                                # tags which have identical parts:
+                                if min_value == 0 or dist2 == 0:
+                                    min_tags_list_zeros.append(tag)
+                                    chimera_tags.append(max_tag)
+                                    # chimeric = True
+                                # else:
+                                    # chimeric = False
+
+                                # if mate_b is False:
+                                #    text = "pos {}: sample tag: {}; HD a = {}; HD b' = {}; similar tag(s): {}; chimeric = {}".format(pos, sample_tag, min_value, dist2, list(max_tag), chimeric)
+                                # else:
+                                #     text = "pos {}: sample tag: {}; HD a' = {}; HD b = {}; similar tag(s): {}; chimeric = {}".format(pos, sample_tag, dist2, min_value, list(max_tag), chimeric)
+                                i += 1
+                            chimera_tags = [x for x in chimera_tags if x != []]
+                            chimera_tags_new = []
+                            for i in chimera_tags:
+                                if len(i) > 1:
+                                    for t in i:
+                                        chimera_tags_new.append(t)
+                                else:
+                                    chimera_tags_new.extend(i)
+                            chimera_tags_new = np.asarray(chimera_tags_new)
+                            chimera = ", ".join(chimera_tags_new)
+                        else:
+                            chimera_tags_new = []
+                            chimera = ""
+
+                        trimmed = False
+                        contradictory = False
+                        chimeric_dcs = False
+
+                        if chimera_correction and len(chimera_tags_new) > 0: # chimeras
+                            alt1ff = 0
+                            alt4ff = 0
+                            alt2ff = 0
+                            alt3ff = 0
+                            chimeric_dcs = True
+                            trimmed = False
+                            contradictory = False
+                        elif ((all(float(ij) >= 0.5 for ij in [alt1ff, alt4ff]) & # contradictory variant
+                            all(float(ij) == 0. for ij in [alt2ff, alt3ff])) |
+                            (all(float(ij) >= 0.5 for ij in [alt2ff, alt3ff]) &
+                            all(float(ij) == 0. for ij in [alt1ff, alt4ff]))):
+                            alt1ff = 0
+                            alt4ff = 0
+                            alt2ff = 0
+                            alt3ff = 0
+                            trimmed = False
+                            contradictory = True
+                            chimeric_dcs = False
+                        else:
+                            if ((read_pos1 >= 0) and ((read_pos1 <= trim) | (abs(read_len_median1 - read_pos1) <= trim))):
+                                beg1 = total1new
+                                total1new = 0
+                                alt1ff = 0
+                                alt1f = 0
+                                trimmed = True
+
+                            if ((read_pos4 >= 0) and ((read_pos4 <= trim) | (abs(read_len_median4 - read_pos4) <= trim))):
+                                beg4 = total4new
+                                total4new = 0
+                                alt4ff = 0
+                                alt4f = 0
+                                trimmed = True
+                            
+                            if ((read_pos2 >= 0) and ((read_pos2 <= trim) | (abs(read_len_median2 - read_pos2) <= trim))):
+                                beg2 = total2new
+                                total2new = 0
+                                alt2ff = 0
+                                alt2f = 0
+                                trimmed = True
+                                
+                            if ((read_pos3 >= 0) and ((read_pos3 <= trim) | (abs(read_len_median3 - read_pos3) <= trim))):
+                                beg3 = total3new
+                                total3new = 0
+                                alt3ff = 0
+                                alt3f = 0
+                                trimmed = True
+                            details1 = (total1, total4, total1new, total4new, ref1, ref4, alt1, alt4, ref1f, ref4f, alt1f, alt4f, na1, na4, lowq1, lowq4, beg1, beg4)
+                            details2 = (total2, total3, total2new, total3new, ref2, ref3, alt2, alt3, ref2f, ref3f, alt2f, alt3f, na2, na3, lowq2, lowq3, beg2, beg3)
+
+                        # assign tiers
+                        if ((all(int(ij) >= 3 for ij in [total1new, total4new]) &
+                             all(float(ij) >= 0.75 for ij in [alt1ff, alt4ff])) |
+                            (all(int(ij) >= 3 for ij in [total2new, total3new]) &
+                             all(float(ij) >= 0.75 for ij in [alt2ff, alt3ff]))):
+                            tier = "1.1"
+                            counter_tier11 += 1
+                            tier_dict[key1]["tier 1.1"] += 1
+
+                        elif (all(int(ij) >= 1 for ij in [total1new, total2new, total3new, total4new]) &
+                              any(int(ij) >= 3 for ij in [total1new, total4new]) &
+                              any(int(ij) >= 3 for ij in [total2new, total3new]) &
+                              all(float(ij) >= 0.75 for ij in [alt1ff, alt2ff, alt3ff, alt4ff])):
+                            tier = "1.2"
+                            counter_tier12 += 1
+                            tier_dict[key1]["tier 1.2"] += 1
+
+                        elif ((all(int(ij) >= 1 for ij in [total1new, total4new]) &
+                               any(int(ij) >= 3 for ij in [total1new, total4new]) &
+                               all(float(ij) >= 0.75 for ij in [alt1ff, alt4ff])) |
+                              (all(int(ij) >= 1 for ij in [total2new, total3new]) &
+                               any(int(ij) >= 3 for ij in [total2new, total3new]) &
+                               all(float(ij) >= 0.75 for ij in [alt2ff, alt3ff]))):
+                            tier = "2.1"
+                            counter_tier21 += 1
+                            tier_dict[key1]["tier 2.1"] += 1
+
+                        elif (all(int(ij) >= 1 for ij in [total1new, total2new, total3new, total4new]) &
+                              all(float(ij) >= 0.75 for ij in [alt1ff, alt2ff, alt3ff, alt4ff])):
+                            tier = "2.2"
+                            counter_tier22 += 1
+                            tier_dict[key1]["tier 2.2"] += 1
+
+                        elif ((all(int(ij) >= 1 for ij in [total1new, total4new]) &
+                               any(int(ij) >= 3 for ij in [total2new, total3new]) &
+                               all(float(ij) >= 0.75 for ij in [alt1ff, alt4ff]) &
+                               any(float(ij) >= 0.75 for ij in [alt2ff, alt3ff])) |
+                              (all(int(ij) >= 1 for ij in [total2new, total3new]) &
+                               any(int(ij) >= 3 for ij in [total1new, total4new]) &
+                               all(float(ij) >= 0.75 for ij in [alt2ff, alt3ff]) &
+                               any(float(ij) >= 0.75 for ij in [alt1ff, alt4ff]))):
+                            tier = "2.3"
+                            counter_tier23 += 1
+                            tier_dict[key1]["tier 2.3"] += 1
+
+                        elif ((all(int(ij) >= 1 for ij in [total1new, total4new]) &
+                               all(float(ij) >= 0.75 for ij in [alt1ff, alt4ff])) |
+                              (all(int(ij) >= 1 for ij in [total2new, total3new]) &
+                               all(float(ij) >= 0.75 for ij in [alt2ff, alt3ff]))):
+                            tier = "2.4"
+                            counter_tier24 += 1
+                            tier_dict[key1]["tier 2.4"] += 1
+
+                        elif ((len(pure_tags_dict_short[key1]) > 1) &
+                              (all(float(ij) >= 0.5 for ij in [alt1ff, alt4ff]) |
+                               all(float(ij) >= 0.5 for ij in [alt2ff, alt3ff]))):
+                            tier = "3.1"
+                            counter_tier31 += 1
+                            tier_dict[key1]["tier 3.1"] += 1
+
+                        elif ((all(int(ij) >= 1 for ij in [total1new, total4new]) &
+                               all(float(ij) >= 0.5 for ij in [alt1ff, alt4ff])) |
+                              (all(int(ij) >= 1 for ij in [total2new, total3new]) &
+                               all(float(ij) >= 0.5 for ij in [alt2ff, alt3ff]))):
+                            tier = "3.2"
+                            counter_tier32 += 1
+                            tier_dict[key1]["tier 3.2"] += 1
+
+                        elif (trimmed):
+                            tier = "4.1"
+                            counter_tier41 += 1
+                            tier_dict[key1]["tier 4.1"] += 1
+
+                        elif (contradictory):
+                            tier = "4.2"
+                            counter_tier42 += 1
+                            tier_dict[key1]["tier 4.2"] += 1
+
+                        elif chimera_correction and chimeric_dcs:
+                            tier = "4.3"
+                            counter_tier43 += 1
+                            tier_dict[key1]["tier 4.3"] += 1
+
+                        else:
+                            tier = "5"
+                            counter_tier5 += 1
+                            tier_dict[key1]["tier 5"] += 1
+
+                        if (read_pos1 == -1):
+                            read_pos1 = read_len_median1 = None
+                        if (read_pos4 == -1):
+                            read_pos4 = read_len_median4 = None
+                        if (read_pos2 == -1):
+                            read_pos2 = read_len_median2 = None
+                        if (read_pos3 == -1):
+                            read_pos3 = read_len_median3 = None
+                        line = (var_id, tier, key2[:-5], 'ab1.ba2', read_pos1, read_pos4, read_len_median1, read_len_median4, dcs_median) + details1 + (sscs_mut_ab, sscs_mut_ba, sscs_ref_ab, sscs_ref_ba, add_mut14, chimera)
+                        ws1.write_row(row, 0, line)
+                        line = ("", "", key2[:-5], 'ab2.ba1', read_pos2, read_pos3, read_len_median2, read_len_median3, dcs_median) + details2 + (sscs_mut_ab, sscs_mut_ba, sscs_ref_ab, sscs_ref_ba, add_mut23, chimera)
+                        ws1.write_row(row + 1, 0, line)
+
+                        ws1.conditional_format('L{}:M{}'.format(row + 1, row + 2),
+                                               {'type': 'formula',
+                                                'criteria': '=OR($B${}="1.1", $B${}="1.2")'.format(row + 1, row + 1),
+                                                'format': format1,
+                                                'multi_range': 'L{}:M{} T{}:U{} B{}'.format(row + 1, row + 2, row + 1, row + 2, row + 1, row + 2)})
+                        ws1.conditional_format('L{}:M{}'.format(row + 1, row + 2),
+                                               {'type': 'formula',
+                                                'criteria': '=OR($B${}="2.1", $B${}="2.2", $B${}="2.3", $B${}="2.4")'.format(row + 1, row + 1, row + 1, row + 1),
+                                                'format': format3,
+                                                'multi_range': 'L{}:M{} T{}:U{} B{}'.format(row + 1, row + 2, row + 1, row + 2, row + 1, row + 2)})
+                        ws1.conditional_format('L{}:M{}'.format(row + 1, row + 2),
+                                               {'type': 'formula',
+                                                'criteria': '=$B${}>="3"'.format(row + 1),
+                                                'format': format2,
+                                                'multi_range': 'L{}:M{} T{}:U{} B{}'.format(row + 1, row + 2, row + 1, row + 2, row + 1, row + 2)})
+
+                        row += 3
+
+    if chimera_correction:
+        header_line2 = ('variant ID', 'cvrg', 'AC alt (all tiers)', 'AF  (all tiers)', 'cvrg (tiers 1.1-2.4)', 'AC alt (tiers 1.1-2.4)', 'AF (tiers 1.1-2.4)', 'AC alt (orginal DCS)', 'AF (original DCS)',
+                    'tier 1.1', 'tier 1.2', 'tier 2.1', 'tier 2.2', 'tier 2.3', 'tier 2.4',
+                    'tier 3.1', 'tier 3.2', 'tier 4.1', 'tier 4.2', 'tier 4.3', 'tier 5', 'AF 1.1-1.2', 'AF 1.1-2.1', 'AF 1.1-2.2',
+                    'AF 1.1-2.3', 'AF 1.1-2.4', 'AF 1.1-3.1', 'AF 1.1-3.2', 'AF 1.1-4.1', 'AF 1.1-4.2', 'AF 1.1-4.3', 'AF 1.1-5')
+    else:
+        header_line2 = ('variant ID', 'cvrg', 'AC alt (all tiers)', 'AF  (all tiers)', 'cvrg (tiers 1.1-2.4)', 'AC alt (tiers 1.1-2.4)', 'AF (tiers 1.1-2.4)', 'AC alt (orginal DCS)', 'AF (original DCS)',
+                    'tier 1.1', 'tier 1.2', 'tier 2.1', 'tier 2.2', 'tier 2.3', 'tier 2.4',
+                    'tier 3.1', 'tier 3.2', 'tier 4.1', 'tier 4.2', 'tier 5', 'AF 1.1-1.2', 'AF 1.1-2.1', 'AF 1.1-2.2',
+                    'AF 1.1-2.3', 'AF 1.1-2.4', 'AF 1.1-3.1', 'AF 1.1-3.2', 'AF 1.1-4.1', 'AF 1.1-4.2', 'AF 1.1-5')
+
+    ws2.write_row(0, 0, header_line2)
+    row = 0
+
+    for key1, value1 in sorted(tier_dict.items()):
+        if key1 in pure_tags_dict_short.keys():
+            i = np.where(np.array(['#'.join(str(i) for i in z)
+                                   for z in zip(mut_array[:, 1], mut_array[:, 2])]) == key1)[0][0]
+            ref = mut_array[i, 9]
+            alt = mut_array[i, 10]
+            chrom, pos = re.split(r'\#', key1)
+            ref_count = cvrg_dict[key1][0]
+            alt_count = cvrg_dict[key1][1]
+            cvrg = ref_count + alt_count
+
+            var_id = '-'.join([chrom, pos, ref, alt])
+            lst = [var_id, cvrg]
+            used_tiers = []
+            cum_af = []
+            for key2, value2 in sorted(value1.items()):
+                # calculate cummulative AF
+                used_tiers.append(value2)
+                if len(used_tiers) > 1:
+                    cum = safe_div(sum(used_tiers), cvrg)
+                    cum_af.append(cum)
+            lst.extend([sum(used_tiers), safe_div(sum(used_tiers), cvrg)])
+            if chimera_correction:
+                lst.extend([(cvrg - sum(used_tiers[-6:])), sum(used_tiers[0:6]), safe_div(sum(used_tiers[0:6]), (cvrg - sum(used_tiers[-6:]))), alt_count, safe_div(alt_count, cvrg)])
+            else:
+                lst.extend([(cvrg - sum(used_tiers[-5:])), sum(used_tiers[0:6]), safe_div(sum(used_tiers[0:6]), (cvrg - sum(used_tiers[-5:]))), alt_count, safe_div(alt_count, cvrg)])
+            lst.extend(used_tiers)
+            lst.extend(cum_af)
+            lst = tuple(lst)
+            ws2.write_row(row + 1, 0, lst)
+            ws2.conditional_format('J{}:K{}'.format(row + 2, row + 2), {'type': 'formula', 'criteria': '=$J$1="tier 1.1"', 'format': format1, 'multi_range': 'J{}:K{} J1:K1'.format(row + 2, row + 2)})
+            ws2.conditional_format('L{}:O{}'.format(row + 2, row + 2), {'type': 'formula', 'criteria': '=$L$1="tier 2.1"', 'format': format3, 'multi_range': 'L{}:O{} L1:O1'.format(row + 2, row + 2)})
+            if chimera_correction:
+                ws2.conditional_format('P{}:U{}'.format(row + 2, row + 2), {'type': 'formula', 'criteria': '=$P$1="tier 3.1"', 'format': format2, 'multi_range': 'P{}:U{} P1:U1'.format(row + 2, row + 2)})
+            else:
+                ws2.conditional_format('P{}:T{}'.format(row + 2, row + 2), {'type': 'formula', 'criteria': '=$P$1="tier 3.1"', 'format': format2, 'multi_range': 'P{}:T{} P1:T1'.format(row + 2, row + 2)})
+
+            row += 1
+
+    # sheet 3
+    if chimera_correction:
+        sheet3 = [("tier 1.1", counter_tier11), ("tier 1.2", counter_tier12), ("tier 2.1", counter_tier21),
+                  ("tier 2.2", counter_tier22), ("tier 2.3", counter_tier23), ("tier 2.4", counter_tier24),
+                  ("tier 3.1", counter_tier31), ("tier 3.2", counter_tier32), ("tier 4.1", counter_tier41), 
+                  ("tier 4.2", counter_tier42), ("tier 4.3", counter_tier43), ("tier 5", counter_tier5)]
+    else:
+        sheet3 = [("tier 1.1", counter_tier11), ("tier 1.2", counter_tier12), ("tier 2.1", counter_tier21),
+                  ("tier 2.2", counter_tier22), ("tier 2.3", counter_tier23), ("tier 2.4", counter_tier24),
+                  ("tier 3.1", counter_tier31), ("tier 3.2", counter_tier32), ("tier 4.1", counter_tier41), 
+                  ("tier 4.2", counter_tier42), ("tier 5", counter_tier5)]
+
+
+    header = ("tier", "count")
+    ws3.write_row(0, 0, header)
+
+    for i in range(len(sheet3)):
+        ws3.write_row(i + 1, 0, sheet3[i])
+        ws3.conditional_format('A{}:B{}'.format(i + 2, i + 2),
+                               {'type': 'formula',
+                                'criteria': '=OR($A${}="tier 1.1", $A${}="tier 1.2")'.format(i + 2, i + 2),
+                                'format': format1})
+        ws3.conditional_format('A{}:B{}'.format(i + 2, i + 2),
+                               {'type': 'formula',
+                                'criteria': '=OR($A${}="tier 2.1", $A${}="tier 2.2", $A${}="tier 2.3", $A${}="tier 2.4")'.format(i + 2, i + 2, i + 2, i + 2),
+                                'format': format3})
+        ws3.conditional_format('A{}:B{}'.format(i + 2, i + 2),
+                               {'type': 'formula',
+                                'criteria': '=$A${}>="3"'.format(i + 2),
+                                'format': format2})
+    if chimera_correction:
+        description_tiers = [("Tier 1.1", "both ab and ba SSCS present (>75% of the sites with alternative base) and minimal FS>=3 for both SSCS in at least one mate"), ("", ""), ("Tier 1.2", "both ab and ba SSCS present (>75% of the sites with alt. base) and mate pair validation (min. FS=1) and minimal FS>=3 for at least one of the SSCS"), ("Tier 2.1", "both ab and ba SSCS present (>75% of the sites with alt. base) and minimal FS>=3 for at least one of the SSCS in at least one mate"), ("Tier 2.2", "both ab and ba SSCS present (>75% of the sites with alt. base) and mate pair validation (min. FS=1)"), ("Tier 2.3", "both ab and ba SSCS present (>75% of the sites with alt. base) and minimal FS=1 for both SSCS in one mate and minimal FS>=3 for at least one of the SSCS in the other mate"), ("Tier 2.4", "both ab and ba SSCS present (>75% of the sites with alt. base) and minimal FS=1 for both SSCS in at least one mate"), ("Tier 3.1", "both ab and ba SSCS present (>50% of the sites with alt. base) and recurring mutation on this position"), ("Tier 3.2", "both ab and ba SSCS present (>50% of the sites with alt. base) and minimal FS>=1 for both SSCS in at least one mate"), ("Tier 4.1", "variants at the start or end of the reads"), ("Tier 4.2", "mates with contradictory information"), ("Tier 4.3", "variants that are chimeric"), ("Tier 5", "remaining variants")]
+    else:
+        description_tiers = [("Tier 1.1", "both ab and ba SSCS present (>75% of the sites with alternative base) and minimal FS>=3 for both SSCS in at least one mate"), ("", ""), ("Tier 1.2", "both ab and ba SSCS present (>75% of the sites with alt. base) and mate pair validation (min. FS=1) and minimal FS>=3 for at least one of the SSCS"), ("Tier 2.1", "both ab and ba SSCS present (>75% of the sites with alt. base) and minimal FS>=3 for at least one of the SSCS in at least one mate"), ("Tier 2.2", "both ab and ba SSCS present (>75% of the sites with alt. base) and mate pair validation (min. FS=1)"), ("Tier 2.3", "both ab and ba SSCS present (>75% of the sites with alt. base) and minimal FS=1 for both SSCS in one mate and minimal FS>=3 for at least one of the SSCS in the other mate"), ("Tier 2.4", "both ab and ba SSCS present (>75% of the sites with alt. base) and minimal FS=1 for both SSCS in at least one mate"), ("Tier 3.1", "both ab and ba SSCS present (>50% of the sites with alt. base) and recurring mutation on this position"), ("Tier 3.2", "both ab and ba SSCS present (>50% of the sites with alt. base) and minimal FS>=1 for both SSCS in at least one mate"), ("Tier 4.1", "variants at the start or end of the reads"), ("Tier 4.2", "mates with contradictory information"), ("Tier 5", "remaining variants")]
+
+    examples_tiers = [[("Chr5:5-20000-11068-C-G", "1.1", "AAAAAGATGCCGACTACCTT", "ab1.ba2", "254", "228", "287", "288", "289",
+                        "3", "6", "3", "6", "0", "0", "3", "6", "0", "0", "1", "1", "0", "0", "0", "0", "0", "0",
+                        "4081", "4098", "5", "10", "", ""),
+                       ("", "", "AAAAAGATGCCGACTACCTT", "ab2.ba1", None, None, None, None,
+                        "289", "0", "0", "0", "0", "0", "0", "0", "0", None, None, None, None,
+                        "0", "0", "0", "0", "0", "0", "4081", "4098", "5", "10", "", "")],
+                      [("Chr5:5-20000-11068-C-G", "1.1", "AAAAATGCGTAGAAATATGC", "ab1.ba2", "254", "228", "287", "288", "289",
+                        "33", "43", "33", "43", "0", "0", "33", "43", "0", "0", "1", "1", "0", "0", "0", "0", "0",
+                        "0", "4081", "4098", "5", "10", "", ""),
+                       ("", "", "AAAAATGCGTAGAAATATGC", "ab2.ba1", "268", "268", "270", "288", "289",
+                        "11", "34", "10", "27", "0", "0", "10", "27", "0", "0", "1", "1", "0", "0", "1",
+                        "7", "0", "0", "4081", "4098", "5", "10", "", "")],
+                      [("Chr5:5-20000-10776-G-T", "1.2", "CTATGACCCGTGAGCCCATG", "ab1.ba2", "132", "132", "287", "288", "290",
+                        "4", "1", "4", "1", "0", "0", "4", "1", "0", "0", "1", "1", "0", "0", "0", "0",
+                        "0", "0", "1", "6", "47170", "41149", "", ""),
+                       ("", "", "CTATGACCCGTGAGCCCATG", "ab2.ba1", "77", "132", "233", "200", "290",
+                        "4", "1", "4", "1", "0", "0", "4", "1", "0", "0", "1", "1", "0", "0", "0", "0",
+                        "0", "0", "1", "6", "47170", "41149", "", "")],
+                      [("Chr5:5-20000-11068-C-G", "2.1", "AAAAAAACATCATACACCCA", "ab1.ba2", "246", "244", "287", "288", "289",
+                        "2", "8", "2", "8", "0", "0", "2", "8", "0", "0", "1", "1", "0", "0", "0", "0", "0", "0",
+                        "4081", "4098", "5", "10", "", ""),
+                       ("", "", "AAAAAAACATCATACACCCA", "ab2.ba1", None, None, None, None,
+                        "289", "0", "0", "0", "0", "0", "0", "0", "0", None, None, None, None, "0", "0",
+                        "0", "0", "0", "0", "4081", "4098", "5", "10", "", "")],
+                      [("Chr5:5-20000-11068-C-G", "2.2", "ATCAGCCATGGCTATTATTG", "ab1.ba2", "72", "72", "217", "288", "289",
+                        "1", "1", "1", "1", "0", "0", "1", "1", "0", "0", "1", "1", "0", "0", "0", "0", "0", "0",
+                        "4081", "4098", "5", "10", "", ""),
+                       ("", "", "ATCAGCCATGGCTATTATTG", "ab2.ba1", "153", "164", "217", "260", "289",
+                        "1", "1", "1", "1", "0", "0", "1", "1", "0", "0", "1", "1", "0", "0", "0", "0", "0", "0",
+                        "4081", "4098", "5", "10", "", "")],
+                      [("Chr5:5-20000-11068-C-G", "2.3", "ATCAATATGGCCTCGCCACG", "ab1.ba2", None, None, None, None,
+                        "289", "0", "5", "0", "5", "0", "0", "0", "5", None, None, None, "1", "0",
+                        "0", "0", "0", "0", "0", "4081", "4098", "5", "10", "", ""),
+                       ("", "", "ATCAATATGGCCTCGCCACG", "ab2.ba1", "202", "255", "277", "290", "289",
+                        "1", "3", "1", "3", "0", "0", "1", "3", "0", "0", "1", "1", "0", "0", "0", "0",
+                        "0", "0", "4081", "4098", "5", "10", "", "")],
+                      [("Chr5:5-20000-11068-C-G", "2.4", "ATCAGCCATGGCTATTTTTT", "ab1.ba2", "72", "72", "217", "288", "289",
+                        "1", "1", "1", "1", "0", "0", "1", "1", "0", "0", "1", "1", "0", "0", "0", "0", "0", "0", "4081",
+                        "4098", "5", "10", "", ""),
+                       ("", "", "ATCAGCCATGGCTATTTTTT", "ab2.ba1", "153", "164", "217", "260", "289",
+                        "1", "1", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "1", "1", "0", "0", "0", "0", "4081",
+                        "4098", "5", "10", "", "")],
+                      [("Chr5:5-20000-10776-G-T", "3.1", "ATGCCTACCTCATTTGTCGT", "ab1.ba2", "46", "15", "287", "288", "290",
+                        "3", "3", "3", "2", "3", "1", "0", "1", "1", "0.5", "0", "0.5", "0", "0", "0", "1",
+                        "0", "0", "3", "3", "47170", "41149", "", ""),
+                       ("", "", "ATGCCTACCTCATTTGTCGT", "ab2.ba1", None, "274", None,
+                        "288", "290", "0", "3", "0", "2", "0", "1", "0", "1", None, "0.5", None, "0.5",
+                        "0", "0", "0", "1", "0", "0", "3", "3", "47170", "41149", "", "")],
+                      [("Chr5:5-20000-11315-C-T", "3.2", "ACAACATCACGTATTCAGGT", "ab1.ba2", "197", "197", "240", "255", "271",
+                        "2", "3", "2", "3", "0", "1", "2", "2", "0", "0.333333333333333", "1",
+                        "0.666666666666667", "0", "0", "0", "0", "0", "0", "1", "1", "6584", "6482", "", ""),
+                       ("", "", "ACAACATCACGTATTCAGGT", "ab2.ba1", "35", "35", "240", "258", "271",
+                        "2", "3", "2", "3", "0", "1", "2", "2", "0", "0.333333333333333", "1",
+                        "0.666666666666667", "0", "0", "0", "0", "0", "0", "1", "1", "6584", "6482", "", "")],
+                      [("Chr5:5-20000-13983-G-C", "4.1", "AAAAAAAGAATAACCCACAC", "ab1.ba2", "0", "100", "255", "276", "269",
+                        "5", "6", "0", "6", "0", "0", "5", "6", "0", "0", "0", "1", "0", "0", "0", "0", "5", "0", "1", "1", "5348", "5350", "", ""),
+                       ("", "", "AAAAAAAGAATAACCCACAC", "ab2.ba1", None, None, None, None,
+                        "269", "0", "0", "0", "0", "0", "0", "0", "0", None, None, None, None, "0",
+                        "0", "0", "0", "0", "0", "1", "1", "5348", "5350", "", "")],
+                       [("Chr5:5-20000-13963-T-C", "4.2", "TTTTTAAGAATAACCCACAC", "ab1.ba2", "38", "38", "240", "283", "263",
+                        "110", "54", "110", "54", "0", "0", "110", "54", "0", "0", "1", "1", "0", "0", "0",
+                        "0", "0", "0", "1", "1", "5348", "5350", "", ""),
+                       ("", "", "TTTTTAAGAATAACCCACAC", "ab2.ba1", "100", "112", "140", "145", "263",
+                        "7", "12", "7", "12", "7", "12", "0", "0", "1", "1", "0",
+                        "0", "0", "0", "0", "0", "0", "0", "1", "1", "5348", "5350", "", "")]]
+
+    if chimera_correction: 
+        examples_tiers.extend([[("Chr5:5-20000-13963-T-C", "4.3", "TTTTTAAGAAGCTATTTTTT", "ab1.ba2", "72", "72", "217", "288", "289",
+                        "1", "1", "1", "1", "0", "0", "1", "1", "0", "0", "1", "1", "0", "0", "0", "0", "0", "0", "4081",
+                        "4098", "5", "10", "", "TTTTTAAGAATAACCCACAC"),
+                       ("", "", "TTTTTAAGAAGCTATTTTTT", "ab2.ba1", "153", "164", "217", "260", "289",
+                        "1", "1", "0", "0", "0", "0", "0", "0", "0", "0", "0", "0", "1", "1", "0", "0", "0", "0", "4081",
+                        "4098", "5", "10", "", "TTTTTAAGAATAACCCACAC")],
+                      [("Chr5:5-20000-13983-G-C", "5", "ATGTTGTGAATAACCCACAC", "ab1.ba2", None, "186", None, "276", "269",
+                        "0", "6", "0", "6", "0", "0", "0", "6", "0", "0", "0", "1", "0", "0", "0", "0", "0",
+                        "0", "1", "1", "5348", "5350", "", ""),
+                       ("", "", "ATGTTGTGAATAACCCACAC", "ab2.ba1", None, None, None, None,
+                        "269", "0", "0", "0", "0", "0", "0", "0", "0", None, None, None, None, "0",
+                        "0", "0", "0", "0", "0", "1", "1", "5348", "5350", "", "")]])
+    else:
+        examples_tiers.extend([
+                      [("Chr5:5-20000-13983-G-C", "5", "ATGTTGTGAATAACCCACAC", "ab1.ba2", None, "186", None, "276", "269",
+                        "0", "6", "0", "6", "0", "0", "0", "6", "0", "0", "0", "1", "0", "0", "0", "0", "0",
+                        "0", "1", "1", "5348", "5350", "", ""),
+                       ("", "", "ATGTTGTGAATAACCCACAC", "ab2.ba1", None, None, None, None,
+                        "269", "0", "0", "0", "0", "0", "0", "0", "0", None, None, None, None, "0",
+                        "0", "0", "0", "0", "0", "1", "1", "5348", "5350", "", "")]])
+
+    start_row = 15
+    ws3.write(start_row, 0, "Description of tiers with examples")
+    ws3.write_row(start_row + 1, 0, header_line)
+    row = 0
+    for i in range(len(description_tiers)):
+        ws3.write_row(start_row + 2 + row + i + 1, 0, description_tiers[i])
+        ex = examples_tiers[i]
+        for k in range(len(ex)):
+            ws3.write_row(start_row + 2 + row + i + k + 2, 0, ex[k])
+        ws3.conditional_format('L{}:M{}'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3), {'type': 'formula', 'criteria': '=OR($B${}="1.1", $B${}="1.2")'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 2), 'format': format1, 'multi_range': 'L{}:M{} T{}:U{} B{}'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3)})
+        ws3.conditional_format('L{}:M{}'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3),
+                               {'type': 'formula', 'criteria': '=OR($B${}="2.1",$B${}="2.2", $B${}="2.3", $B${}="2.4")'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 2),
+                                'format': format3,
+                                'multi_range': 'L{}:M{} T{}:U{} B{}'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3)})
+        ws3.conditional_format('L{}:M{}'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3),
+                               {'type': 'formula',
+                                'criteria': '=$B${}>="3"'.format(start_row + 2 + row + i + k + 2),
+                                'format': format2,
+                                'multi_range': 'L{}:M{} T{}:U{} B{}'.format(start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3, start_row + 2 + row + i + k + 2, start_row + 2 + row + i + k + 3)})
+        row += 3
+    workbook.close()
+
+
+if __name__ == '__main__':
+    sys.exit(read2mut(sys.argv))
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/read2mut.xml	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,84 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<tool id="read2mut" name="Call specific mutations in reads:" version="1.0.2" profile="19.01">
+    <description>Looks for reads with mutation at known positions and calculates frequencies and stats.</description>
+    <macros>
+        <import>va_macros.xml</import>
+    </macros>
+    <expand macro="requirements">
+        <requirement type="package" version="1.1.0">xlsxwriter</requirement>
+    </expand>
+    <command><![CDATA[
+        ln -s '$file2' bam_input.bam &&
+        ln -s '${file2.metadata.bam_index}' bam_input.bam.bai &&
+        python '$__tool_directory__/read2mut.py' 
+        --mutFile '$file1'
+        --bamFile bam_input.bam
+        --inputJson '$file3'
+        --sscsJson '$file4'
+        --thresh '$thresh'
+        --phred '$phred'
+        --trim '$trim'
+        $chimera_correction
+        --outputFile '$output_xlsx'
+    ]]>
+    </command>
+    <inputs>
+        <param name="file1" type="data" format="tabular" label="DCS Mutation File" optional="false" help="TABULAR file with DCS mutations. See Help section below for a detailed explanation."/>
+        <param name="file2" type="data" format="bam" label="BAM File of raw reads" optional="false" help="BAM file with aligned raw reads of selected tags."/>
+        <param name="file3" type="data" format="json" label="JSON File with DCS tag stats" optional="false" help="JSON file generated by DCS mutations to tags/reads"/>
+        <param name="file4" type="data" format="json" label="JSON File with SSCS tag stats" optional="false" help="JSON file generated by DCS mutations to SSCS stats."/>
+        <param name="thresh" type="integer" label="Tag count threshold" value="0" help="Integer threshold for displaying mutations. Only mutations occuring in DCS of less than thresh tags are displayed. Default of 0 displays all."/>
+        <param name="phred" type="integer" label="Phred quality score threshold" min="0" max="41" value="20" help="Integer threshold for Phred quality score. Only reads higher than this threshold are considered. Default = 20."/>
+        <param name="trim" type="integer" label="Trimming threshold" value="10" help="Integer threshold for assigning mutations at start and end of reads to lower tier. Default 10."/>
+        <param name="chimera_correction" type="boolean" label="Apply chimera correction?" truevalue="--chimera_correction" falsevalue="" checked="False" help="Add additional tier for chimeric variants and correct the variant frequencies."/>
+    </inputs>
+    <outputs>
+        <data name="output_xlsx" format="xlsx" label="${tool.name} on ${on_string}: XLSX"/>
+    </outputs>
+    <tests>
+        <test>
+            <param name="file1" value="DCS_Mutations_test_data_VA.tabular"/>
+            <param name="file2" value="Interesting_Reads_test_data_VA.trim.bam"/>
+            <param name="file3" value="tag_count_dict_test_data_VA.json"/>
+            <param name="file4" value="SSCS_counts_test_data_VA.json"/>
+            <param name="thresh" value="0"/>
+            <param name="phred" value="20"/>
+            <param name="trim" value="10"/>
+            <param name="chimera_correction" value="False"/>
+            <output name="output_xlsx" file="mutant_reads_summary_short_trim_test_data_VA.xlsx" decompress="true" lines_diff="10"/>
+        </test>
+    </tests>
+    <help> <![CDATA[
+**What it does**
+
+Takes a tabular file with mutations, a BAM file of aligned raw reads, and JSON files 
+created by the tools **DCS mutations to tags/reads** and **DCS mutations to SSCS stats** 
+as input and calculates frequencies and stats for DCS mutations based on information 
+from the raw reads.
+
+**Input** 
+
+**Dataset 1:** Tabular file with duplex consesus sequence (DCS) mutations as 
+generated by the **Variant Annotator** tool.
+
+**Dataset 2:** BAM file of aligned raw reads. This file can be obtained by the 
+tool `Map with BWA-MEM <https://arxiv.org/abs/1303.3997>`_.
+
+**Dataset 3:** JSON file generated by the **DCS mutations to tags/reads** tool 
+containing dictonaries of the tags of reads containing mutations 
+in the DCS.
+
+**Dataset 4:** JSON file generated by the **DCS mutations to SSCS stats** tool 
+stats of tags that carry a mutation in the SSCS at the same position a mutation 
+is called in the DCS.
+
+**Output**
+
+The output is an XLSX file containing frequencies stats for DCS mutations based 
+on information from the raw reads. In addition to that a tier based 
+classification is provided based on the amout of support for a true variant call.
+
+    ]]> 
+    </help>
+    <expand macro="citation" />
+</tool>
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/Aligned_Families_test_data_VA.tabular	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,41 @@
+GATAACCTTGCTTCGTGATTAATC	ab	1	M01897:257:000000000-AYB6W:1:2112:28792:17250 2:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAATGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	GGGGGGGGGGGGGGFGFFGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGAFGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGFGG8<FED@FGEGGGGGGGGGGGGGGGGGGGGGGGGFEFGGGGGGGEEGGGGGGGGGCECE8EGGGEFFGGGGGGFGGGDEDD5EGGGGGFGGGDCFCFGFGFGCFAGGFDFFFEEGGFF3><:>>FD>FFC=4=:0<;DD>6461992<)892<AFBFFFFFF244:-1:1>:=0306(4)-.42((44(667(449?0,
+GATAACCTTGCTTCGTGATTAATC	ba	2	M01897:257:000000000-AYB6W:1:1108:16316:3620 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGGGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGG	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGFGGGGGGGGGGFGGGGGGGGGGGFGGGGGGFCFGGGGFGGGGEGGFGGGGGGGFFGGGGGGGGGGGGGGGGBB?CFGGGGGGCCGGFGGGGGGGCFGECC79CEECGDDD99CFGGGGF4>>GG>BDE@BBFG5=B?AF98::A42<EF?;B::((7:?7???<)/:?91;1,6?F?29902(
+GATAACCTTGCTTCGTGATTAATC	ba	2	M01897:257:000000000-AYB6W:1:1118:22651:3876 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGG	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGDGGGGEEGGGGGGGGGGGGGGGGGGGCEFDFGGGFGFGGEGEBFGG<FGFCFEFDFDF<DGGGGGGFGFGFGGGGGFGGGCE:FFBFGGGGDGGGGGGGGDD@FGGFC6AE7E1CGGFGCCFGGGGCEGFGGGCCFG9A*59@FGD><?9=CFF6>3BBDFFF392?G)-96<2<:<:44<232:B3:F>6??F0(34A248:>?1,(.404-,((4(
+GATAACCTTGCTTCGTGATTAATC	ba	2	M01897:257:000000000-AYB6W:1:1118:5518:20674 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACCCCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCCAGAAGCGGGACGGCCGTAAGTCCCAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	FGEF9FGGCCG@FGGGGGGCFCC@EGGG<FFGGFGGCGGGGGEG?FGGGGDGGGGGGGGFEFGGGGGGGGG7==F:F,=FEGGFFFFG<<F<=8FG>CAF9E8CFCEFFFFF,?F=F8FFD=,DFFE+@CGGGFF7D<FGGEFEGGFG2DCCFGECC*=CGEGGGGGGGGCCFFFGF7FFGFGGGGGGGFFG=E56C55CEGF3:F*9*./>FG***27)?::D)5557@>BFD@)/))).(().9<2((-29BF>F4(83,:12-)4)2,3??<<1:(7>((,
+GATAACCTTGCTTCGTGATTAATC	ab	2	M01897:257:000000000-AYB6W:1:2112:28792:17250 1:N:0:1	GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTGCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACG	FGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGFGGGGGFFGGGGGGGGGGGGGGGGCFGGFFGGGCGGGGGGGGGCEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDFGGGGGGGGGGGGGGGGGGGGGGGGEGFGGFGFFGFGFGFFBEFGFFFF7F?FFB?DF>FDFFB:)9>FBFFF?F099E>;<?1:<?0>F0;BB1
+GATAACCTTGCTTCGTGATTAATC	ba	1	M01897:257:000000000-AYB6W:1:1108:16316:3620 2:N:0:1	GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTGCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCGCGGGACACG	GGGGGGGGGGGGGGGGGFGGGGGGGGGGGFGGGGGGGGGGGGGGDGGGGGGGGGGGFGGGGGGGGFGGGGGGGGGGGGDGCEEFGGGGGCGGGGGGGGGGGGGFGGGGGGGGFGGGGGGGGGGFGGGCGFFGGGGGGGGGGFFFFGGGGCEG==CECFFCDGGGGGGGGGGGGG597*<FGGFGDFC35>+*:=6FDFF4CFFF9B204>G?FE)5FAF?:7>FBB<A?FB(9?<AFFF<0?B?F4:BFF2>B69>;B))6<<?(,(46((4,42(7>926(82
+GATAACCTTGCTTCGTGATTAATC	ba	1	M01897:257:000000000-AYB6W:1:1118:22651:3876 2:N:0:1	GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTTCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCCCGGGACACG	FDCCF9FFDFGGGGGGGGGGGGFGGGGGGFGGGGGGGGGGGEGFECFGGGGGGGGGGGGG8EGFGEGFG,F@GGGGGEGGEG7FGGGGG@BFFGCCGGEGGGGGGGFGGGGGGFGGG@FFF9CFGGGGGGGGGGGGGGGGGGGFFGF5E*CECC>EGFGG7EGD==?E8:E7CCE3C+?:C?FFG@D3B5:>78)/C6=FFF<>B>>0:@EBF3))14>B?20>?A<2:>99>F<<AD7??BF0??)8<0<BF?>>>FAA9A:,403>BF?2;B(46(4(((((
+GATTGGATAACGTTGTGGCAATTG	ab	1	M01897:257:000000000-AYB6W:1:1101:14310:2734 2:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGGAAGTCACCGGAATCCGGGACGTCCTGGCAGCTAGGGCGGGCCCCGAGCCAGG	GGGGCCFGGGGGGGGGFFGGGGGGGFGGGGEGGGGGGGFGGGGGGGGFGGGGGGGGGGGGGGGGGGGG<@FGGCCGGGGGGCFFGGGGGGGGC,@ECFBFDDGGG@FGGGGG9CFG@CCFF@DCDFC>=CEGGDEGCC@CDC*=CC*=5>FGCFEGGDFGG?<EGGGFFFD49FD=6>:CGFD>5)/)47C@4),85:B:DF?(8)448:D:,5?7**430;>01661(-4((74,94:)(,-(-18(--(-(-(-(=01,442))(.(,(2((-(-(-((,((
+GATTGGATAACGTTGTGGCAATTG	ab	1	M01897:257:000000000-AYB6W:1:1112:4840:17845 2:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCCGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCATACTGGGCACAGGGCCAGGCGTGAGGGCTCAAGAAGCGGGACCGCCGTCAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTCGGCGTGTCCCTAGCTAGC	GGGC6<,CFDG8FG,CF6C<FGGDD<FGGGFGGGFG<E8,6FCFC,77BF,CFGGGFFFGGFEF:<CECC7:F:@?DFCEGFFE8?EGGF<,+,44=7,B,C@DGC@7+?F8,7=D,A=>,9=FFFG:@=BC7CCEFGFDGGGG788CEF66EFGGG7CF*:**=C5=FGG5AC=+:C*2:EFF*7*2/97DD>FC)7>@G@5(704(255005FFFB??FFB39((,--32()(./6>B<(())9))-38>0,43(-((((33<)-,((--(.4)).43)).(
+GATTGGATAACGTTGTGGCAATTG	ab	1	M01897:257:000000000-AYB6W:1:1118:8154:20084 2:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	GGGGGGGGGGGGGGFGGGGGGGGGGGGFGGGG@EDFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG=FGGGGGFFBFGGGGGGGFGGFGGGGFGGGFGGGGGEGGGGGGG:FGGGG5EGGGG:FEFGEGGGGGGGGGGGGGGGD:EG?FGFFGFCEGG>GFFFCGGFFFGEFE:>?(7.()44>B>G*=F<7:F9>D>9>F03;26:6)6>B<9(38<7A?FB2>>?FF(=:?(((.2:A:)-4((.63-,49>:?0
+GATTGGATAACGTTGTGGCAATTG	ab	1	M01897:257:000000000-AYB6W:1:2110:10849:23965 2:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG@CGGGGGGGFGGGGECGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGEGGGCCCGGGGGGGGGGGGGGGGGDGGEGGGDGGGGDFGGGGGFDGDDFFFGGGFFFGFGFE@:?GFFFFFGFFFFFD2?BFFFF09>B9>F(7)2.9A2)6:44<@A7BF?>BF?>6>:((,(,5AF?F91(-:B<>,(3>00(
+GATTGGATAACGTTGTGGCAATTG	ba	2	M01897:257:000000000-AYB6W:1:1113:13084:11145 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGTACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGGCFGGGGFFFFBGFFFBEFFFFGFF?@AFFB?FFFFFFFFFFFFFFFB00:?FFFAFFFFFFFFF66>FFF<1
+GATTGGATAACGTTGTGGCAATTG	ba	2	M01897:257:000000000-AYB6W:1:1115:14952:6061 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	GFGGGGGEGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGGGGGGGGGGGGGGFGGGGGGGGGGGGG,FGFFFFFGAFGGGGGGGGGGGGFGFGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGEGFGGGGGEGFGGGDCBEGGGGGGGGGGGGGGGGFEGGFGFGGGGGGFDGGGCDGD9DFFFGD4>FGG4FF@9DD>DD>>FFDBGFFDFFFFFFFFF=?:8?F><F>?F?FBFF?7>F>DB<>?B9>>?9>9?F1
+GATTGGATAACGTTGTGGCAATTG	ba	2	M01897:257:000000000-AYB6W:1:2101:12835:23979 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGACGTAAGTCCAAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGGGTCCCGAGCCAGC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG8FFGGGGFFGGGGGGGGGGFGGGGGGDFGGGGGFGGGGGGGGGGGGGEGG7FFGGGG@FGGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGFGGFGGGGC6CFFFGGCFFGGGGGG:C:47FFD3*1<677<6<;EGB@><)-3:>55-9))).:12<6)4430;>3>0(*4??F1(.:7?>(,((-.8B1999B?1
+GATTGGATAACGTTGTGGCAATTG	ba	2	M01897:257:000000000-AYB6W:1:2102:25716:13556 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGCTTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGCGCCAGC	,C@FCGGGGGE8C@FCFFGEECGCFGGDFFEFFGFFGGFGGGFFGFCF@FFAFGGGGGCGEG,EFGGG?=F@EGGGC,<=FF?DEGFFFGFFGFFDGCGG?FDDG>EFGFGA9?EFE@FGGGDFGFCFFGC+@EEE@:F:E7C1:FBCF@7<2CFFF**::8CFEEGE7C8BFF?CFGGC<9CFGCEG+CCC8:CFDCDC=:**202:65*CF5CGFD)6?5))).753>><:5>9@466-.((.9::0)4B>)8><>:0(80:2))501(--3:FF(,4>02(
+GATTGGATAACGTTGTGGCAATTG	ba	2	M01897:257:000000000-AYB6W:1:2107:10919:21008 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACCGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGCEGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG>EGEGGGFGGGGGGGGGGGGGGGGGFFGGGGGGGGGGGGGGCGGGGGGGGFGGFGGGFG@EGGBDGGFFFFGFFFFF)*4<:B@F?G6<>9BFFFFB?BFF?DF?BFF00:BFFFB;BFBB2
+GATTGGATAACGTTGTGGCAATTG	ab	2	M01897:257:000000000-AYB6W:1:1101:14310:2734 1:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCGCTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGCGTCTCCCGGGTGCCTGGTGGCGGTGTGGGGG	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGFFEGGGGGFFGGGGGGGGGGGGGGFFEE:=FED=FFFFFFEEEGGGGGG:FC:FFGGGGCFGGGGGGGFGFGGGCCG9FGGGGGGGGGFFEFBFGGEEE<EGGGGGGGGGGGGGF@FCCCEEGGCFFGGFGGGGGG:>EEGDGFDC4>EDDFGGEBEFGE5>CGGFFF*)<FF<<FF:61:<BFFB7??9::?FF:07<7)(.,,2<1(11(,3:>7:773(-766223:((
+GATTGGATAACGTTGTGGCAATTG	ab	2	M01897:257:000000000-AYB6W:1:1112:4840:17845 1:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCTTGGCCCGGTCCTGGCCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCCGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGCGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC	FCEFGG@@FGGGGEFGGFFGGDCF8CG86ECEGGGFGF,C,C@@FFFGGE:,CFFDG7CFFGECE=<FF9<F,C<7++=7+4@+@=FFFGGFA,CB,EF9@7F::3@F@BFCC7E=FC@FCCCFF<=FGA7:FGFG,37F9FCG7:3?>7:FCGGG:@FC6B,=EE7FFE>EGG9C?*=5CC7887*/=*:?C5E76C*::*//C>D*8D4377CF*;:?)055.547;FF?4*7F)<2@<AF:))766:23)(731F>?>(41((4))8>7:0(--,-338((
+GATTGGATAACGTTGTGGCAATTG	ab	2	M01897:257:000000000-AYB6W:1:1118:8154:20084 1:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC	GGGGGGGGGGGGGGGGGGGDDGGFFGGDFGGGGGGGGGGGGFEECFGGCFGEGGGGGGGGGGGGGGGFGGFCCFFFFGGGDGGGDGGGFGGGGGCFFFGGFGGGGGGGGFGC7EGG<=BFFGEG<DCFGGFGGGGGGFGGGGEGGGFGGGGGCGCGGGFCCFGGGGGGGGEGGGGGGGGGGGGGGFEEGEGGG5>FGGGGGDGGGDGFFFGFGGGFGGFFF6@FFFFFFFF<FFFFFF???FFA>?B2>B<0<?78AFF1706>9B?AF:?:0(3139:FF<?0
+GATTGGATAACGTTGTGGCAATTG	ab	2	M01897:257:000000000-AYB6W:1:2110:10849:23965 1:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGEFGFGGGGGGGGGGGGGGGEGGGGGGGFGGFDGGGGGGGGGGDGGGGGGGGGGGGFGFF=EF9E7E?FFGFDFC?>GGFDBGFGFC?FFGABFFFB>G?AFFBFB@>>?BB>BEF?AB:??0:B;71??F?(.
+GATTGGATAACGTTGTGGCAATTG	ba	1	M01897:257:000000000-AYB6W:1:1113:13084:11145 2:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTACCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGFGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGFCDFBGGGFFFFFFFGF:?D>FFBF>?F@FFFFFBFFFF??FFFF?62>:?FF>?FDFFFF?0:?F?ABFF><AFB>09BF?AFFF?:?>>9>?FF::0
+GATTGGATAACGTTGTGGCAATTG	ba	1	M01897:257:000000000-AYB6W:1:1115:14952:6061 2:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC	GFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDGGGCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGEDGF@FGGGGGGGGFGGGGGGCGGGFGGGGFGFFFCD@7DGF58FFFFFGF3:>D:6>>GBBFF474<?FFFF?B?B0(:1:F?068>:79?28508?>>4<04>AA<09>0>F?:6<B?F0(4969<F:?:(
+GATTGGATAACGTTGTGGCAATTG	ba	1	M01897:257:000000000-AYB6W:1:2101:12835:23979 2:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCCCTTCCCATCTGGGTCCCCAACGGCCTCTCCTGCGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGTC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGEEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGD7FEFBFGGGGGGGCFGGGFGGGGGCGEGGG=CFGGECBFCFACC7DDFGGC9FFGCCFFFFF<@F?FFEGGGGGG8:57ACC=@60;CC:7,CEEGCGG4CFCA<?C<<FGGAFCDC5:>6C?3C.76)0319:*4)57?F*5<?2=FFF:?328()7395?(,22)((79>?:7:93:B)1)21>99(489<1(((,751681(8-(
+GATTGGATAACGTTGTGGCAATTG	ba	1	M01897:257:000000000-AYB6W:1:2102:25716:13556 2:N:0:1	GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTTGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCCTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGCGCCCACTTCCCATCTGCGTCCCCACAGGCCTCTCCTGTTGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGTCC	;,EACFFEGD6FFGGG8F<@<FEGGGGGGFCE7F@,EDG9<@:FCFCEF,CFGF>BF6+,9BFA+4EDFFFEFGCC,BFCB:FFGD==C?CFG,,CEE9E7CGGE@FCCGD+8:CC<9DFGG,@:B9:F9BC@5DD5FFG;@DFBCECEE7EC7,,?CFDCDGGGFFFFC://=?=4CF4C+C558DDF5EDGBB5/*.<DEE:))4()-*-)*948*6:74).4;9?F:1((>0?@B2?04)4)1699<<:>247667<340;>B?A71((-49@(-((4())
+CCTCCCGGCAGTGCGAAAATGTCA	ab	1	M01897:257:000000000-AYB6W:1:1108:17396:7377 2:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGGCCAGGCCGGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGTGGCTGCCCAGGCGGCCTGTTTTTTTGCAGGCTCCCTACGCTACGGGGTGGGCTTTTTCCGTTTCATCTTGGTGTTGCCGGCTGGGACGCCTTGCGCC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG7FGGGGGGGGGGGC<FGGGDGGGGGGGGGGGGGGGGGGGGDGGGGGGGGGDEC8C>C5*:/C:*:<+2/>C:*:*+*<>?+*+0<5:/>E5<35***<6293*935=DC)))1707C5)(1*))())()*06)(((0,(*(,(,(-4(9),4D6(4((5)4*(,).2))-).5)5:228))-1(-(((((-((,()5(-(
+CCTCCCGGCAGTGCGAAAATGTCA	ba	2	M01897:257:000000000-AYB6W:1:1106:22053:22582 1:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGCCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGACGCGGGCAGTGTGTATGCAGTCATCCTCAGCTACGGGCTGGGCTTCTTCCTGTTTATCCTGGTGGTGGCGGCTGTGTCGCTCTGCCGTC	CFGGAFCFCFGGGFDGDDDGGDGGGG;F:BFGEGFGGGGFF<FFDECG@CFDGGF@FECFAEGFGGGGGAFFEGGGEGF<?E@FFGFEFGEGG+BEF=<FGGCFCFGGGGGGGG8FDFGGDF@FFGGGEEG*88:C88AFEC>8A:@;EFG8>:EEGE0<CCF+<E:CE/C8C*8C*;;C:0*;=EFEDG*/0*7*:7*18*27:CFGD?>>7+CGG>?F:?4*7?FG6).-))7)/<BF0)6.)/--/)67.:F209304(((493(,:5-)(2;:<2).4((
+CCTCCCGGCAGTGCGAAAATGTCA	ba	2	M01897:257:000000000-AYB6W:1:1109:21342:10199 1:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGCCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGAGGCGGGCAGTGTGTATGCAGGCATCCTCAGCTACGGGGTGGGCTTCTTCCTGTTCATCCTGGTGGTGGCGGCTGTGACGCTCTGCCGCC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDGC?>DFGGFGGGGGGFFFFFFFFFF@FFFFCDFGF?FFAFFFDAAFFBFB9?FFD08<<6?BFFF;F?2<??6??<7>B>9
+CCTCCCGGCAGTGCGAAAATGTCA	ba	2	M01897:257:000000000-AYB6W:1:1111:28216:18792 1:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGCCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGAGGCGGGCAGTGTGTATGCAGGCATCCTCAGCTACGGGGTGGGCTTCTTCCTGTTCATCCTGGTGGTGGCGGCTGTGACGCTCTGCCGCC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG=8FGGEGDGGGGGGGGCFFGGGGGGGGGFFFFFFFFFGFFFFB5<BEFB>8AABAFF<9<5FBF?):F:B?:2@FFFF1.54<?:.323<?FF9
+CCTCCCGGCAGTGCGAAAATGTCA	ab	2	M01897:257:000000000-AYB6W:1:1108:17396:7377 1:N:0:1	CTAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGTGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCGCGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGAAGCCCACCCCGT	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGFFGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGBCGGGGGGGGGGGGGGCFGGGGGGGGEGGEGGGGGGGEGGGGGGGGGGGGGGGDGDDDEFDGGFFGFFFFFGFFFF>EFBFFFGFFFFF:BFFF?F?FFFFFF?F<BBF??BBFFFFBBFF
+CCTCCCGGCAGTGCGAAAATGTCA	ba	1	M01897:257:000000000-AYB6W:1:1106:22053:22582 2:N:0:1	CTAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCACCAGGAGGCCTGGCGGGCCGGCAGCTCAGAACCTGATATCTACTTTCTGTTAGCTGTCGCTCGAGCGGGAAGCGGGAGATCTTGTGCGCGGTGTGGGAGCCTAGCCCTTTCTTGGGGTGGCTGCGCAGGCGGCAGAGCGTCACAGCTGCTACAACCAGGATGAACAGGAAGAGCCCCACCCCGTC	FCF<9C@F8E9@FGC,,,,<,CF<,C@B@CC@<F,,@F::FD+FC@@F,CFFEEDFGD:C=<<B?FF:E8,B,B,AC<FA8C44++B=>F7F?+A7FF+==<F+:+@7+AFB,8C:F**>CC@F?CCFFCFC@C,26,3224@C@C,,?CG+<+2CFC*:*:);C7E*21*9CE**>DDFC7+:0=/))5C)1)(*)00>*9:(.4(,577:*=47)721),,),(-(4(47()((43460(.)(0..).))).4(()(,(,)6)((((,4((((4(-((((((
+CCTCCCGGCAGTGCGAAAATGTCA	ba	1	M01897:257:000000000-AYB6W:1:1109:21342:10199 2:N:0:1	CTAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCGCGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGCAGCCCAGCCCGT-	GGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGFGGGGGGGGGFGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGGGGGGGGGGGGGFEGDFGFGGGGFGFGGFGGGEG?FGCDGGEGGGGGGGGG6>FEGFDFGGFFGGGEE3DFF@=@FFGF2?>FB9FFFFFBFFFBFFFFFF9>>F>F68?>>?:BABFFFFF6B??:?BF5<>BB<49?:?:?(4?:0:0(.3399 
+CCTAGTCTTTGATTGGCCACTTTT	ab	1	M01897:257:000000000-AYB6W:1:1106:12553:14962 2:N:0:1	TAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCACGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGAAGCCCACCCCGT--	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDFFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGEEEGGGGGGGGFGGGGGGFGGFFGGD@9DFFFGGFFFFF7AFFGFFFFFFFFFFFFFFBFF6>FFFFBBB>FBBADBFFF((428?F<?F>:?:DFA:1:?FF7??F10(3<F96  
+CCTAGTCTTTGATTGGCCACTTTT	ab	1	M01897:257:000000000-AYB6W:1:1106:15615:18803 2:N:0:1	TAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCTCCAGGCGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCACGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGAAGCCCACCCCGTAC	GGCGGGGGGGGFFFGGCFGFGFGGGD@@F7FGGDCC:EFCFFGG:@DGGGGFGGGGGGGEGGG9EFGGGGGGGGGGG,CEEGDGEGGGGGGFFGGGGCGGGFG<FFDGGFEG@EG7FDEGCEEDECGGCFGGGGGFCFF2CCEFCFGCF;FGFGFF5BCFFCFGGGGGG5CGC=EE<CGGGGG6CBFFEFF53@6CED>755>:DF>AFFGA=>FF>>09>09@6>>BFBF(3.116>F:0)-:1)37(:<?BFFF<)57?AA2:637?0,(4.10(:39(31(
+CCTAGTCTTTGATTGGCCACTTTT	ab	1	M01897:257:000000000-AYB6W:1:1110:11692:17499 2:N:0:1	TAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCACGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGACGCCCACCCCGT--	GGGGG@EGGGGGGGGGFGGGGFGGGGCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDEGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGFFGGDFFFGFFFGFFF<FFBFFGFFFFFDFFFFFFFFF6((49>>)4:D>F6:AFBBFFFFFF0(8,1<6?7(6<B:22::?(4999<F?BBF?  
+CCTAGTCTTTGATTGGCCACTTTT	ab	1	M01897:257:000000000-AYB6W:1:1110:21292:16434 2:N:0:1	TAGGCTCTACATGGTGAGCAGAGACGCGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGTGGGAAGCGGGAGATCTTGTGCACGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGTTGAACAGGAAGAAGCCCACCCCTT--	GFGGGGD8C@FGGFCFCFFGGGGGEE+@CCFFEDCDFCF@CFCFGGGGGGGFEDFFGGGDGF@FFGEGF<FGGCF7<+4,=BFCFGGGDGG:?<=EE:FGFD8F+:C+>FGACCCCFGGE7EGGGGGGGFGCCFF;@A<;E;<BFCCGG;>CFGFFFFG5:*=+588C8>57EEGGDEFGGGCFC?*9@FGGGB>)97)<=@F?/))3:BGFFFE:)5925?>326909>>6>54,7((-8)-8-71(--24641:B)47445270,(3124(.,(,:<>(.  
+CCTAGTCTTTGATTGGCCACTTTT	ab	2	M01897:257:000000000-AYB6W:1:1106:12553:14962 1:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGG--CCAGGCAGGGCCCCAAGCCCCTTGTCT-TGCAGCCGGGGGGGGGGCGGTGGGAGCCTAACAAGCGGGGCGGGGGGTTGGAGGCCTCCCCAAGTTCGGGGGTGGCTTCTTCCTGTTCATCCTTGGTGTGGGGGCTGTGACGCCTTTGCGGC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG*  ?8*;*;**:*;***2***2A***00+< C++0++;***:**:*****://:**;**0++*++2*/:E/*1**)))/)1)+*1**))9))**)/**)03>))8D)(8().5<*)-7))1)67)6/.8118((4(-,.()-(()(-)).(,-
+CCTAGTCTTTGATTGGCCACTTTT	ab	2	M01897:257:000000000-AYB6W:1:1106:15615:18803 1:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGA--CCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGAGGCGGGCAGTGTGTATGCAGGCATCCTCAGCTACGGGGTGGGCTTCTTCCTGTTCATCCTGGTGGTGGCGGCTGTGACGCTCTGCC-GC	GGGGGFGGGGFGGGGGGGGGFGGGGGGGGGGEGGGGGGGGGDEFGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGCFGGGGGGGCFGEGFGEEGGGGGGGDGGDGGCFGGGGGGGDGGGGGGFGGGGGGGGGGG  GDGGGGGGGFGGGGFGGGGGGGGGFGG9FFGGGGGGGGGCEGECGGCFCEEGGGGGGGGGGGCGGCEC3*:C>DG=FC<?CGGFFGFFFGFFGFFFFFG:BDFFF<:AFFFFFFFF4<4?BAFB:<BB??09>B0?B><:D243847 10
+CCTAGTCTTTGATTGGCCACTTTT	ab	2	M01897:257:000000000-AYB6W:1:1110:11692:17499 1:N:0:1	CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGGGGCCAGGCCAGGCCCCAACGCCCATGTCTTTGCAGCCGAGGGGGAGCTGGTTGGGGCTGACGAGGCGGGCAGTGGGTATGCAGGCATCCTCAGCTACGGGGTGGGCTTCTTCCTGTTCATCCTGGGGGTGGCGGCTTGTACCCTCTTCC---	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG8>5*;:8C8E@;?*:;88*2CE*8*++<C99+9+@C**88*:C*:?6<C*+1*:858C*;7/9E*1CGCCC0*)1)*+<C7C.5766<69=<))9*05>3/4;<31<2)9:4=).0))/69?<((213(7:960(,1.-))))(()-)   
+GATAAGCCAACTGCCATCTAGAAT	ab	1	M01897:257:000000000-AYB6W:1:1105:25798:19415 2:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTACGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTTTCCCGAGCCAGC	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGEGGGGGFGGGFGGGGGGGGGGGGGGEFCFFGGDEGG8EGGGGFGGGGEGFEGFFG<EDGGGGGGGGGCC7FGGGGGGGFGG=>FFCDBF)7:>7FF:EF?<?FEE:@F@?6??F6>B01>BF;FFF*4(,2:24?FBBF>?F?FFBF0;B2:0(:??FF7:BF?03:2<BBFBFB?0
+GATAAGCCAACTGCCATCTAGAAT	ba	2	M01897:257:000000000-AYB6W:1:2104:15100:19675 1:N:0:1	CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTCGCGTGTCCCGAGCCAGC	GGFGFGGDGGGGGGCGGGD@EFGGFFGEGGGGGGGFECDEDFFGGGFGGGGFGDEGDGGEGGGGGG:FGGG@CFECFGGGGGFFGGGGGGF<CGGEGGGEFGGGEE7FGFF,=B=DBBFFDDFFFGGGGGGGEGGG:><FEGCF:FGEFGGFFFGGGDGGFEGDFGGGFGGGFFGGGCGGGEGFGGGGFFCFFEDG57CGGCFFC6C*CEGG6:CGGG:6<C>>CEFDGB7B5/<<:9<>>><F279?FG<>>>:>:D(47:6<26)402346>2<>(-49??0
+GATAAGCCAACTGCCATCTAGAAT	ab	2	M01897:257:000000000-AYB6W:1:1105:25798:19415 1:N:0:1	GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTGCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACG	GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFFFFGGFFGDB9FFGFFFF0F?:?FFFFFFFFFFFFBF@FFBA?B9;9B9>BB>FFF>FF>><?4
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/DCS_Mutations_test_data_VA.tabular	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,5 @@
+#SAMPLE	CHR	POS	A	C	G	T	CVRG	ALLELES	MAJOR	MINOR	MAF	BIAS
+__NONE__	ACH_TDII_5regions	505	1	2208	0	0	2209	1	C	A	0.00045	1.09465
+__NONE__	ACH_TDII_5regions	571	0	2817	0	1	2818	1	C	T	0.00035	1.04139
+__NONE__	ACH_TDII_5regions	958	0	1	0	14667	14668	1	T	C	7e-05	1.03624
+
Binary file test-data/DCS_test_data_VA.bam has changed
Binary file test-data/DCS_test_data_VA.bam.bai has changed
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/Interesting_Reads_test_data_VA.fastq	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,124 @@
+@GATAACCTTGCTTCGTGATTAATC.ab.1
+CTAGAGGGCCAGACCCTGGAGAGAATGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+GGGGGGGGGGGGGGFGFFGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGAFGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGFGG8<FED@FGEGGGGGGGGGGGGGGGGGGGGGGGGFEFGGGGGGGEEGGGGGGGGGCECE8EGGGEFFGGGGGGFGGGDEDD5EGGGGGFGGGDCFCFGFGFGCFAGGFDFFFEEGGFF3><:>>FD>FFC=4=:0<;DD>6461992<)892<AFBFFFFFF244:-1:1>:=0306(4)-.42((44(667(449?0,
+@GATAACCTTGCTTCGTGATTAATC.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGGGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGG
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGFGGGGGGGGGGFGGGGGGGGGGGFGGGGGGFCFGGGGFGGGGEGGFGGGGGGGFFGGGGGGGGGGGGGGGGBB?CFGGGGGGCCGGFGGGGGGGCFGECC79CEECGDDD99CFGGGGF4>>GG>BDE@BBFG5=B?AF98::A42<EF?;B::((7:?7???<)/:?91;1,6?F?29902(
+@GATAACCTTGCTTCGTGATTAATC.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGG
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGDGGGGEEGGGGGGGGGGGGGGGGGGGCEFDFGGGFGFGGEGEBFGG<FGFCFEFDFDF<DGGGGGGFGFGFGGGGGFGGGCE:FFBFGGGGDGGGGGGGGDD@FGGFC6AE7E1CGGFGCCFGGGGCEGFGGGCCFG9A*59@FGD><?9=CFF6>3BBDFFF392?G)-96<2<:<:44<232:B3:F>6??F0(34A248:>?1,(.404-,((4(
+@GATAACCTTGCTTCGTGATTAATC.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCGCCCGGGAGACCCCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGTCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCCAGAAGCGGGACGGCCGTAAGTCCCAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+FGEF9FGGCCG@FGGGGGGCFCC@EGGG<FFGGFGGCGGGGGEG?FGGGGDGGGGGGGGFEFGGGGGGGGG7==F:F,=FEGGFFFFG<<F<=8FG>CAF9E8CFCEFFFFF,?F=F8FFD=,DFFE+@CGGGFF7D<FGGEFEGGFG2DCCFGECC*=CGEGGGGGGGGCCFFFGF7FFGFGGGGGGGFFG=E56C55CEGF3:F*9*./>FG***27)?::D)5557@>BFD@)/))).(().9<2((-29BF>F4(83,:12-)4)2,3??<<1:(7>((,
+@GATAACCTTGCTTCGTGATTAATC.ab.2
+GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTGCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACG
++
+FGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGFGGGGGFFGGGGGGGGGGGGGGGGCFGGFFGGGCGGGGGGGGGCEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDFGGGGGGGGGGGGGGGGGGGGGGGGEGFGGFGFFGFGFGFFBEFGFFFF7F?FFB?DF>FDFFB:)9>FBFFF?F099E>;<?1:<?0>F0;BB1
+@GATAACCTTGCTTCGTGATTAATC.ba.1
+GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTGCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCGCGGGACACG
++
+GGGGGGGGGGGGGGGGGFGGGGGGGGGGGFGGGGGGGGGGGGGGDGGGGGGGGGGGFGGGGGGGGFGGGGGGGGGGGGDGCEEFGGGGGCGGGGGGGGGGGGGFGGGGGGGGFGGGGGGGGGGFGGGCGFFGGGGGGGGGGFFFFGGGGCEG==CECFFCDGGGGGGGGGGGGG597*<FGGFGDFC35>+*:=6FDFF4CFFF9B204>G?FE)5FAF?:7>FBB<A?FB(9?<AFFF<0?B?F4:BFF2>B69>;B))6<<?(,(46((4,42(7>926(82
+@GATAACCTTGCTTCGTGATTAATC.ba.1
+GATCCTGCCGTGTGGACTCTGTGCGGTGCCCGCAGGGCGGTGCTGGCGCTCGCCTATCGCTCTGCTCTCTCTTTGTAGACGGCGGGCGCTAACACCACCGACAAGGAGCTAGAGGTTCTCTCCTTGCACAACGTCACCTTTGAGGACGCCGGGGAGTACACCTGCCTGGCGGGCAATTCTATTGGGTTTTCTCATCACTCTGCGTGGCTGGTGGTTCTGCCAGGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCCCGGGACACG
++
+FDCCF9FFDFGGGGGGGGGGGGFGGGGGGFGGGGGGGGGGGEGFECFGGGGGGGGGGGGG8EGFGEGFG,F@GGGGGEGGEG7FGGGGG@BFFGCCGGEGGGGGGGFGGGGGGFGGG@FFF9CFGGGGGGGGGGGGGGGGGGGFFGF5E*CECC>EGFGG7EGD==?E8:E7CCE3C+?:C?FFG@D3B5:>78)/C6=FFF<>B>>0:@EBF3))14>B?20>?A<2:>99>F<<AD7??BF0??)8<0<BF?>>>FAA9A:,403>BF?2;B(46(4(((((
+@GATTGGATAACGTTGTGGCAATTG.ab.1
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGGAAGTCACCGGAATCCGGGACGTCCTGGCAGCTAGGGCGGGCCCCGAGCCAGG
++
+GGGGCCFGGGGGGGGGFFGGGGGGGFGGGGEGGGGGGGFGGGGGGGGFGGGGGGGGGGGGGGGGGGGG<@FGGCCGGGGGGCFFGGGGGGGGC,@ECFBFDDGGG@FGGGGG9CFG@CCFF@DCDFC>=CEGGDEGCC@CDC*=CC*=5>FGCFEGGDFGG?<EGGGFFFD49FD=6>:CGFD>5)/)47C@4),85:B:DF?(8)448:D:,5?7**430;>01661(-4((74,94:)(,-(-18(--(-(-(-(=01,442))(.(,(2((-(-(-((,((
+@GATTGGATAACGTTGTGGCAATTG.ab.1
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCCGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCATACTGGGCACAGGGCCAGGCGTGAGGGCTCAAGAAGCGGGACCGCCGTCAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTCGGCGTGTCCCTAGCTAGC
++
+GGGC6<,CFDG8FG,CF6C<FGGDD<FGGGFGGGFG<E8,6FCFC,77BF,CFGGGFFFGGFEF:<CECC7:F:@?DFCEGFFE8?EGGF<,+,44=7,B,C@DGC@7+?F8,7=D,A=>,9=FFFG:@=BC7CCEFGFDGGGG788CEF66EFGGG7CF*:**=C5=FGG5AC=+:C*2:EFF*7*2/97DD>FC)7>@G@5(704(255005FFFB??FFB39((,--32()(./6>B<(())9))-38>0,43(-((((33<)-,((--(.4)).43)).(
+@GATTGGATAACGTTGTGGCAATTG.ab.1
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+GGGGGGGGGGGGGGFGGGGGGGGGGGGFGGGG@EDFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG=FGGGGGFFBFGGGGGGGFGGFGGGGFGGGFGGGGGEGGGGGGG:FGGGG5EGGGG:FEFGEGGGGGGGGGGGGGGGD:EG?FGFFGFCEGG>GFFFCGGFFFGEFE:>?(7.()44>B>G*=F<7:F9>D>9>F03;26:6)6>B<9(38<7A?FB2>>?FF(=:?(((.2:A:)-4((.63-,49>:?0
+@GATTGGATAACGTTGTGGCAATTG.ab.1
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG@CGGGGGGGFGGGGECGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGEGGGCCCGGGGGGGGGGGGGGGGGDGGEGGGDGGGGDFGGGGGFDGDDFFFGGGFFFGFGFE@:?GFFFFFGFFFFFD2?BFFFF09>B9>F(7)2.9A2)6:44<@A7BF?>BF?>6>:((,(,5AF?F91(-:B<>,(3>00(
+@GATTGGATAACGTTGTGGCAATTG.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGTACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGGCFGGGGFFFFBGFFFBEFFFFGFF?@AFFB?FFFFFFFFFFFFFFFB00:?FFFAFFFFFFFFF66>FFF<1
+@GATTGGATAACGTTGTGGCAATTG.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+GFGGGGGEGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGGGGGGGGGGGGGGFGGGGGGGGGGGGG,FGFFFFFGAFGGGGGGGGGGGGFGFGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGEGFGGGGGEGFGGGDCBEGGGGGGGGGGGGGGGGFEGGFGFGGGGGGFDGGGCDGD9DFFFGD4>FGG4FF@9DD>DD>>FFDBGFFDFFFFFFFFF=?:8?F><F>?F?FBFF?7>F>DB<>?B9>>?9>9?F1
+@GATTGGATAACGTTGTGGCAATTG.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGACGTAAGTCCAAGGATTCCCGTCCGTCCTGGCAGCTTTGGCGGGTCCCGAGCCAGC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG8FFGGGGFFGGGGGGGGGGFGGGGGGDFGGGGGFGGGGGGGGGGGGGEGG7FFGGGG@FGGGGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGFGGFGGGGC6CFFFGGCFFGGGGGG:C:47FFD3*1<677<6<;EGB@><)-3:>55-9))).:12<6)4430;>3>0(*4??F1(.:7?>(,((-.8B1999B?1
+@GATTGGATAACGTTGTGGCAATTG.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACAGGCTTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGCGCCAGC
++
+,C@FCGGGGGE8C@FCFFGEECGCFGGDFFEFFGFFGGFGGGFFGFCF@FFAFGGGGGCGEG,EFGGG?=F@EGGGC,<=FF?DEGFFFGFFGFFDGCGG?FDDG>EFGFGA9?EFE@FGGGDFGFCFFGC+@EEE@:F:E7C1:FBCF@7<2CFFF**::8CFEEGE7C8BFF?CFGGC<9CFGCEG+CCC8:CFDCDC=:**202:65*CF5CGFD)6?5))).753>><:5>9@466-.((.9::0)4B>)8><>:0(80:2))501(--3:FF(,4>02(
+@GATTGGATAACGTTGTGGCAATTG.ba.2
+CTAGAGGGCCAGACCCTGGAGAGAAGGAGCCCAGCAGAGCCAGCCAGTCCCACACCGCCACCAGGCACCCGGGAGACACCAGAGCCACAGGAGAGGCCTTTGGGGACCCAGATGGGAAGTGGGCTCGAGGGGGCTGAGGGGGCCCCTCTGGGACCAGGACCGGGCCAGGCCAACTTTGTCCCCACACTGGGCACAGGGCCAGGAGTGAGGGCTCAAGAAGCGGGACGGCCGTAAGTCACCGGATTCCCGTCCGTCCTGGCAGCTTTGGCGTGTCCCGAGCCAGC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGCEGGGGGGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG>EGEGGGFGGGGGGGGGGGGGGGGGFFGGGGGGGGGGGGGGCGGGGGGGGFGGFGGGFG@EGGBDGGFFFFGFFFFF)*4<:B@F?G6<>9BFFFFB?BFF?DF?BFF00:BFFFB;BFBB2
+@GATTGGATAACGTTGTGGCAATTG.ab.2
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCGCTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGCGTCTCCCGGGTGCCTGGTGGCGGTGTGGGGG
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGFFEGGGGGFFGGGGGGGGGGGGGGFFEE:=FED=FFFFFFEEEGGGGGG:FC:FFGGGGCFGGGGGGGFGFGGGCCG9FGGGGGGGGGFFEFBFGGEEE<EGGGGGGGGGGGGGF@FCCCEEGGCFFGGFGGGGGG:>EEGDGFDC4>EDDFGGEBEFGE5>CGGFFF*)<FF<<FF:61:<BFFB7??9::?FF:07<7)(.,,2<1(11(,3:>7:773(-766223:((
+@GATTGGATAACGTTGTGGCAATTG.ab.2
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCTTGGCCCGGTCCTGGCCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCCGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGCGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC
++
+FCEFGG@@FGGGGEFGGFFGGDCF8CG86ECEGGGFGF,C,C@@FFFGGE:,CFFDG7CFFGECE=<FF9<F,C<7++=7+4@+@=FFFGGFA,CB,EF9@7F::3@F@BFCC7E=FC@FCCCFF<=FGA7:FGFG,37F9FCG7:3?>7:FCGGG:@FC6B,=EE7FFE>EGG9C?*=5CC7887*/=*:?C5E76C*::*//C>D*8D4377CF*;:?)055.547;FF?4*7F)<2@<AF:))766:23)(731F>?>(41((4))8>7:0(--,-338((
+@GATTGGATAACGTTGTGGCAATTG.ab.2
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC
++
+GGGGGGGGGGGGGGGGGGGDDGGFFGGDFGGGGGGGGGGGGFEECFGGCFGEGGGGGGGGGGGGGGGFGGFCCFFFFGGGDGGGDGGGFGGGGGCFFFGGFGGGGGGGGFGC7EGG<=BFFGEG<DCFGGFGGGGGGFGGGGEGGGFGGGGGCGCGGGFCCFGGGGGGGGEGGGGGGGGGGGGGGFEEGEGGG5>FGGGGGDGGGDGFFFGFGGGFGGFFF6@FFFFFFFF<FFFFFF???FFA>?B2>B<0<?78AFF1706>9B?AF:?:0(3139:FF<?0
+@GATTGGATAACGTTGTGGCAATTG.ab.2
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGEFGFGGGGGGGGGGGGGGGEGGGGGGGFGGFDGGGGGGGGGGDGGGGGGGGGGGGFGFF=EF9E7E?FFGFDFC?>GGFDBGFGFC?FFGABFFFB>G?AFFBFB@>>?BB>BEF?AB:??0:B;71??F?(.
+@GATTGGATAACGTTGTGGCAATTG.ba.1
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTACCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGFGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGFCDFBGGGFFFFFFFGF:?D>FFBF>?F@FFFFFBFFFF??FFFF?62>:?FF>?FDFFFF?0:?F?ABFF><AFB>09BF?AFFF?:?>>9>?FF::0
+@GATTGGATAACGTTGTGGCAATTG.ba.1
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCACTTCCCATCTGGGTCCCCAAAGGCCTCTCCTGTGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGAC
++
+GFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDGGGCGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGEDGF@FGGGGGGGGFGGGGGGCGGGFGGGGFGFFFCD@7DGF58FFFFFGF3:>D:6>>GBBFF474<?FFFF?B?B0(:1:F?068>:79?28508?>>4<04>AA<09>0>F?:6<B?F0(4969<F:?:(
+@GATTGGATAACGTTGTGGCAATTG.ba.1
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTGGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCTTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGAGCCCCCTTCCCATCTGGGTCCCCAACGGCCTCTCCTGCGGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGGTC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGEEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGD7FEFBFGGGGGGGCFGGGFGGGGGCGEGGG=CFGGECBFCFACC7DDFGGC9FFGCCFFFFF<@F?FFEGGGGGG8:57ACC=@60;CC:7,CEEGCGG4CFCA<?C<<FGGAFCDC5:>6C?3C.76)0319:*4)57?F*5<?2=FFF:?328()7395?(,22)((79>?:7:93:B)1)21>99(489<1(((,751681(8-(
+@GATTGGATAACGTTGTGGCAATTG.ba.1
+GGTACCGGCTTCTGCTGCTGCTGCTGCTCCGCACTGTCTGGGGGACGCTTGCTCGGGACACGCCAAAGCTGCCAGGACGGACGGGAATCCTGTGACTTACGGCCGTCCCGCCTCTTGAGCCCTCACTCCTGGCCCTGTGCCCAGTGTGGGGACAAAGTTGGCCTGGCCCGGTCCTGGTCCCAGAGGGGCCCCCTCAGCCCCCTCGCGCCCACTTCCCATCTGCGTCCCCACAGGCCTCTCCTGTTGCTCTGGTGTCTCCCGGGTGCCTGGTGGCGGTGTGGTCC
++
+;,EACFFEGD6FFGGG8F<@<FEGGGGGGFCE7F@,EDG9<@:FCFCEF,CFGF>BF6+,9BFA+4EDFFFEFGCC,BFCB:FFGD==C?CFG,,CEE9E7CGGE@FCCGD+8:CC<9DFGG,@:B9:F9BC@5DD5FFG;@DFBCECEE7EC7,,?CFDCDGGGFFFFC://=?=4CF4C+C558DDF5EDGBB5/*.<DEE:))4()-*-)*948*6:74).4;9?F:1((>0?@B2?04)4)1699<<:>247667<340;>B?A71((-49@(-((4())
+@CCTCCCGGCAGTGCGAAAATGTCA.ab.1
+CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGGCCAGGCCGGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGTGGCTGCCCAGGCGGCCTGTTTTTTTGCAGGCTCCCTACGCTACGGGGTGGGCTTTTTCCGTTTCATCTTGGTGTTGCCGGCTGGGACGCCTTGCGCC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG7FGGGGGGGGGGGC<FGGGDGGGGGGGGGGGGGGGGGGGGDGGGGGGGGGDEC8C>C5*:/C:*:<+2/>C:*:*+*<>?+*+0<5:/>E5<35***<6293*935=DC)))1707C5)(1*))())()*06)(((0,(*(,(,(-4(9),4D6(4((5)4*(,).2))-).5)5:228))-1(-(((((-((,()5(-(
+@CCTCCCGGCAGTGCGAAAATGTCA.ba.2
+CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGCCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGACGCGGGCAGTGTGTATGCAGTCATCCTCAGCTACGGGCTGGGCTTCTTCCTGTTTATCCTGGTGGTGGCGGCTGTGTCGCTCTGCCGTC
++
+CFGGAFCFCFGGGFDGDDDGGDGGGG;F:BFGEGFGGGGFF<FFDECG@CFDGGF@FECFAEGFGGGGGAFFEGGGEGF<?E@FFGFEFGEGG+BEF=<FGGCFCFGGGGGGGG8FDFGGDF@FFGGGEEG*88:C88AFEC>8A:@;EFG8>:EEGE0<CCF+<E:CE/C8C*8C*;;C:0*;=EFEDG*/0*7*:7*18*27:CFGD?>>7+CGG>?F:?4*7?FG6).-))7)/<BF0)6.)/--/)67.:F209304(((493(,:5-)(2;:<2).4((
+@CCTCCCGGCAGTGCGAAAATGTCA.ba.2
+CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGCCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGAGGCGGGCAGTGTGTATGCAGGCATCCTCAGCTACGGGGTGGGCTTCTTCCTGTTCATCCTGGTGGTGGCGGCTGTGACGCTCTGCCGCC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGDGC?>DFGGFGGGGGGFFFFFFFFFF@FFFFCDFGF?FFAFFFDAAFFBFB9?FFD08<<6?BFFF;F?2<??6??<7>B>9
+@CCTCCCGGCAGTGCGAAAATGTCA.ba.2
+CTAGACTCACTGGCGTTACTGACTGCGAGACCCTCCAGACAAGGCGCGTGCTGAGGTTCTGAGCCCCCTTCCGCTCCCAGTGGTGCCTGCGGCTCTGGGCCAGGGGCATCCATGGGAGCCCCGTGGGGGGGGGGCCAGGCCAGGCCTCAACGCCCATGTCTTTGCAGCCGAGGAGGAGCTGGTGGAGGCTGACGAGGCGGGCAGTGTGTATGCAGGCATCCTCAGCTACGGGGTGGGCTTCTTCCTGTTCATCCTGGTGGTGGCGGCTGTGACGCTCTGCCGCC
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG=8FGGEGDGGGGGGGGCFFGGGGGGGGGFFFFFFFFFGFFFFB5<BEFB>8AABAFF<9<5FBF?):F:B?:2@FFFF1.54<?:.323<?FF9
+@CCTCCCGGCAGTGCGAAAATGTCA.ab.2
+CTAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGTGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCGCGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGAAGCCCACCCCGT
++
+GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGFFGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGBCGGGGGGGGGGGGGGCFGGGGGGGGEGGEGGGGGGGEGGGGGGGGGGGGGGGDGDDDEFDGGFFGFFFFFGFFFF>EFBFFFGFFFFF:BFFF?F?FFFFFF?F<BBF??BBFFFFBBFF
+@CCTCCCGGCAGTGCGAAAATGTCA.ba.1
+CTAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCACCAGGAGGCCTGGCGGGCCGGCAGCTCAGAACCTGATATCTACTTTCTGTTAGCTGTCGCTCGAGCGGGAAGCGGGAGATCTTGTGCGCGGTGTGGGAGCCTAGCCCTTTCTTGGGGTGGCTGCGCAGGCGGCAGAGCGTCACAGCTGCTACAACCAGGATGAACAGGAAGAGCCCCACCCCGTC
++
+FCF<9C@F8E9@FGC,,,,<,CF<,C@B@CC@<F,,@F::FD+FC@@F,CFFEEDFGD:C=<<B?FF:E8,B,B,AC<FA8C44++B=>F7F?+A7FF+==<F+:+@7+AFB,8C:F**>CC@F?CCFFCFC@C,26,3224@C@C,,?CG+<+2CFC*:*:);C7E*21*9CE**>DDFC7+:0=/))5C)1)(*)00>*9:(.4(,577:*=47)721),,),(-(4(47()((43460(.)(0..).))).4(()(,(,)6)((((,4((((4(-((((((
+@CCTCCCGGCAGTGCGAAAATGTCA.ba.1
+CTAGGCTCTACATGGTGAGCAGAGACGAGGAGAGGGGAGCCCGCCTGGCTGCAGAGAGGGCTCACACAGCCCAGGACCAGCGTGGGCCGAGGTGGGGCTCCAGGAGGCCTGGCGGGCAGGCAGCTCAGAACCTGGTATCTACTTTCTGTTACCTGTCGCTTGAGCGGGAAGCGGGAGATCTTGTGCGCGGTGGGGGAGCCCAGGCCTTTCTTGGGGGGGCTGCGCAGGCGGCAGAGCGTCACAGCCGCCACCACCAGGATGAACAGGAAGCAGCCCAGCCCGT
++
+GGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGFGGGGGGGGGFGEGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFFGGGGGGGGGGGGGFEGDFGFGGGGFGFGGFGGGEG?FGCDGGEGGGGGGGGG6>FEGFDFGGFFGGGEE3DFF@=@FFGF2?>FB9FFFFFBFFFBFFFFFF9>>F>F68?>>?:BABFFFFF6B??:?BF5<>BB<49?:?:?(4?:0:0(.3399
Binary file test-data/Interesting_Reads_test_data_VA.trim.bam has changed
Binary file test-data/Interesting_Reads_test_data_VA.trim.bam.bai has changed
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/SSCS_counts_test_data_VA.json	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,1 @@
+[{"ACH_TDII_5regions#505": {"ab": 2, "ba": 1}, "ACH_TDII_5regions#571": {"ab": 1, "ba": 1}, "ACH_TDII_5regions#958": {"ab": 1, "ba": 1}}, {"ACH_TDII_5regions#505": {"ab": 1, "ba": 1}, "ACH_TDII_5regions#571": {"ab": 2, "ba": 1}, "ACH_TDII_5regions#958": {"ab": 1}}]
\ No newline at end of file
Binary file test-data/SSCS_test_data_VA.bam has changed
Binary file test-data/SSCS_test_data_VA.bam.bai has changed
Binary file test-data/mutant_reads_summary_short_trim_test_data_VA.xlsx has changed
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/test-data/tag_count_dict_test_data_VA.json	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,1 @@
+[{"GATAACCTTGCTTCGTGATTAATC": {"ACH_TDII_5regions#505": "A"}, "GATTGGATAACGTTGTGGCAATTG": {"ACH_TDII_5regions#571": "T"}, "CCTCCCGGCAGTGCGAAAATGTCA": {"ACH_TDII_5regions#958": "C"}}, {"ACH_TDII_5regions#505": [1, 1, 173.0], "ACH_TDII_5regions#571": [1, 1, 143.0], "ACH_TDII_5regions#958": [0, 1, 195.0]}]
\ No newline at end of file
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/va_macros.xml	Sun Oct 04 17:19:39 2020 +0000
@@ -0,0 +1,20 @@
+<macros>
+    <xml name="citation">
+        <citations>
+            <citation type="bibtex">
+@misc{duplex,
+    author = {Povysil, Gundula and Heinzl, Monika and Salazar, Renato and Stoler, Nicholas and Nekrutenko, Anton and Tiemann-Boege, Irene},
+    year = {2019},
+    title = {{Variant Analyzer: a quality control for variant calling in duplex sequencing data (manuscript)}}
+ }
+           </citation>
+        </citations>
+    </xml>
+    <xml name="requirements">
+        <requirements>
+            <requirement type="package" version="3.1.2">matplotlib</requirement>        
+            <requirement type="package" version="0.15">pysam</requirement>
+            <yield/>
+        </requirements>
+    </xml>
+</macros>