shell:
    "abricate {params} {input} > {output.csv}"
    " && python3 {scripts_directory}get_abricate_info_list.py {output.csv} {output.yaml}"

shell:
    "abricate {params} {input} > {output.csv}"
    " && python3 {scripts_directory}get_abricate_info_list.py {output.csv} {output.yaml}"

shell:
    "fastqc {input} -o {output.dir} {params} && python3 {scripts_directory}move_fastqc_output.py {wildcards.sample}"

shell:
    "fastqc {input.read_1} {input.read_2} -o {output.dir} {params} -t {threads}"

shell:
    "fastqc {input.read} -o {output.dir} {params} -t {threads}"

shell:
    "fasttree {params} {input} > {output.tree} && cp {output.tree} {output.nwk}"

shell:
    "fasttree {params} {input} > {output}"

shell:
    "cp {input.tree} {output.tree} &&"
    "cp {input.tree} {output.nwk}"

shell:
  "python {coverage_script} -i {input.depth} -r {REFERENCE_FASTA}"
   " -o {output.coverage} -t {params.threshold} -g {REFERENCE_GB}"

shell:
    "python {scripts_directory}merge_coverage.py '{input}' {output.coverage_regular} {output.coverage_translate} {REFERENCE_GB}"

shell:
    "mkdir align_samples/{wildcards.sample}/reference -p && "
    "cp {REFERENCE_FASTA} align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta && "
    "bwa index align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta && "
    "bwa mem align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta {input.reads_1} | "
    "samtools view -u -T align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta -q {params.minqual} | "
    "samtools sort --reference align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta > {output}"

shell:
    "samtools mpileup -A -d 600000 -B -Q 0 {input.bam} --reference align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta | "
    "ivar variants -p align_samples/{wildcards.sample}/iVar/snps -q {params.mapqual} -t {params.minfrac} -m {params.mincov}  align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta "

shell:
    "samtools mpileup -d 0 -A -Q {params.mapqual} {input.bam} --reference align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta | "
    "ivar consensus -p align_samples/{wildcards.sample}/iVar/{params.consensus} -q {params.mapqual} -t {params.minfrac} -n N -m {params.mincov} "

shell:
    "python ../workflow/scripts/split_consensus.py {input.consensus} {input.depth} {REFERENCE_GB} {output.consensus}"

shell:
    "python ../workflow/scripts/split_depth_file.py {input.zipped} {REFERENCE_GB}"

shell:
    "python ../workflow/scripts/mask_consensus_by_deep.py align_samples/{wildcards.sample}/reference/{REFERENCE_NAME}.fasta {input.first_consensus} {output.align_file} {wildcards.seg}"

shell:
    "python ../workflow/scripts/msa_masker.py -i {input.align_file} -df {input.depth} -o {output} {params}"

shell:
    "python ../workflow/scripts/get_consensus_medaka.py '{input}' {output}"

shell:
    "python {scripts_directory}mask_regions.py {input.consensus} {output.final_consensus} {params} "

shell:
    "python {scripts_directory}convert_vcf.py {input.tsv_file} {output.vcf_file} "

shell:
    "mafft --thread {threads} {params} {input.align_file} > {output.aligned_file}"

shell:
    "mafft --thread {threads} {params} {input.align_file} > {output.aligned_file}"

shell:
    "mafft --thread {threads} {params} {input} > {output}"

shell:
    "mafft --thread {threads} {params} {input} > {output}"

shell:
    "mafft --thread {threads} {params} {input} > {output}"

shell:
    "mafft --thread {threads} {params} {input} > {output}"

shell:
    "cp {input} {output}"

shell:
    "mkdir {output} && cp {software_parameters_path} projects/{wildcards.project}/"

shell:
    "mkdir projects/{wildcards.project}/sample_{wildcards.sample}/ -p && "
    " cp -r {params} projects/{wildcards.project}/sample_{wildcards.sample}/ "

shell:
    "python {scripts_directory}generate_AllConsensus.py {input.coverage} {REFERENCE_GB} '{input.every_consensus}' {REFERENCE_FASTA} {output.AllConsensus} {output.all_consensus_no_ref} "
    "&& python {scripts_directory}concat_segments.py '{input.every_consensus}' {REFERENCE_GB} {output.All_nt} {input.coverage} {REFERENCE_FASTA} {output.All_nt_only_90plus}"

shell:
    "python {scripts_directory}split_files_by_locus.py {input} projects/{PROJECT_NAME}/main_result {REFERENCE_GB}"

shell:
    "samtools depth {params} {input.i} | bgzip -c > {output.only_depth} "
    "&& tabix -p vcf {output.only_depth} && gunzip -c {output.only_depth}  > {output.depth} "

shell:
    "python3 {scripts_directory}split_depth_file.py {input} {REFERENCE_GB}"

shell:
    "medaka variant --verbose {input.ref} {input.hdf} {output.vcf}"

shell:
    "medaka tools annotate {input.vcf} {input.ref} {input.bam} {output.snps} && "
    "bgzip {output.snps} -c > {output.vcf_zipped} && tabix {output.vcf_zipped}"

shell:
    "touch {output.temp_file} && "
    "python {scripts_directory}add_freq_medaka_consensus.py {input.normal_reference_fasta} {input.vcf_file} {input.file_coverage} "
    "{output.vcf_file_out} {output.vcf_file_out_removed_by_filter} {output.temp_file} {output.normal_reference_fasta_removed} {output.temp_file} {params} &&"
    "bgzip -c -f {output.vcf_file_out} > {output.vcf_file_out_compr} && tabix {output.vcf_file_out_compr}"

shell:
    "python {scripts_directory}mask_consensus_by_deep.py {REFERENCE_FASTA} {input.first_consensus} {output.align_file} {wildcards.seg}"

shell:
    "python {scripts_directory}msa_masker.py -i {input.align_file} -df {input.depth} -o {output} {params}"

shell:
    "python {scripts_directory}get_consensus_medaka.py '{input}' {output}"

shell:
    "python {scripts_directory}mask_regions.py {input} {output} {params}"

shell:
    "NanoStat --fastq {input} --outdir {output.dir} -n {wildcards.sample}_stats.txt -t {threads}"

shell:
    "abricate {params} {REFERENCE_FASTA} > {output.csv}"
    " && python3 {scripts_directory}get_abricate_info.py {output.csv} {output.yaml}"

shell:
    "echo 'This is not SARS-CoV-2' > {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "seqret {params} -sequence {input} -outseq {output}"

shell:
    "rm -r align_samples/{wildcards.sample}/snippy/ && "
    "snippy --cpus {threads} --se {input.read} --ref {REFERENCE_FASTA} --outdir align_samples/{wildcards.sample}/snippy/ {params}"

shell:
    "gunzip -c {input.zipped} > {output.unzipped}"

shell:
    "python3 {scripts_directory}split_depth_file.py align_samples/{wildcards.sample}/snippy/depth/snps.depth {REFERENCE_GB} "
    "&& touch {output.unzipped}"

shell:
    "python {scripts_directory}mask_consensus_by_deep.py {REFERENCE_FASTA} {input.first_consensus} {output.align_file} {wildcards.seg}"

shell:
    "python {scripts_directory}msa_masker.py -i {input.align_file} -df {input.depth} -o {output} {params}"

shell:
    "python {scripts_directory}get_consensus_medaka.py '{input}' {output}"

shell:
    "python {scripts_directory}mask_regions.py {input} {output} {params}"

shell:
    "python {scripts_directory}create_snpeff_text.py $CONDA_PREFIX {input.ref_gb} {input.ref_fa} {REFERENCE_NAME} {output} "

shell:
    "{replace} {input.i1} && "
    "snpEff {params} -v {REFERENCE_NAME} {input.i1}  > {output}"

shell:
    "spades.py -t {threads}  {params} -1 {input.read_1} -2 {input.read_2} -o {output.dir}"

shell:
    "python {scripts_directory}translate.py {REFERENCE_GB} {input.Alignment_nt} {output} '{wildcards.locus}' '{wildcards.gene}' {input.coverage}"

shell:
    "trimmomatic SE "
    "-threads {threads} "
    "{input} "
    "{output} "
    "{params}"

shell:
    "python {scripts_directory}variants.py '{input}' '{output}' validated_variants"

shell:
    "python {scripts_directory}variants.py '{input}' {output} minor_iSNVs"

shell:
    "python {scripts_directory}variants.py '{input}' {output} minor_iSNVs_inc_indels"

shell:
    "python {scripts_directory}proportions_iSNVs_graph.py '{input}' {output.out_file}"

shell:
    "echo 'No sample with all segments above 90% coverage for the value given in the user/parameters.yaml file.' > {output}"

from Bio import SeqIO


def is_fasta(filename: str) -> bool:
    """
    The is_fasta function checks if a file is in FASTA format.

    :param filename: str: Specify the file to be checked
    :return: True if the file is a fasta file, and false otherwise
    :doc-author: Trelent
    """
    with open(filename, "r", encoding="utf-8") as handle:
        fasta = SeqIO.parse(handle, "fasta")
        return any(fasta)


def is_genbank(filename: str) -> bool:
    """
    The is_gbk function checks if a file is in GenBank format.

    :param filename: str: Specify the file to be checked
    :return: True if the file is a genbank file, and false otherwise
    :doc-author: Trelent
    """
    with open(filename, "r", encoding="utf-8") as handle:
        gbk = SeqIO.parse(handle, "genbank")
        return any(gbk)


def same_identifiers(fasta_file: str, gbk_file: str) -> bool:
    """
    The same_identifiers function checks whether the IDs in a FASTA file match those
    in a GenBank file. It does this by opening both files and parsing them with
    the SeqIO module from Biopython. The function returns True if the lists of IDs
    are identical, and False otherwise.

    :param fasta_file: str: Specify the path to the fasta file
    :param gbk_file: str: Specify the path to the genbank file
    :return: True if the order of the ids in both files are equal
    :doc-author: Trelent
    """
    with open(fasta_file, "r", encoding="utf-8") as handle:
        fasta = SeqIO.parse(handle, "fasta")  # type: ignore
        fasta_ids: list[str] = [record.id for record in fasta]
    with open(gbk_file, "r", encoding="utf-8") as handle:
        gbk = SeqIO.parse(handle, "genbank")
        gbk_locus: list[str] = [record.name for record in gbk]
    return sorted(fasta_ids) == sorted(gbk_locus)

import sys
import csv
import re
from Bio import SeqIO
from extract_gb_info import get_locus
from yaml_io import read_yaml


def get_fasta_reference_concat(fasta_ref):
    """
    The get_fasta_reference_concat function takes a fasta file as inputand returns the concatenated
    reference sequence.
    The function is used to create a single, concatenated reference sequence from multiple contigs
    in the assembly.

    :param fasta_ref: Specify the reference genome to be used for mapping
    :return: The concatenated reference sequence in the form of a seqrecord
    """
    concat_reference = []
    count = 0
    for record in SeqIO.parse(fasta_ref, "fasta"):
        if count == 0:
            concat_reference.append(record)
            start = fasta_ref.index("reference")
            end = fasta_ref.index(".fasta")
            # why + 11 => 16/05
            concat_reference[0].id = fasta_ref[start + 11 : end]
        else:
            concat_reference[0].seq += record.seq
        count += 1
    return concat_reference[0]


def get_id_sequence_from_consensus(consensus, sample_name, value_list):
    """
    The get_id_sequence_from_consensus function takes a fasta file of consensus sequences and
    returns the sequence of the locus that is specified by the value_list. The first record in this
    fasta file will be returned with its id changed to sample_name.
    This function assumes that there are no duplicate records in consensus, which may not be true.

    :param consensus: Define the consensus sequence file
    :param sample_name: Name the sequence in the seqrecord object
    :param value_list: Specify which segments to concatenate
    :return: A seqrecord object that contains the consensus sequence for a given sample
    :doc-author: Trelent
    """

    concat_consensus = []
    count = 1

    starting_segment = value_list[0]
    for record in SeqIO.parse(consensus, "fasta"):
        if count == starting_segment:
            concat_consensus.append(record)
            concat_consensus[0].id = sample_name
        elif count in value_list:
            concat_consensus[0].seq += record.seq
        count += 1
    return concat_consensus[0]


def get_id_sequence_from_consensus_strict(consensus, sample_name, value_list):
    """
    The get_id_sequence_from_consensus_strict function takes a consensus sequence file, the sample name, and a list of integers.
    It then parses through the consensus sequence file to find all sequences that match the starting segment in the list.
    The function concatenates these sequences together into one record and assigns it an ID based on its sample name.  It returns this record.

    :param consensus: Specify the path to a fasta file containing the consensus sequence for each sample
    :param sample_name: Name the sequence
    :param value_list: Specify the segments that will be concatenated
    :return: A seqrecord object from a fasta file
    :doc-author: Trelent
    """
    concat_consensus = []
    count = 1

    starting_segment = value_list[0]
    for record in SeqIO.parse(consensus, "fasta"):
        if count == starting_segment:
            concat_consensus.append(record)
            concat_consensus[0].id = ""
            concat_consensus[0].id = sample_name
            concat_consensus[0].description = ""
        elif count in value_list:
            concat_consensus[0].seq += record.seq
        count += 1
    return concat_consensus[0]


def create_consensus_file_for_alignment(
    consensus_list: list,
    reference_gb: str,
    output: str,
    coverage_file: str,
    reference_fasta: str,
    output_only_90_plus: str,
):
    """
    The create_consensus_file_for_alignment function takes a list of files, the reference genome in genbank format,
    the species name and an output file name as input. It then creates a fasta file with all the sequences from the
    consensus files that have coverage above 90%. The function also takes into account if it is flu or not. If it is flu
    it will only include samples that have 8 loci covered at 90% or more in one file and in the other it will concatonate all
    segments that have 90% coverage or more.

    :param consensus: Specify the path to the directory containing all of your consensus sequences
    :param reference_gb: Get the name of the reference sequence in a genbank file
    :param species: Determine which reference file to use
    :param output: Name the output file
    :param coverage_file: Specify the file containing all the coverage values for each sample
    :param fasta_file: Get the reference sequence from the fasta file
    :param output_only_90_plus: Create a fasta file with only the sequences that have at least 90% coverage
    :return: A fasta file with the consensus sequences of all samples that have a coverage over 90%
    :doc-author: Trelent
    """
    with open(coverage_file, encoding="UTF8", newline="") as csvfile:
        csv_reader = csv.reader(csvfile, delimiter=",")
        coverage_list = list(csv_reader)
    dic_directory = read_yaml("../config/constants.yaml")

    software_parameters = read_yaml(dic_directory["software_parameters"])
    coverage_value = software_parameters["min_coverage_consensus"]
    for row in coverage_list:
        for value in range(len(row) - 1, 0, -1):
            # print(row[value])
            if float(row[value]) > coverage_value:
                row[value] = value
            else:
                row.pop(value)

    coverage_dic = {}
    for i in coverage_list:
        coverage_dic[i[0]] = i[1:]

    fasta_reference_concat = get_fasta_reference_concat(reference_fasta)
    new_concat_file = [fasta_reference_concat]
    new_concat_file_strict = [fasta_reference_concat]
    if len(get_locus(reference_gb)) > 1:
        for file in consensus_list:
            sample_name = re.findall("(?<=sample_)(.*?)(?=/)", file)[0]
            if len(coverage_dic[sample_name]) > 0:
                new_sequence = get_id_sequence_from_consensus(
                    file, sample_name, coverage_dic[sample_name]
                )
                new_concat_file.append(new_sequence)
            if len(coverage_dic[sample_name]) == 8:
                new_sequence = get_id_sequence_from_consensus_strict(
                    file, sample_name, coverage_dic[sample_name]
                )
                new_concat_file_strict.append(new_sequence)
        SeqIO.write(new_concat_file, output_only_90_plus, "fasta")
        SeqIO.write(new_concat_file_strict, output, "fasta")

    else:
        for file in consensus_list:
            sample_name = re.findall("(?<=sample_)(.*?)(?=/)", file)[0]
            if len(coverage_dic[sample_name]) == 1:
                new_sequence = get_id_sequence_from_consensus_strict(
                    file, sample_name, coverage_dic[sample_name]
                )
                new_concat_file_strict.append(new_sequence)
        SeqIO.write(new_concat_file_strict, output_only_90_plus, "fasta")
        SeqIO.write(new_concat_file_strict, output, "fasta")


if __name__ == "__main__":
    CONSENSUS_LIST = sys.argv[1].split(" ")
    REFERENCE_GB = sys.argv[2]
    OUTPUT_FILE = sys.argv[3]
    COVERAGE_FILE = sys.argv[4]
    REFERENCE_FASTA = sys.argv[5]
    OUTPUT_FILE_WITH_90_PLUS_COVERAGE = sys.argv[6]

    create_consensus_file_for_alignment(
        CONSENSUS_LIST,
        REFERENCE_GB,
        OUTPUT_FILE,
        COVERAGE_FILE,
        REFERENCE_FASTA,
        OUTPUT_FILE_WITH_90_PLUS_COVERAGE,
    )

import os
import sys
import re
import errno
import argparse


def parse_args(args=None):
    Description = "Convert iVar variants tsv file to vcf format."
    Epilog = """Example usage: python ivar_variants_to_vcf.py <FILE_IN> <FILE_OUT>"""

    parser = argparse.ArgumentParser(description=Description, epilog=Epilog)
    parser.add_argument("FILE_IN", help="Input tsv file.")
    parser.add_argument("FILE_OUT", help="Full path to output vcf file.")
    parser.add_argument(
        "-po",
        "--pass_only",
        dest="PASS_ONLY",
        help="Only output variants that PASS all filters.",
        action="store_true",
    )
    parser.add_argument(
        "-af",
        "--allele_freq_thresh",
        type=float,
        dest="ALLELE_FREQ_THRESH",
        default=0,
        help="Only output variants where allele frequency greater than this number (default: 0).",
    )

    return parser.parse_args(args)


def make_dir(path):
    if not len(path) == 0:
        try:
            os.makedirs(path)
        except OSError as exception:
            if exception.errno != errno.EEXIST:
                raise


def ivar_variants_to_vcf(FileIn, FileOut, passOnly=False, minAF=0):
    filename = os.path.splitext(FileIn)[0]
    header = (
        "##fileformat=VCFv4.2\n"
        "##source=iVar\n"
        '##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">\n'
        '##FILTER=<ID=PASS,Description="Result of p-value <= 0.05">\n'
        '##FILTER=<ID=FAIL,Description="Result of p-value > 0.05">\n'
        '##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">\n'
        '##FORMAT=<ID=REF_DP,Number=1,Type=Integer,Description="Depth of reference base">\n'
        '##FORMAT=<ID=REF_RV,Number=1,Type=Integer,Description="Depth of reference base on reverse reads">\n'
        '##FORMAT=<ID=REF_QUAL,Number=1,Type=Integer,Description="Mean quality of reference base">\n'
        '##FORMAT=<ID=ALT_DP,Number=1,Type=Integer,Description="Depth of alternate base">\n'
        '##FORMAT=<ID=ALT_RV,Number=1,Type=Integer,Description="Deapth of alternate base on reverse reads">\n'
        '##FORMAT=<ID=ALT_QUAL,Number=1,Type=String,Description="Mean quality of alternate base">\n'
        '##FORMAT=<ID=ALT_FREQ,Number=1,Type=String,Description="Frequency of alternate base">\n'
    )
    header += (
        "#CHROM\tPOS\tID\tREF\tALT\tQUAL\tFILTER\tINFO\tFORMAT\t" + filename + "\n"
    )

    varList = []
    varCountDict = {"SNP": 0, "INS": 0, "DEL": 0}
    OutDir = os.path.dirname(FileOut)
    make_dir(OutDir)
    fout = open(FileOut, "w")
    fout.write(header)
    with open(FileIn) as f:
        for line in f:
            if not re.match("REGION", line):
                line = re.split("\t", line)
                CHROM = line[0]
                POS = line[1]
                ID = "."
                REF = line[2]
                ALT = line[3]
                var_type = "SNP"
                if ALT[0] == "+":
                    ALT = REF + ALT[1:]
                    var_type = "INS"
                elif ALT[0] == "-":
                    REF += ALT[1:]
                    ALT = line[2]
                    var_type = "DEL"
                QUAL = line[9]
                pass_test = line[13]
                if pass_test == "TRUE":
                    FILTER = "PASS"
                else:
                    FILTER = "FAIL"
                INFO = "DP=" + line[11]
                FORMAT = "GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ"
                SAMPLE = (
                    "1:"
                    + line[4]
                    + ":"
                    + line[5]
                    + ":"
                    + line[6]
                    + ":"
                    + line[7]
                    + ":"
                    + line[8]
                    + ":"
                    + line[9]
                    + ":"
                    + line[10]
                )
                oline = (
                    CHROM
                    + "\t"
                    + POS
                    + "\t"
                    + ID
                    + "\t"
                    + REF
                    + "\t"
                    + ALT
                    + "\t"
                    + QUAL
                    + "\t"
                    + FILTER
                    + "\t"
                    + INFO
                    + "\t"
                    + FORMAT
                    + "\t"
                    + SAMPLE
                    + "\n"
                )
                writeLine = True
                if passOnly and FILTER != "PASS":
                    writeLine = False
                if float(line[10]) < minAF:
                    writeLine = False
                if (CHROM, POS, REF, ALT) in varList:
                    writeLine = False
                else:
                    varList.append((CHROM, POS, REF, ALT))
                if writeLine:
                    varCountDict[var_type] += 1
                    fout.write(oline)
    fout.close()

    ## Print variant counts to pass to MultiQC
    # varCountList = [(k, str(v)) for k, v in sorted(varCountDict.items())]
    # print("\t".join(["sample"] + [x[0] for x in varCountList]))
    # print("\t".join([filename] + [x[1] for x in varCountList]))


def main(args=None):
    args = parse_args(args)
    ivar_variants_to_vcf(
        args.FILE_IN, args.FILE_OUT, args.PASS_ONLY, args.ALLELE_FREQ_THRESH
    )


if __name__ == "__main__":
    sys.exit(main())

import subprocess
import sys
from extract_gb_info import get_id_version, get_locus

CONFIG_FILE = "../workflow/db/snpeff.config"
DATA_DIRECTORY = "../workflow/db/data/"


def prepare_snpeff_run(ref_path_gb, locus, ref_path_fa, reference_name, output):
    """
    The prepare_snpeff_run function creates a snpEFF database for the reference genome.
    It takes as input:
        - The path to the reference genome in GenBank format, and
        - The name of the output file that will be created by this function.

       This function creates a directory called 'config/data' if it does not already exist,
       and then creates three files inside that directory: genes.gbk, sequences.fa,
       and snpeff_genes_prod_snpsiftdb-passed-only-geneidmap
       (the last one is created by running SnpSift).

    :param snpeff_path: Specify the path to snpeff
    :param ref_path_gb: Specify the path to the genbank file of a reference genome
    :param locus: Determine the reference genome to use for snpeff
    :param ref_path_fa: Specify the path to the reference genome fasta file
    :param reference_name: Name the reference files in the snpeff config file
    :param output: Create a file that indicates the database is ready
    :return: A text file with a message "database ready!"
    :doc-author: Trelent
    """
    if len(locus) > 1:
        text = [
            f"{reference_name}.chromosome : {' '.join(str(locus))}\n",
            f"{reference_name}.genome : flu\n",
        ]
        for i in locus:
            text.append(
                f"{reference_name}.{i}.codonTable : Bacterial_and_Plant_Plastid\n"
            )
    else:
        text = [
            f"\n{reference_name}.genome: {reference_name}\n",
            f"{reference_name}.{get_id_version(ref_path_gb)}.codonTable : Bacterial_and_Plant_Plastid\n",
        ]
    with open(CONFIG_FILE, "r", encoding="UTF8") as snpeff_config:
        lines = snpeff_config.readlines()
        joined = "".join(lines)
        find = joined.find("".join(text))

    if find == -1:
        subprocess.run(
            f"mkdir  {DATA_DIRECTORY}{reference_name} -p", shell=True, check=False
        )
        subprocess.run(
            f"cat {ref_path_gb} > {DATA_DIRECTORY}{reference_name}/genes.gbk",
            shell=True,
            check=False,
        )
        subprocess.run(
            f"cat {ref_path_fa} > {DATA_DIRECTORY}{reference_name}/sequences.fa",
            shell=True,
            check=False,
        )
        with open(CONFIG_FILE, "a", encoding="UTF8") as snpeff:
            snpeff.writelines(text)

        subprocess.run(
            f"snpEff build -genbank {reference_name} -c {CONFIG_FILE} > {output}",
            shell=True,
            check=False,
        )
    else:
        with open(output, "w", encoding="UTF8") as out:
            out.write("Database Ready!")


if __name__ == "__main__":
    SNPEFF_PATH = sys.argv[1]
    REFERENCE_GB = sys.argv[2]
    REFERENCE_FASTA = sys.argv[3]
    REFERENCE_NAME = sys.argv[4]
    OUTPUT = sys.argv[5]
    SEGMENTS_LIST = get_locus(REFERENCE_GB)
    prepare_snpeff_run(
        REFERENCE_GB, SEGMENTS_LIST, REFERENCE_FASTA, REFERENCE_NAME, OUTPUT
    )