#!python


"""A tool for bacterial/archaeal RNA-Seq based genome annotations"""
import argparse
import os
from annogesiclib.controller import Controller

__author__ = "Sung-Huan Yu <sung-huan.yu@uni-wuerzburg.de>"
__email__ = "sung-huan.yu@uni-wuerzburg.de"
__version__ = "0.7.14"

def main():
    print("""
       ___    _   ___   ______                  _     
      /   |  / | / / | / / __ \____ ____  _____(_)____ \\
  __ / /| | /  |/ /  |/ / / / / __ `/ _ \/ ___/ / ___/__\\
 |  / ___ |/ /|  / /|  / /_/ / /_/ /  __(__  ) / /__    /
 | /_/  |_/_/ |_/_/ |_/\____/\__, /\___/____/_/\___/   /
 |                          /____/ 
 |__________________
 |_____________________
 |________________________________________________
 |                                                \\
 |________________________________________________/
""")
    home_path = os.environ["HOME"]
    parser = argparse.ArgumentParser()
    parser.add_argument(
        "--version", "-v", default=False, action="store_true",
        help="show version")
    subparsers = parser.add_subparsers(help="commands")
    # Arguments for project creation
    create_project_parser = subparsers.add_parser(
        "create", help="Create a project")
    create_project_basic = create_project_parser.add_argument_group(
            'basic arguments')
    create_project_basic.add_argument(
        "--project_path", "-pj", required=True, help="Name/path of the project.")
    create_project_parser.set_defaults(func=create_project)
    
    # Parameters for get input files
    get_input_parser = subparsers.add_parser(
        "get_input_files", help="Get required files. "
        "(i.e. annotation files, fasta files)")
    get_input_basic = get_input_parser.add_argument_group(
            'basic arguments')
    get_input_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    get_input_basic.add_argument(
        "--ftp_path", "-F", required=True,
        help="Path of folder on the NCBI FTP server where the required files "
        "are located.")
    get_input_basic.add_argument(
        "--ref_fasta", "-f", default=False, action="store_true",
        help="Download fasta files of the reference. Default is False.")
    get_input_basic.add_argument(
        "--ref_gff", "-g", default=False, action="store_true",
        help="Download gff files of the reference. Default is False.")
    get_input_basic.add_argument(
        "--ref_gbk", "-k", default=False, action="store_true",
        help="Download genbank files of the reference. Default is False.")
    get_input_add = get_input_parser.add_argument_group('additional arguments')
    get_input_add.add_argument(
        "--ref_ptt", "-p", default=False, action="store_true",
        help="Download ptt files of the reference. Default is False.")
    get_input_add.add_argument(
        "--ref_rnt", "-r", default=False, action="store_true",
        help="Download rnt files of the reference. Default is False.")
    get_input_add.add_argument(
        "--convert_embl", "-e", default=False, action="store_true",
        help="Convert gbk to embl files of the reference. Default is False.")
    get_input_parser.set_defaults(func=get_input)
    # get target fasta
    get_target_fasta_parser = subparsers.add_parser(
        "update_genome_fasta", help="Get fasta files of reference genomes "
        "if the reference sequences do not exist.")
    get_target_fasta_basic = get_target_fasta_parser.add_argument_group(
            'basic arguments')
    get_target_fasta_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    get_target_fasta_basic.add_argument(
        "--related_fasta_files", "-c", required=True, nargs='+',
        help="Path of the genome fasta files of the closely related species.")
    get_target_fasta_basic.add_argument(
        "--mutation_table", "-m", required=True,
        help="Path of the mutation table which stores the mutation "
        "information between the reference genome and genome of the closely "
        "related species. For an example check "
        "https://github.com/Sung-Huan/ANNOgesic/"
        "blob/master/tutorial_data/mutation.csv")
    get_target_fasta_basic.add_argument(
        "--combine_to_one_fasta", "-cm", default=False, action="store_true",
        help="Combinine all updated sequences in --mutation_table to "
        "one fasta file. Default is False.")
    get_target_fasta_parser.set_defaults(func=get_target_fasta)    

    # run RATT
    RATT_parser = subparsers.add_parser(
        "annotation_transfer", help="Transfer the annotations "
        "from a closely related species genome to a target genome.")
    RATT_basic = RATT_parser.add_argument_group(
        'basic arguments')
    RATT_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    RATT_basic.add_argument(
        "--compare_pair", "-p", required=True, nargs='+',
        help="Please assign the sequence identifier of genome pairs, e.g. "
        "NC_007795:NEW_NC_007795. The related genome is NC_007795 and "
        "the target genome is NEW_NC_007795. The assigned names are the "
        "headers of the fasta file (start with \">\"), not the "
        "filename of fasta file. If the headers contain space or '|', only "
        "the string from '>' to the first space or '|' will be considered "
        "as the name for --compare_pair (normally this part is the "
        "accession number). If multiple sequences need to be assigned, "
        "please use spaces to separate them.")
    RATT_basic.add_argument(
        "--related_embl_files", "-ce", default=None, nargs='+',
        help="The paths of the embl files of the related species. "
        "If --related_embl_files is assigned, --related_gbk_files is not "
        "needed.")
    RATT_basic.add_argument(
        "--related_gbk_files", "-cg", default=None, nargs='+',
        help="The paths of the genbank files of the related species. "
        "The genbank can be ended by .gbk, .gbff or .gb. If "
        "--related_gbk_files is assigned, --related_embl_files is not "
        "needed.")
    RATT_basic.add_argument(
        "--related_fasta_files", "-cf", required=True, nargs='+',
        help="The paths of the fasta files of the related species.")
    RATT_basic.add_argument(
        "--target_fasta_files", "-tf", required=True, nargs='+',
        help="The paths of the target fasta files.")
    RATT_add = RATT_parser.add_argument_group(
        'additional arguments')
    RATT_add.add_argument(
        "--ratt_path", default="start.ratt.sh",
        help="Path of the start.ratt.sh file of RATT folder. "
        "Default is start.ratt.sh.")
    RATT_add.add_argument(
        "--element", "-e", default="chromosome",
        help="--element will become the prefix of all output file."
        "Default is \"chromosome\".")
    RATT_add.add_argument(
        "--transfer_type", "-t", default="Strain",
        help="The transfer type for running RATT. (For the details, "
        "please refer to the manual of RATT.) Default is Strain.")
    RATT_parser.set_defaults(func=run_RATT)

    # Parameters of TSSpredator
    TSSpredator_parser = subparsers.add_parser(
        "tss_ps", help="Detect TSSs or "
        "processing sites.")
    TSSpredator_basic = TSSpredator_parser.add_argument_group(
        'basic arguments')
    TSSpredator_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    TSSpredator_basic.add_argument(
        "--program", "-p", required=True, choices=['TSS', 'PS'],
        help="The feature to predict. Please assign \"TSS\" or "
        "\"PS\". ")
    TSSpredator_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    TSSpredator_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs, etc.")
    TSSpredator_basic.add_argument(
        "--tex_notex_libs", "-tl", required=True, nargs='+',
        help="TEX+/- wig files for TSSpredator. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_tex_forward.wig:tex:1:a:+ "
        "my_lib_tex_reverse.wig:tex:1:a:-.")
    TSSpredator_basic.add_argument(
        "--replicate_tex", "-rt", default=["all_1"], nargs='+',
        help="This value is the minimal number of replicates that a TSS has "
        "to be detected in. The format is "
        "$NUMBERofCONDITION_$NUMBERofREPLICATE. "
        "If different --replicate_tex values need to be assigned to "
        "different conditions, please use spaces to separate them. "
        "For example, 1_2 2_2 3_3 means that --replicate_tex is 2 in "
        "number 1 and number 2 conditions. In number 3 condition, "
        "--replicate_tex is 3. For assigning the same --replicate_tex to all "
        "conditions, just use like all_1 (--replicate_tex is 1 in all "
        "conditions). Default is all_1.")
    TSSpredator_basic.add_argument(
        "--condition_names", "-cn", required=True, nargs='+',
        help="The output prefix of all conditions. If multiple conditions "
        "need to be assigned, please use spaces to separate them. "
        "For an example, prefix_condition1 prefix_condition2.")
    TSSpredator_add = TSSpredator_parser.add_argument_group(
        'additional arguments')
    TSSpredator_add.add_argument(
        "--tsspredator_path", default="/usr/local/bin/TSSpredator.jar",
        help="Path of TSSpredator. Default is /usr/local/bin/TSSpredator.jar")
    TSSpredator_add.add_argument(
        "--height", "-he", default=0.3,
        help="This value relates to the minimal number of read starts at a "
	"certain genomic position to be considered as a TSS candidate. "
        "Default is 0.3.")
    TSSpredator_add.add_argument(
        "--height_reduction", "-rh", default=0.2,
        help="When comparing different genomes/conditions and "
	"the step height threshold is reached in at least one "
        "genome/condition, the threshold is reduced for the other "
        "genomes/conditions by the value set here. "
        "This value must be smaller than the step height threshold. "
        "Default is 0.2.")
    TSSpredator_add.add_argument(
        "--factor", "-fa", default=2.0,
        help="The minimal factor by which the TSS height has to "
	"exceed the local expression background. Default is 2.0.")
    TSSpredator_add.add_argument(
        "--factor_reduction", "-rf", default=0.5,
        help="When comparing different genomes/conditions and "
        "the step factor threshold is reached in at least one "
        "genome/condition, the threshold is reduced for the other "
        "genomes/conditions by the value set here. "
        "This value must be smaller than the step factor threshold. "
        "Default is 0.5.")
    TSSpredator_add.add_argument(
        "--enrichment_factor", "-ef", default=2.0,
        help="The minimal enrichment factor. Default is 2.0.")
    TSSpredator_add.add_argument(
        "--processing_factor", "-pf", default=1.5,
        help="The minimal processing factor. If the value for the "
        "untreated library is higher than the treated library the positions"
        "is considered as a processing site and not annotated as detected. "
        "Default is 1.5.")
    TSSpredator_add.add_argument(
        "--base_height", "-bh", default=0,
        help="The minimal number of reads should be mapped on TSS. "
        "Default is 0.0.")
    TSSpredator_add.add_argument(
        "--utr_length", "-u", default=300, type=int,
        help="The length of UTRs. Default is 300.")
    TSSpredator_add.add_argument(
        "--tolerance", "-to", default=5, type=int,
        help="The 5'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated TSSs ("
        "--compare_transcript_files is provided). Default is 5.")
    TSSpredator_add.add_argument(
        "--cluster", "-c", default=2, type=int,
        help="This value defines the maximal distance (nucleotides) between "
        "TSS candidates have to be clustered together. If the distance "
        "between these multiple TSSs is smaller or equal to this value, "
        "only one of them will be printed out. Default is 2.")
    TSSpredator_add.add_argument(
        "--manual_files", "-m", default=None, nargs='+',
        help="If gff files of the manual checked TSS are provided, this "
        "function will merge manual checked ones and TSSpredator predicted "
        "ones. Please assign the path of manual-checked TSS gff files.")
    TSSpredator_add.add_argument(
        "--curated_sequence_length", "-le", default=None, nargs='+',
        help="The length of the sequence used for the manual set of TSS/PS. "
        "This value is required to calculate the accurracy. If the whole "
        "genome was used write \"all\". Otherwise use the name "
        "of the reference sequence in the folowing format: $GENOME:SLENGTH. "
        "Multiple entries are accepted. For an example, "
        "test.gff contains two sequences s1 and s2. For s1 100 kb were "
        "checked while for s2 the whole sequence was curated. The value of "
        "this argument would be s1:100000 s2:all. Per default all the full "
        "length of all sequences will be used.")
    TSSpredator_add.add_argument(
        "--validate_gene", "-v", default=False, action="store_true",
        help="Using TSS candidates to validate genes in annotation file. "
        "it will be store in statistics folder. Default is False.")
    TSSpredator_add.add_argument(
        "--compare_transcript_files", "-ta", default=None, nargs='+',
        help="If the paths of transcript gff files are provided, this "
        "function will compare TSS and transcript to obtain the overlap "
        "information. Default is False. ")
    TSSpredator_add.add_argument(
        "--re_check_orphan", "-ro", default=False, action="store_true",
        help="If there is no information of gene or locus_tag in genome "
        "annotation gff file, all TSSs will be assigned to orphan TSSs by "
        "TSSpredator. The function can compare TSSs with CDSs to classify "
        "the TSS correctly. Default is False.")
    TSSpredator_add.add_argument(
        "--remove_overlap_feature", "-of", default=False, action="store_true",
        help="If a processing site and a TSS are overlaping, keep \"TSS\", "
        "The predicted feature (based on --program) will be removed. "
        "Default is False.")
    TSSpredator_add.add_argument(
        "--compare_overlap_gff", "-rg", default=None, nargs='+',
        help="If --overlap_feature is \"TSS\" or \"PS\", "
        "--reference_gff_files need to be assigned. For TSS prediction, "
        "please assign the path of processing site. For processing site "
        "prediction, please assign the path of TSS. Don't use this flag "
        "if --overlap_feature is \"both\".")
    TSSpredator_add.add_argument(
        "--remove_low_expression", "-rl", default=None,
        help="For removing the low expressed TSS by comparing the manual "
        "detected TSSs and predicted ones. Please assign the manual-checked "
        "TSS in gff format.")
    TSSpredator_parser.set_defaults(func=run_TSSpredator)
    
    # Parameter of optimization of TSSpredator
    op_TSSpredator_parser = subparsers.add_parser(
        "optimize_tss_ps", help="Optimize TSSs or processing "
        "sites based on manual detected ones.")
    op_TSSpredator_basic = op_TSSpredator_parser.add_argument_group(
        'basic arguments')
    op_TSSpredator_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    op_TSSpredator_basic.add_argument(
        "--program", "-p", default="TSS", choices=['TSS', 'PS'],
        help="The feature for optimization. Please assign \"TSS\" or "
        "\"PS\". Default is TSS.")
    op_TSSpredator_basic.add_argument(
        "--fasta_files", "-f", nargs='+', required=True,
        help="Paths of the fasta file of the reference genome.")
    op_TSSpredator_basic.add_argument(
        "--annotation_files", "-g", nargs='+', required=True,
        help="Paths of the genome annotation gff file containing CDSs, "
        "tRNAs, rRNAs, etc.")
    op_TSSpredator_basic.add_argument(
        "--manual_files", "-m", nargs='+', required=True,
        help="Paths of the manual-checked TSS or PS in gff format. It is used "
        "for benchmarking the prediction during the optimization. The set "
        "should comprise roughly 50 TSS/PS or more.")
    op_TSSpredator_basic.add_argument(
        "--curated_sequence_length", "-le", default=None, nargs='+',
        help="The length of the sequence used for the manual set of TSS/PS. "
        "This value is required to calculate the accurracy. If the whole "
        "genome was used write \"all\". Otherwise use the name "
        "of the reference sequence in the folowing format: $GENOME:SLENGTH. "
        "Multiple entries are accepted. For an example, "
        "test.gff contains two sequences s1 and s2. For s1 100 kb were "
        "checked while for s2 the whole sequence was curated. The value of "
        "this argument would be s1:100000 s2:all. Per default all the full "
        "length of all sequences will be used.")
    op_TSSpredator_basic.add_argument(
        "--tex_notex_libs", "-tl", required=True, nargs='+',
        help="TEX+/- wig files for TSSpredator. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_tex_forward.wig:tex:1:a:+ "
        "my_lib_tex_reverse.wig:tex:1:a:-.")
    op_TSSpredator_basic.add_argument(
        "--replicate_tex", "-rt", default=["all_1"], nargs='+',
        help="This value is the minimal number of replicates that a TSS has "
        "to be detected in. The format is "
        "$NUMBERofCONDITION_$NUMBERofREPLICATE. "
        "If different --replicate_tex values need to be assigned to "
        "different conditions, please use spaces to separate them. "
        "For example, 1_2 2_2 3_3 means that --replicate_tex is 2 in "
        "number 1 and number 2 conditions. In number 3 condition, "
        "--replicate_tex is 3. For assigning the same --replicate_tex to all "
        "conditions, just use like all_1 (--replicate_tex is 1 in all "
        "conditions). Default is all_1.")
    op_TSSpredator_basic.add_argument(
        "--condition_names", "-cn", required=True, nargs='+',
        help="The prefixs of output filename."
        "If multiple conditions need to be assigned, please use spaces to "
        "separate them. For an example, prefix_condition1 prefix_condition2.")
    op_TSSpredator_add = op_TSSpredator_parser.add_argument_group(
        'additional arguments')
    op_TSSpredator_add.add_argument(
        "--tsspredator_path", default="/usr/local/bin/TSSpredator.jar",
        help="Path of TSSpredator. Default is /usr/local/bin/TSSpredator.jar")
    op_TSSpredator_add.add_argument(
        "--max_height", "-he", default=2.5, type=float,
        help="This value relates to the minimum number of read starts at a "
        "certain genomic position to be considered as a TSS candidate. "
        "During optimization, --max_height will be never larger than this "
        "value. Default is 2.5.")
    op_TSSpredator_add.add_argument(
        "--max_height_reduction", "-rh", default=2.4, type=float,
        help="When comparing different genomes/conditions and the step "
        "height threshold is reached in at least one genome/condition, "
        "the threshold is reduced for the other genomes/conditions "
        "by the value set here. This value must be smaller than the step "
        "height threshold. During optimization, --max_height_reduction will "
        "be never larger than this value. Default is 2.4.")
    op_TSSpredator_add.add_argument(
        "--max_factor", "-fa", default=10, type=float,
        help="The minimum factor by which the TSS height has to "
        "exceed the local expression background. During optimization, "
        "--max_factor will be never larger than this value. Default is 10.")
    op_TSSpredator_add.add_argument(
        "--max_factor_reduction", "-rf", default=9.9, type=float,
        help="When comparing different genomes/conditions and the step factor "
        "threshold is reached in at least one genome/condition, "
        "the threshold is reduced for the other genomes/conditions "
        "by the value set here. This value must be smaller than the step "
        "factor threshold. During optimization, --max_factor_reduction will "
        "be never larger than this value. Default is 9.9.")
    op_TSSpredator_add.add_argument(
        "--max_base_height", "-bh", default=0.06, type=float,
        help="The minimum number of reads should be mapped on TSS. "
        "During optimization, --max_base_height will be never larger than "
        "this value. Default is 0.06.")
    op_TSSpredator_add.add_argument(
        "--max_enrichment_factor", "-ef", default=6.0, type=float,
        help="The minimum enrichment factor. "
        "During optimization, --max_enrichment_factor will be never larger "
        "than this value. Default is 6.0.")
    op_TSSpredator_add.add_argument(
        "--max_processing_factor", "-pf", default=6.0, type=float,
        help="The minimum processing factor. If the value for the "
        "untreated library is higher than the treated library the positions"
        "is considered as a processing site and not annotated as detected. "
        "During optimization, --max_processing_factor will be never larger "
        "than this value. Default is 6.0")
    op_TSSpredator_add.add_argument(
        "--utr_length", "-u", default=300, type=int,
        help="The length of UTR. Default is 300.")
    op_TSSpredator_add.add_argument(
        "--cluster", "-cu", default=2, type=int,
        help="This value defines the maximal distance (nucleotides) between "
        "TSS candidates have to be clustered together. Default is 2.")
    op_TSSpredator_add.add_argument(
        "--parallels", "-c", type=int, default=4,
        help="Number of parallel threats for optimization. Default is 4.")
    op_TSSpredator_add.add_argument(
        "--steps", "-s", default=4000, type=int,
        help="Number of tota runs for the optimization. Default is 4000 runs.")
    op_TSSpredator_parser.set_defaults(func=optimize_TSSpredator)
    
    # Parameter of Terminator
    Terminator_parser = subparsers.add_parser(
        "terminator", help="Detect rho-independent terminators.")
    Terminator_basic = Terminator_parser.add_argument_group(
        'basic arguments')
    Terminator_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    Terminator_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    Terminator_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs,"
        " tRNAs, rRNAs, etc. Transcript gff files and sRNA gff files need to "
        "be separately assigned to --transcript_files and --srna_files, "
        "respectively.")
    Terminator_basic.add_argument(
        "--transcript_files", "-a", required=True, nargs='+',
        help="Paths of the transcript gff files.")
    Terminator_basic.add_argument(
        "--tex_notex_libs", "-tl", default=None, nargs='+',
        help="TEX+/- wig files. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_tex_forward.wig:tex:1:a:+ "
        "my_lib_tex_reverse.wig:tex:1:a:-.")
    Terminator_basic.add_argument(
        "--frag_libs", "-fl", default=None, nargs='+',
        help="Wig files of RNA-Seq with transcript fragmented. The format is: "
        "wig_file_path:frag:condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_frag_forward.wig:frag:1:a:+ "
        "my_lib_frag_reverse.wig:frag:1:a:-.")
    Terminator_basic.add_argument(
        "--tex_notex", "-te", default=1, type=int, choices=[1, 2],
        help="The value is for TEX+/- libraries to decide the terminator "
        "should be detected in both (TEX+ and TEX-) or can be detected in "
        "only one library (TEX+ or TEX-). Please assign 1 or 2. Default is 1.")
    Terminator_basic.add_argument(
        "--replicate_tex", "-rt", default=None, nargs='+',
        help="This value is the minimal number of replicates that a TSS has "
        "to be detected in. The format is "
        "$NUMBERofCONDITION_$NUMBERofREPLICATE. "
        "If different --replicate_tex values need to be assigned to "
        "different conditions, please use spaces to separate them. "
        "For example, 1_2 2_2 3_3 means that --replicate_tex is 2 in "
        "number 1 and number 2 conditions. In number 3 condition, "
        "--replicate_tex is 3. For assigning the same --replicate_tex to all "
        "conditions, just use like all_1 (--replicate_tex is 1 in all "
        "conditions). Default is all_1.")
    Terminator_basic.add_argument(
        "--replicate_frag", "-rf", default=None, nargs='+',
        help="Similar to --replicates_tex. "
        "This value is for fragmented (or conventional) libraries.")
    Terminator_add = Terminator_parser.add_argument_group(
        'additional arguments')
    Terminator_add.add_argument(
        "--transterm_path", default="transterm",
        help="Path of \"transterm\" in TransTermHP.")
    Terminator_add.add_argument(
        "--expterm_path", default="/usr/local/bin/expterm.dat",
        help="Path of expterm.dat for TransTermHP. "
        "Default is /usr/local/bin/expterm.dat")
    Terminator_add.add_argument(
        "--rnafold_path", default="RNAfold",
        help="Path of RNAfold of Vienna package.")
    Terminator_add.add_argument(
        "--srna_files", "-sr", default=None, nargs='+',
        help="Paths of sRNA gff files if sRNA information "
        "need to be considered as well.")
    Terminator_add.add_argument(
        "--decrease", "-d", default=0.5, type=float,
        help="The maximum ratio -- (lowest coverage / highest "
        "coverage) within (or nearby) the terminator. If the ratio is smaller "
        "than --decrease, the candidate will be considered as highly-confidence "
        "terminator. Default is 0.5.")
    Terminator_add.add_argument(
        "--tolerance_detect_coverage", "-tc", default=30, type=int,
        help="The extended region (nucleotides) of the "
        "terminators for detecting coverage significant drop. "
        "For example, the location of terminator is 300-400, and "
        "--tolerance_detect_coverage is 30. If the coverage decrease is detected "
        "within 270-430, this candidate is still considered as the "
        "terminator which have coverage dramatic decrease. Default is 30.")
    Terminator_add.add_argument(
        "--tolerance_within_transcript", "-tut", default=30, type=int,
        help="If the candidates are within transcript and the distance "
        "(nucleotides) between the end of transcript and terminator is "
        "within this value, the candidate will be considered as a terminator. "
        "Otherwise, it will be removed. Default is 30.")
    Terminator_add.add_argument(
        "--tolerance_downstream_transcript", "-tdt", default=30, type=int,
        help="The meaning is similar to --tolerance_within_transcript. "
        "This value is for the candidates which are at the downstream of "
        "transcript. Default is 30.")
    Terminator_add.add_argument(
        "--tolerance_within_gene", "-twg", default=10, type=int,
        help="The meaning is similar to --tolerance_within_transcript. "
        "This value is for gene in stead of transcript. Default is 10.")
    Terminator_add.add_argument(
        "--tolerance_downstream_gene", "-tdg", default=310, type=int,
        help="The meaning is similar to --tolerance_downstream_transcript. "
        "This value is for gene in stead of transcript. Default is 310.")
    Terminator_add.add_argument(
        "--highest_coverage", "-hc", default=10, type=float,
        help="The minimum value of the highest coverage of terminator. "
        "The low expressed terminator are not included in "
        "\"best_candidates\", but are still in \"all_candidates\". Default is 10.")
    Terminator_add.add_argument(
        "--window_size", "-wz", default=100, type=int,
        help="Window size for searching secondary structure in intergenic "
        "region. Default is 100 nts.")
    Terminator_add.add_argument(
        "--window_shift", "-ws", default=20, type=int,
        help="The number of nucleotides for window shift. Default is 20 nts.")
    Terminator_add.add_argument(
        "--min_loop_length", "-ml", default=3, type=int,
        help="The minimum loop length of terminator. Default is 3 nts.")
    Terminator_add.add_argument(
        "--max_loop_length", "-Ml", default=10, type=int,
        help="The maximum loop length of terminator. Default is 10 nts.")
    Terminator_add.add_argument(
        "--min_stem_length", "-ms", default=4, type=int,
        help="The minimum stem length of terminator. Default is 4 nts.")
    Terminator_add.add_argument(
        "--max_stem_length", "-Ms", default=20, type=int,
        help="The maximum stem length of terminator. Default is 20 nts.")
    Terminator_add.add_argument(
        "--miss_rate", "-mr", default=0.25, type=float,
        help="The percentage of nucleotides which can be no pair in the stem. "
        "Default is 0.25.")
    Terminator_add.add_argument(
        "--min_u_tail", "-mu", default=5, type=int,
        help="The minimum number of U in poly U-tail of terminator. Default is 5.")
    Terminator_add.add_argument(
        "--mutation_u_tail", "-uu", default=2, type=int,
        help="The number of nts which are not U can be tolerated. Default is 2.")
    Terminator_add.add_argument(
        "--keep_multi_term", "-kp", default=False, action="store_true",
        help="Sometimes, one gene is associated with multiple terminators "
        "In default, it will only keep the highly-confidence one. "
        "This flag can keep all terminators which are associated with the "
        "same gene. Default is False.")
    Terminator_parser.set_defaults(func=run_Terminator)

    # Parameter of Transcript
    Transcript_parser = subparsers.add_parser(
        "transcript", help="Detect transcripts based on "
        "coverage file.")
    Transcript_basic = Transcript_parser.add_argument_group(
        'basic arguments')
    Transcript_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    Transcript_basic.add_argument(
        "--annotation_files", "-g", default=None, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs, etc. TSS gff files and terminator gff files need to be "
        "separately assigned to --tss_files and --terminator_files, respectively.")
    Transcript_basic.add_argument(
        "--modify_transcript", "-mt", default=["merge_overlap"], nargs='+',
        choices=["merge_overlap", "extend_3end", "extend_5end",
                 "within_extend_ends", "none"],
        help="If --annotation_files is provided, the post-modification of "
        "transcript based on genome annotations can be assigned. There are "
        "five opetions. 1. \"merge_overlap\": if multiple transcripts overlap "
        "with the same gene, they will be merged as one complete transcript. "
        "2. \"extend_3end\": if a transcript starts at the "
        "upstream of a gene and ends within the gene, the end point of the "
        "transcript will be extended to the end point of the gene. "
        "3. \"extend_5end\": if a transcript starts within a gene and "
        "ends at the downstream of gene, the starting point of the "
        "transcript will be extended to the starting point of the gene. "
        "4. \"within_extend_ends\": if a transcript is within a gene, "
        "the two ends of the transcript will be extended to the two ends of "
        "gene. 5. \"none\": the transcript will not be modified by the genome "
        "annotations. For using mutliple modifications, please separated them "
        "by spaces. Default is merge_overlapped.")
    Transcript_basic.add_argument(
        "--tex_notex_libs", "-tl", default=None, nargs='+',
        help="TEX+/- wig files. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_tex_forward.wig:tex:1:a:+ "
        "my_lib_tex_reverse.wig:tex:1:a:-.")
    Transcript_basic.add_argument(
        "--frag_libs", "-fl", default=None, nargs='+',
        help="Wig files of RNA-Seq with transcript fragmented. The format is: "
        "wig_file_path:frag:condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_frag_forward.wig:frag:1:a:+ "
        "my_lib_frag_reverse.wig:frag:1:a:-.")
    Transcript_basic.add_argument(
        "--replicate_tex", "-rt", default=None, nargs='+',
        help="This value is the minimal number of replicates that a TSS has "
        "to be detected in. The format is "
        "$NUMBERofCONDITION_$NUMBERofREPLICATE. "
        "If different --replicate_tex values need to be assigned to "
        "different conditions, please use spaces to separate them. "
        "For example, 1_2 2_2 3_3 means that --replicate_tex is 2 in "
        "number 1 and number 2 conditions. In number 3 condition, "
        "--replicate_tex is 3. For assigning the same --replicate_tex to all "
        "conditions, just use like all_1 (--replicate_tex is 1 in all "
        "conditions). Default is all_1.")
    Transcript_basic.add_argument(
        "--replicate_frag", "-rf", default=None, nargs='+',
        help="Similar to --replicates_tex. "
        "This value is for fragmented (or conventional) libraries.")
    Transcript_basic.add_argument(
        "--tex_notex", "-te", default=1, type=int, choices=[1, 2],
        help="The value is for TEX+/- libraries to decide the transcript "
        "should be detected in both (TEX+ and TEX-) or can be detected in "
        "only one library (TEX+ or TEX-). Please assign 1 or 2. Default is 1.")
    Transcript_add = Transcript_parser.add_argument_group(
        'additional arguments')
    Transcript_add.add_argument(
        "--length", "-l", default=20, type=int,
        help="The minimum length of the transcript after modifying base on "
        "genome annotation. Default is 20.")
    Transcript_add.add_argument(
        "--height", "-he", default=10, type=int,
        help="The minimum coverage of the transcript. If --tex_notex is 2, "
        "the average coverage of TEX+ and TEX- libraries should be higher "
        "than this value. The default is 10.")
    Transcript_add.add_argument(
        "--width", "-w", default=20, type=int,
        help="The minimum length of the transcript without modifying by "
        "genome annotation. The default is 20.")
    Transcript_add.add_argument(
        "--tolerance", "-t", default=5, type=int,
        help="The number of nucleotides which coveraes can drop below the "
        "--height in a transcript. The default is 5.")
    Transcript_add.add_argument(
        "--tolerance_coverage", "-tc", default=0, type=int,
        help="The minimum coverage of the nucleotides which match the "
        "situation of --tolerance, Default is 0.")
    Transcript_add.add_argument(
        "--tss_files", "-ct", default=None, nargs='+',
        help="Paths of TSS files for comparing transcripts and TSSs.")
    Transcript_add.add_argument(
        "--compare_feature_genome", "-cf", default=None, nargs='+',
        help="If --compare_genome_annotation is provided, please assign the "
        "feature for comparing. Multiple features can be separated by spaces. "
        "Default is None.")
    Transcript_add.add_argument(
        "--tss_tolerance", "-tt", default=5, type=int,
        help="The 5'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated TSSs (--tss_files is "
        "provided). Default is 5.")
    Transcript_add.add_argument(
        "--terminator_files", "-e", default=None, nargs='+',
        help="Paths of terminator gff files for comparing transcripts and "
        "terminators. Default is None.")
    Transcript_add.add_argument(
        "--terminator_tolerance", "-et", default=30, type=int,
        help="The 3'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated terminators. "
        "Default is 30.")
    Transcript_add.add_argument(
        "--max_length_distribution", "-mb", default=2000, type=int,
        help="For generating the figure of distribution of transcript length, "
        "please assign the maximum length. Default is 2000.")
    Transcript_parser.set_defaults(func=run_Transcript_Assembly)

    # Parameter of UTR detection
    UTR_parser = subparsers.add_parser(
        "utr", help="Detect 5'UTRs and 3'UTRs.")
    UTR_basic = UTR_parser.add_argument_group(
        'basic arguments')
    UTR_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    UTR_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs, etc. Gff files of TSSs, terminators, and transcripts need "
        "to be separately assigned to --tss_files, terminator_files, and "
        "transcript_files, respectively.")
    UTR_basic.add_argument(
        "--tss_files", "-t", required=True, nargs='+',
        help="Paths of the TSS files.")
    UTR_basic.add_argument(
        "--transcript_files", "-a", required=True, nargs='+',
        help="Paths of the transcript gff files.")
    UTR_basic.add_argument(
        "--terminator_files", "-e", default=None, nargs='+',
        help="If the paths of terminator files are assigned, "
        "terminator will be used to detect 3'UTR.")
    UTR_add = UTR_parser.add_argument_group(
        'additional arguments')
    UTR_add.add_argument(
        "--tss_source", "-s", default=True, action="store_false",
        help="The TSS gff file is generated by ANNOgesic or not. "
        "Default is True (from ANNOgesic).")
    UTR_add.add_argument(
        "--base_5utr", "-b5", default="both",
        choices=['both', 'transcript', 'TSS'],
        help="The information for detection of 5'UTR. "
        "It can be \"TSS\" or \"transcript\" or \"both\". Default is both.")
    UTR_add.add_argument(
        "--utr_length", "-l", default=300, type=int,
        help="The maximum UTR length. Default is 300.")
    UTR_add.add_argument(
        "--base_3utr", "-b3", default="transcript",
        choices=['both', 'transcript', 'terminator'],
        help="please assign the information for detection of 3'UTR. It can "
        "be \"transcript\" or \"terminator\" or \"both\". "
        "Default is transcript.")
    UTR_add.add_argument(
        "--terminator_tolerance", "-et", default=30, type=int,
        help="The 3'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated terminators. "
        "Default is 30.")
    UTR_add.add_argument(
        "--tolerance_3utr", "-t3", default=10, type=int,
        help="The length of 3'UTR can be extended or withdrew by this "
        "value (nucleotides). It only works when transcript information "
        "is provided. Default is 10.")
    UTR_add.add_argument(
        "--tolerance_5utr", "-t5", default=5, type=int,
        help="The length of 5'UTR can be extended or withdrew by this "
        "value (nucleotides). It only works when transcript information "
        "is provided. Default is 5.")
    UTR_parser.set_defaults(func=run_UTR_detection)

    # Parameter of sRNA detection
    sRNA_parser = subparsers.add_parser(
        "srna", help="Detect intergenic, antisense and UTR-derived sRNAs.")
    sRNA_basic = sRNA_parser.add_argument_group(
        'basic arguments')
    sRNA_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    sRNA_basic.add_argument(
        "--utr_derived_srna", "-u", default=False, action="store_true",
        help="Assign to detect UTR-derived sRNA. Default is False.")
    sRNA_basic.add_argument(
        "--filter_info", "-d",
        default=["tss", "sec_str", "blast_nr", "blast_srna"], nargs='+',
        choices=['tss', 'sec_str', 'blast_nr', 'blast_srna', 'sorf', 'term',
                 'promoter', 'none'],
        help="The filters for improving the sRNA "
        "detection: 1. tss (sRNA has to start with a TSS), 2. sec_str "
        "(free energy change of secondary structure (normalized by length) "
        "has to be smaller than --cutoff_energy), 3. blast_nr (the number "
        "of the homologs in the non-redundant database has to be below the "
        "--cutoff_nr_hit ), 4. blast_srna (as long as the homologs "
        "can be found in the sRNA database, the candidates will be included to "
        "the best candidtes without considering other filters), 5. sorf (sRNA "
        "must not overlap with sORFs), 6. term (sRNA has to be associated with "
        "a terminator), 7. promoter (sRNA has to be associated with a promoter "
        "motif). For using multiple filters, please separated them by spaces. "
        "If blast_srna was assigned, the headers of sequences in sRNA "
        "database should be $ID|$GENOME|$SRNANAME. "
        "\"tss sec_str blast_nr blast_srna\" are recommended to be used. "
        "If \"none\" is assigned, no filters "
        "are applied. Default is tss sec_str blast_nr blast_srna.")
    sRNA_basic.add_argument(
        "--transcript_files", "-a", required=True, nargs='+',
        help="Paths of the transcript files.")
    sRNA_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs, etc. Gff files of transcripts, TSSs, terminators, "
        "processing sites, and sORFs need to be separately assigned to "
        "--transcript_files, --tss_files, --terminator_files, "
        "--processing_site_files, and --sorf_files, respectively.")
    sRNA_basic.add_argument(
        "--tss_files", "-t", default=None, nargs='+',
        help="Paths of TSS gff files. For detecting UTR-derived sRNA or "
        "\"tss\" in --filter_info, TSS gff files MUST be provided.")
    sRNA_basic.add_argument(
        "--fasta_files", "-f", default=None, nargs='+',
        help="paths of fasta files of reference genome, If \"sec_str\", "
        "\"blast_nr\" or \"blast_srna\" is assigned to "
        "--filter_info or --search_poly_u is not 0, fasta files are required. ")
    sRNA_basic.add_argument(
        "--compute_sec_structures", "-cs", default=False, action="store_true",
        help="Computing secondary structures of sRNAs. Default is False.")
    sRNA_basic.add_argument(
        "--srna_database_path", "-sd", default=None,
        help="Path of sRNA database with proper headers of sequences. "
        "Format of the header should be $ID|$GENOME|$NAME. "
        "Please check https://github.com/Sung-Huan/ANNOgesic/blob/"
        "master/database/sRNA_database_BSRD.fa")
    sRNA_basic.add_argument(
        "--nr_database_path", "-nd", default=None,
        help="Path of nr database")
    sRNA_basic.add_argument(
        "--tex_notex_libs", "-tl", default=None, nargs='+',
        help="TEX+/- wig files. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_tex_forward.wig:tex:1:a:+ "
        "my_lib_tex_reverse.wig:tex:1:a:-.")
    sRNA_basic.add_argument(
        "--frag_libs", "-fl", default=None, nargs='+',
        help="Wig files of RNA-Seq with fragmented transcripts. The format is: "
        "wig_file_path:frag:condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_frag_forward.wig:frag:1:a:+ "
        "my_lib_frag_reverse.wig:frag:1:a:-.")
    sRNA_basic.add_argument(
        "--tex_notex", "-te", default=2, type=int, choices=[1, 2],
        help="If TEX+/- libraries are assigned, a sRNA "
        "should be detected in both (TEX+ and TEX-) or needs to be detected "
        " in only one library (TEX+ or TEX-). Default is 2.")
    sRNA_basic.add_argument(
        "--replicate_tex", "-rt", default=None, nargs='+',
        help="This value is the minimal number of replicates that a TSS has "
        "to be detected in. The format is "
        "$NUMBERofCONDITION_$NUMBERofREPLICATE. "
        "If different --replicate_tex values need to be assigned to "
        "different conditions, please use spaces to separate them. "
        "For example, 1_2 2_2 3_3 means that --replicate_tex is 2 in "
        "number 1 and number 2 conditions. In number 3 condition, "
        "--replicate_tex is 3. For assigning the same --replicate_tex to all "
        "conditions, use all_1 (--replicate_tex is 1 in all "
        "conditions). Default is all_1.")
    sRNA_basic.add_argument(
        "--replicate_frag", "-rf", default=None, nargs='+',
        help="Similar to --replicates_tex. "
        "This value is for libraries with fragmented transcripts.")
    sRNA_add = sRNA_parser.add_argument_group(
        'additional arguments')
    sRNA_add.add_argument(
        "--rnafold_path", default="RNAfold",
        help="Path of RNAfold of the Vienna package")
    sRNA_add.add_argument(
        "--relplot_path", default="relplot.pl",
        help="Path of relplot.pl of the Vienna package.")
    sRNA_add.add_argument(
        "--mountain_path", default="mountain.pl",
        help="Path of mountain.pl of the  Vienna package.")
    sRNA_add.add_argument(
        "--blastn_path", default="blastn",
        help="Path of blastn of the  BLAST+ package.")
    sRNA_add.add_argument(
        "--blastx_path", default="blastx",
        help="Path of blastx of the BLAST+ package.")
    sRNA_add.add_argument(
        "--makeblastdb_path", default="makeblastdb",
        help="Path of makeblastdb of the BLAST+ package.")
    sRNA_add.add_argument(
        "--processing_site_files", "-p", default=None, nargs='+',
        help="Paths of the processing site gff files. It can improve the "
        "detection of UTR-derived sRNAs.")
    sRNA_add.add_argument(
        "--terminator_files", "-e", default=None, nargs='+',
        help="Paths of the terminator gff files.")
    sRNA_add.add_argument(
        "--promoter_tables", "-pt", default=None, nargs='+',
        help="If promoter tables can be provided, please assign the paths of "
        "promoter tables, The format of the table is "
        "$GENOME\t$TSS_POSITION\t$TSS_STRAND\t$PROMOTER_NAME.")
    sRNA_add.add_argument(
        "--promoter_names", "-pn", default=None, nargs='+',
        help="If --promoter_tables is provided, please assign the promoter "
        "name (the last column of promoter table). For multiple promoters, "
        "please put spaces between the promoters. Default is None.")
    sRNA_add.add_argument(
        "--sorf_files", "-O", default=None, nargs='+',
        help="If comparison between sRNAs and sORFs needs to be done, "
        "Please assign the paths of sORF gff files")
    sRNA_add.add_argument(
        "--parallel_blast", "-pb", type=int, default=10,
        help="The number of parallel jobs. Default is 10.")
    sRNA_add.add_argument(
        "--tss_source", "-ts", default=True, action="store_false", 
        help="The TSS gff files are generated from ANNOgesic or not. "
        "Default is True (from ANNOgesic).")
    sRNA_add.add_argument(
        "--tss_intergenic_antisense_tolerance", "-tit", default=3, type=int,
        help="The 5'ends of intergenic and antisense sRNA candidates will be "
        "extended or withdrew by this value (nucleotides) for searching the "
        "associated TSSs. Default is 3.")
    sRNA_add.add_argument(
        "--tss_5utr_tolerance", "-t5", default="n_3", type=str,
        help="The 5'ends of 5'UTR-derived sRNAs will be extended or withdrew "
        "by this value (nucleotides) for searching the associated TSSs. The input "
        "type can be percentage (\"p\") or the real amount of reads "
        "(\"n\"). For example, p_0.05 means this value is 5 percent of the "
        "length of 5'UTR. n_10 means this value is 10 nts. Default is n_3.")
    sRNA_add.add_argument(
        "--tss_3utr_tolerance", "-t3", default="p_0.04", type=str,
        help="Similar to --tss_5utr_tolerance. This value is for "
        "3'UTR-derived sRNAs. Default is p_0.04.")
    sRNA_add.add_argument(
        "--tss_intercds_tolerance", "-tc", default="p_0.04", type=str,
        help="Similar to --tss_5utr_tolerance. This value is for "
        "interCDS-derived sRNAs. Default is p_0.04.")
    sRNA_add.add_argument(
        "--terminator_tolerance_in_srna", "-eti", default=30, type=int,
        help="The 3'ends of sRNA candidates will be withdrew by this "
        "value (nucleotides) for searching the associated terminators which "
        "are within sRNAs. Default is 30.")
    sRNA_add.add_argument(
        "--terminator_tolerance_out_srna", "-eto", default=30, type=int,
        help="The 3'ends of sRNA candidates will be extended by this "
        "value (nucleotides) for searching the associated terminators which "
        "are behind of sRNAs. Default is 30.")
    sRNA_add.add_argument(
        "--min_length", "-lm",default=30, type=int,
        help="The minimum sRNA length. Default is 30.")
    sRNA_add.add_argument(
        "--max_length", "-lM",default=500, type=int,
        help="The maximum sRNA length. Default is 500.")
    sRNA_add.add_argument(
        "--min_intergenic_tex_coverage", "-it",default="0,0,0,40,20",
        help="The minimum average coverage of intergenic sRNAs for "
        "TEX+ libraries. This value is based on different types of TSSs. "
        "The order of numbers is \"Primary,Secondary,Internal,Antisense,"
        "Orphan\". For example, 0,0,0,50,10 means that antisense TSS "
        "(minimum coverage is 50) and orphan TSS (minimum coverage is 10) "
        "are used for sRNA prediction. The other types of TSSs will not be "
        "used for the detection (assign to 0). If TSS information is not "
        "provided, the lowest value would be the general cutoff for the "
        "prediction. Default is 0,0,0,40,20.")
    sRNA_add.add_argument(
        "--min_intergenic_notex_coverage", "-in",default="0,0,0,30,10",
        help="Similar to --min_intergenic_tex_coverage. "
        "This value is for TEX- libraries. Default is 0,0,0,30,10.")
    sRNA_add.add_argument(
        "--min_intergenic_fragmented_coverage", "-if",default="400,200,0,50,20",
        help="Similar to --min_intergenic_tex_coverage. "
        "This value is for fragmented (or conventional) libraries. "
        "Default is 400,200,0,50,20.")
    sRNA_add.add_argument(
        "--min_complete_5utr_transcript_coverage", "-ib",default="30,20,30",
        help="Several primary/secondary TSSs are also associated with a "
        "complete transcript containing no CDSs/tRNA/rRNA in 5'UTR of the "
        "following CDS which is located in another transcript. "
        "In order to detect the sRNA candidates in these transcripts, please "
        "assign the minimum average coverage of the sRNA candidates. "
        "The format is $TEX,$NOTEX,$FRAG. For example, 200,100,100 means that "
        "the minimum average coverage is 200 for TEX+ libraries, 100 for TEX- "
        "and fragmented (or conventional) libraries. Default is 30,20,30.")
    sRNA_add.add_argument(
        "--min_antisense_tex_coverage", "-at",default="0,0,0,40,20",
        help="Similar to --min_intergenic_tex_coverage. "
        "This value is for antisense in stead of intergenic. "
        "Default is 0,0,0,40,20.")
    sRNA_add.add_argument(
        "--min_antisense_notex_coverage", "-an",default="0,0,0,30,10",
        help="Similar to --min_intergenic_notex_coverage. "
        "This value is for antisense in stead of intergenic. "
        "Default is 0,0,0,30,10.")
    sRNA_add.add_argument(
        "--min_antisense_fragmented_coverage", "-af",default="400,200,0,50,20",
        help="Similar to --min_intergenic_fragmented_coverage. "
        "This value is for antisense in stead of intergenic. "
        "Default is 400,200,0,50,20.")
    sRNA_add.add_argument(
        "--min_utr_tex_coverage", "-ut",default="p_0.8,p_0.6,p_0.7",
        help="The minimum average coverage of UTR-derived sRNA candidates in "
        "TEX+ libraries. The input can be assigned by the percentile (\"p\") "
        "or real number of coverage (\"n\"). The order of numbers is "
        "\"5'UTR,3'UTR,interCDS\". For example, \"p_0.7,p_0.5,p_0.5\" means "
        "that 70 percentile of all 5'UTR coverages is used for the minimum "
        "coverage of 5'UTR-derived sRNA, median of all 3'UTR and interCDS "
        "coverages is used for minimum coverage of 3'UTR and interCDS-derived "
        "sRNA. Default is p_0.8,p_0.6,p_0.7.")
    sRNA_add.add_argument(
        "--min_utr_notex_coverage", "-un",default="p_0.7,p_0.5,p_0.6",
        help="Similar to --min_utr_tex_coverage. "
        "This value is for TEX- libraries. Default is p_0.7,p_0.5,p_0.6.")
    sRNA_add.add_argument(
        "--min_utr_fragmented_coverage", "-uf",default="p_0.7,p_0.5,p_0.6",
        help="Similar to --min_utr_tex_coverage. "
        "This value is for fragmented (or conventional) libraries. "
        "Default is p_0.7,p_0.5,p_0.6.")
    sRNA_add.add_argument(
        "--min_all_utr_coverage", "-mu", default=50, type=float,
        help="The minimum coverage of UTR-derived sRNAs. "
        "The coverage of UTR-derived sRNAs should not only exceed the "
        "--min_utr_TEX_coverage, --min_utr_noTEX_coverage and "
        "--min_utr_fragmented_coverage, "
        "but also this value. Default is 50.")
    sRNA_add.add_argument(
        "--cutoff_energy", "-ce", default=-0.05, type=float,
        help="If \"sec_str\" is included in --filter_info, please assign the "
        "maximum folding energy change (normalized by length of gene). "
        "Default is -0.05.")
    sRNA_add.add_argument(
        "--mountain_plot", "-m", default=False, action="store_true",
        help="Generating mountain plot of sRNA candidate. "
        "Default is False.")
    sRNA_add.add_argument(
        "--nr_format", "-nf", default=False, action="store_true",
        help="Format nr database. Default is False.")
    sRNA_add.add_argument(
        "--srna_format", "-sf", default=False, action="store_true",
        help="Format sRNA database. Default is False.")
    sRNA_add.add_argument(
        "--decrease_intergenic_antisense","-di", default=0.1, type=float,
        help="This value is for detecting the coverage decrease in "
        "intergenic/antisense transcript. "
        "If the length of intergenic transcript is longer than the --max_length, "
        "it will searching the coverages of the transcript. If the ratio "
        "-- (the lowest coverage / the highest coverage) of the transcript "
        "is smaller than this value and the length is within a given range, "
        "the transcript will be considered as a sRNA as well."
        "Default is 0.1.")
    sRNA_add.add_argument(
        "--decrease_utr","-du", default=0.05, type=float,
        help="Similar to --decrease_intergenic_antisense. This value is "
        "for UTR-derived sRNA. Default is 0.05.")
    sRNA_add.add_argument(
        "--tolerance_intergenic_antisense", "-ti", default=10, type=int,
        help="The 5'ends and 3'ends of intergenic and antisense sRNAs will be "
        "extended by this value (nucleotides) for detecting the significant "
        "coverage decrease. (please check --decrease_intergenic_antisense). "
        "For example, the location of intergenic sRNA is 300-400, and "
        "--tolerance_intergenic_antisense is 30. The searching region is 270-430. "
        "Default is 10.")
    sRNA_add.add_argument(
        "--tolerance_utr", "-tu", default=10, type=int,
        help="Similar to --tolerance_intergenic_antisense. "
        "This is for UTR-derived sRNAs. Default is 10.")
    sRNA_add.add_argument(
        "--cutoff_nr_hit", "-cn", default=0, type=int,
        help="The maximum hits number in nr database. Default is 0.")
    sRNA_add.add_argument(
        "--blast_e_nr", "-en", default=0.0001, type=float,
        help="The maximum e-value for searching in nr database. Default is 0.0001.")
    sRNA_add.add_argument(
        "--blast_e_srna", "-es", default=0.0001, type=float,
        help="The maximum e-value for searching in sRNA database. "
        "Default is 0.0001.")
    sRNA_add.add_argument(
        "--blast_score_srna", "-bs", default=None, type=float,
        help="The minimum score for searching in sRNA database. "
        "Default is None.")
    sRNA_add.add_argument(
        "--blast_score_nr", "-bn", default=None, type=float,
        help="The minimum score for searching in nr database. "
        "Default is None.")
    sRNA_add.add_argument(
        "--detect_srna_in_cds", "-ds", default=False, action="store_true",
        help="Searching sRNA in CDS (e.g. the genome "
        "annotation is not correct). More sRNA candidates which overlap "
        "with CDS will be detected. Default is False.")
    sRNA_add.add_argument(
        "--overlap_percent_cds", "-oc", default=0.5,
        help="The maximum ratio of overlapping between CDS and sRNA candidates. "
        "It only works if --detect_srna_in_cds is true. Default is 0.5")
    sRNA_add.add_argument(
        "--search_poly_u", "-sp", default=15, type=int,
        help="The tolerance length for searching poly U tail of sRNA. "
        "If this value is assigned by 0, the 3'end of sRNA will not be "
        "extended by searching poly U tail. Default is 15.")
    sRNA_add.add_argument(
        "--min_u_poly_u", "-np", default=5, type=int,
        help="The minimum number of U that poly U tail should contain. "
        "Default is 5.")
    sRNA_add.add_argument(
        "--mutation_poly_u", "-mp", default=2, type=int,
        help="The minimum number of nts which are not U can be tolerated. "
        "Default is 2.")
    sRNA_add.add_argument(
        "--ignore_hypothetical_protein", "-ih", default=False,
        action="store_true",
        help="For ignoring hypothetical proteins in "
        "the genome annotation file. Default is False.")
    sRNA_add.add_argument(
        "--ranking_time_promoter", "-rp", default=2, type=float,
        help="If --promoter_tables is provided, the information of promoter "
        "can be use for ranking sRNA candidates. "
        "The ranking score is --ranking_time_promoter * average coverage. "
        "For example, a sRNA candidate which is associated with a promoter and "
        "its average coverage is 10. If --ranking_time_promoter is 2, "
        "the ranking score will be 20 (2*10). "
        "For the candidate which are not associated with a promoter, the "
        "--ranking_time_promoter will be 1. Therefore, --ranking_time_promoter "
        "can not be smaller than 1. Default is 2.")
    sRNA_add.add_argument(
        "--exclude_srna_in_annotation_file", "-ea", default=False,
        action="store_true",
        help="For excluding the sRNAs which are already annotated in "
        "--annotation_files. Default is False.")
    sRNA_parser.set_defaults(func=run_sRNA_detection)
    
    # Parameters of small ORF
    sORF_parser = subparsers.add_parser(
        "sorf", help="Detect expressed sORFs.")
    sORF_basic = sORF_parser.add_argument_group(
        'basic arguments')
    sORF_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    sORF_basic.add_argument(
        "--utr_derived_sorf", "-u", default=False, action="store_true",
        help="Detect UTR-derived sORF. Default is False.")
    sORF_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the fasta files of the reference genome. ")
    sORF_basic.add_argument(
        "--transcript_files", "-a", required=True, nargs='+',
        help="Paths of the transcript gff files.")
    sORF_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs, etc. Gff files of transcripts, sRNAs, and TSSs need to "
        "be separately assigned to --transcript_files, --srna_files, and "
        "--tss_files, respectively.")
    sORF_basic.add_argument(
        "--tss_files", "-t", default=None, nargs='+',
        help="Paths of TSS gff files.")
    sORF_basic.add_argument(
        "--tex_notex_libs", "-tl", default=None, nargs='+',
        help="TEX+/- wig files. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_tex_forward.wig:tex:1:a:+ "
        "my_lib_tex_reverse.wig:tex:1:a:-.")
    sORF_basic.add_argument(
        "--frag_libs", "-fl", default=None, nargs='+',
        help="Wig files of RNA-Seq with fragmented transcripts. "
        "The format is: wig_file_path:frag:condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_frag_forward.wig:frag:1:a:+ "
        "my_lib_frag_reverse.wig:frag:1:a:-.")
    sORF_basic.add_argument(
        "--tex_notex", "-te", default=2, type=int, choices=[1, 2],
        help="If the TEX+/- libraries are provided, this value is that a "
        "sORF should be detected in both (TEX+ and TEX-) or "
        "can be detected in only one library (TEX+ or TEX-). "
        "Please assign 1 or 2. Default is 2.")
    sORF_basic.add_argument(
        "--replicate_tex", "-rt", default=None, nargs='+',
        help="This value is the minimal number of replicates that a TSS has "
        "to be detected in. The format is "
        "$NUMBERofCONDITION_$NUMBERofREPLICATE. "
        "If different --replicate_tex values need to be assigned to "
        "different conditions, please use spaces to separate them. "
        "For example, 1_2 2_2 3_3 means that --replicate_tex is 2 in "
        "number 1 and number 2 conditions. In number 3 condition, "
        "--replicate_tex is 3. For assigning the same --replicate_tex to all "
        "conditions, just use like all_1 (--replicate_tex is 1 in all "
        "conditions). Default is all_1.")
    sORF_basic.add_argument(
        "--replicate_frag", "-rf", default=None, nargs='+',
        help="Similar to --replicates_tex. "
        "This value is for fragmented (or conventional) libraries.")
    sORF_add = sORF_parser.add_argument_group(
        'additional arguments')
    sORF_add.add_argument(
        "--srna_files", "-s", default=None, nargs='+',
        help="Paths of the sRNA gff files for comparing sORF and sRNA to "
        "detect the overlapping.")
    sORF_add.add_argument(
        "--utr_length", "-ul", default=300, type=int,
        help="The utr length for comparing TSS with sORF. "
        "The default number is 300.")
    sORF_add.add_argument(
        "--min_length", "-ml",default=30, type=int,
        help="The minimum nucleotide length of sORF. Default is 30.")
    sORF_add.add_argument(
        "--max_length", "-Ml",default=300, type=int,
        help="The maximum nucleotide length of sORF. Default is 300.")
    sORF_add.add_argument(
        "--cutoff_intergenic_coverage", "-ci", default=10, type=float,
        help="The minimum coverage of intergenic sORF candidates. "
        "Default is 10.")
    sORF_add.add_argument(
        "--cutoff_antisense_coverage", "-ai", default=10, type=float,
        help="The minimum coverage of antisense sORF candidates. "
        "Default is 10.")
    sORF_add.add_argument(
        "--cutoff_5utr_coverage", "-cu5", default="p_0.5",
        help="The minimum coverage for 5'UTR derived sORF candidates. "
        "This value can be assigned by percentile (\"p\") or the amount of "
        "reads (\"n\"). For example, p_0.5 means that the coverage of sORF "
        "candidates should be higher than the 50 percentile of all 5'UTR "
        "transcripts. n_10 means that the coverage of sORF candidates "
        "should be higher than 10 reads. Default is p_0.5.")
    sORF_add.add_argument(
        "--cutoff_3utr_coverage", "-cu3", default="p_0.5",
        help="Similar to --cutoff_5utr_coverage. "
        "This value is for 3'UTRs. Default is p_0.5.")
    sORF_add.add_argument(
        "--cutoff_intercds_coverage", "-cuf", default="p_0.5",
        help="Similar to --cutoff_5utr_coverage. "
        "This value is for interCDS. Default is p_0.5.")
    sORF_add.add_argument(
        "--cutoff_base_coverage", "-cub", default=10, type=float,
        help="The general minimum coverage of all sORF candidates. "
        "All candidates should exceed this value. Default is 10.")
    sORF_add.add_argument(
        "--start_codon", "-ac", default=["ATG"], nargs='+',
        help="The types of start coden. If multiple types of start codon "
        "need to be assigned, please use spaces to separate them. "
        "Default is ATG.")
    sORF_add.add_argument(
        "--stop_codon", "-oc", default=["TAA", "TAG", "TGA"], nargs='+',
        help="The types of stop codon. If multiple types of stop codon need "
        "to be assigned, please use spaces to separate them. "
        "Default is TTA TAG TGA.")
    sORF_add.add_argument(
        "--min_rbs_distance", "-mr", default=3, type=int,
        help="The minimum distance (nucleotides) between the ribosome binding "
        "site (Shine-Dalgarno sequence) and the start codon. Default is 3.")
    sORF_add.add_argument(
        "--max_rbs_distance", "-Mr", default=15, type=int,
        help="The maximum distance (nucleotides) between the ribosome binding "
        "site (Shine-Dalgarno sequence) and the start codon. Default is 15.")
    sORF_add.add_argument(
        "--rbs_not_after_tss", "-at", default=False, action="store_true",
        help="Include the sORFs which are not associated with ribosome "
        "binding site to the high-confidence sORF list. Default is False.")
    sORF_add.add_argument(
        "--tolerance_rbs", "-tr", default=2, type=int,
        help="The number of nucleotides of ribosome binding sites allowed "
        "to be different from AGGAGG. Default is 2.")
    sORF_add.add_argument(
        "--tolerance_3end", "-t3", default=30, type=int,
        help="The number of nucleotides can be extended from the end of "
        "transcript for searching stop codon. Default is 30.")
    sORF_add.add_argument(
        "--tolerance_5end", "-t5", default=5, type=int,
        help="The number of nucleotides can be extended from the starting "
        "point of transcript for searching start codon. Default is 5.")
    sORF_add.add_argument(
        "--print_all_combination", "-pa", default=False, action="store_true",
        help="For printing all combinations of multiple start and stop codons. "
        "Default is False.")
    sORF_add.add_argument(
        "--best_no_srna", "-bs", default=False, action="store_true",
        help="Excluding the sORFs which overlap with sRNAs to highly "
        "confidence sORF list. Default is False.")
    sORF_add.add_argument(
        "--best_no_tss", "-bt", default=False, action="store_true",
        help="Excluding the sORFs which do not start with TSS to highly "
        "confidence sORF list. Default is False.")
    sORF_add.add_argument(
        "--ignore_hypothetical_protein", "-ih", default=False,
        help="For ignoring hypothetical protein in the genome "
        "annotation file. Default is False.")
    sORF_parser.set_defaults(func=run_sORF_detection)
    
    # Parameters of promoter detection
    promoter_parser = subparsers.add_parser(
        "promoter", help="Discover promoter motifs.")
    promoter_basic = promoter_parser.add_argument_group(
        'basic arguments')
    promoter_basic.add_argument(
        "--program", "-p", default="both", choices=['meme', 'glam2', 'both'],
        help="Please choose the program -- meme, glam2 or both. Default is both")
    promoter_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    promoter_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    promoter_basic.add_argument(
        "--tss_files", "-t", required=True, nargs='+',
        help="Paths of the TSS gff files.")
    promoter_basic.add_argument(
        "--motif_width", "-w", default=[50], nargs='+',
        help="Length of the motifs to detects. "
        "For a range insert \"-\" between "
        "two values. Moreover, if multiple lengths need to be assigned, "
        "please use spaces to separate them. For an example, 50 2-10 means "
        "that the lengths of motifs are 50 and within 2 to 10. "
        "The number should be less or equal than --nt_before_TSS. "
        "Default is 50.")
    promoter_basic.add_argument(
        "--num_motifs", "-n", default=10,
        help="The number of motifs. Default is 10.")
    promoter_basic.add_argument(
        "--nt_before_tss", "-b", default=50, type=int,
        help="The number of nucleotides upstream of the TSSs "
        "for promoter prediction. Default is 50.")
    promoter_add = promoter_parser.add_argument_group(
        'additional arguments')
    promoter_add.add_argument(
        "--meme_path", default="meme",
        help="path of MEME.")
    promoter_add.add_argument(
        "--glam2_path", default="glam2",
        help="path of GLAM2.")
    promoter_add.add_argument(
        "--e_value", "-e", default=0.05, type=int,
        help="The maximum e value for running MEME. Default is 0.05.")
    promoter_add.add_argument(
        "--end_run", "-er", default=10000,
        help="If the result of GLAM2 is not improved after running this "
        "number of iteration, it will be ended. Default is 10000.")
    promoter_add.add_argument(
        "--parallels", "-pl", default=None,
        help="The number of  parallel jobs.")
    promoter_add.add_argument(
        "--tss_source", "-s", default=True, action="store_false",
        help="The TSS gff files are generated from ANNOgesic or not. "
        "Default is True (from ANNOgesic)")
    promoter_add.add_argument(
        "--tex_libs", "-tl", default=None, nargs='+',
        help="If --tss_source is False, please assign the name of the "
        "TEX+/- library. The format is: wig_file_path:TEX+/-(tex or notex):"
        "condition_id(integer):replicate_id(alphabet):strand(+ or -). "
        "If multiple wig files need to be assigned, please use spaces to "
        "separate the wig files. For an example, "
        "$WIG_PATH_1:tex:1:a:+ $WIG_PATH_2:tex:1:a:-.")
    promoter_add.add_argument(
        "--annotation_files", "-g", default=None, nargs='+',
        help="If --tss_source is False, please assign the paths of the "
        "genome annotation gff files containing CDSs, tRNAs, rRNAs, etc.")
    promoter_add.add_argument(
        "--combine_all", "-c", default=False, action="store_true",
        help="Generate global promoter motifs across all reference "
        "sequences in the TSS files. Default is False.")
    promoter_parser.set_defaults(func=run_MEME)
    
    # Parameters of operon detection
    operon_parser = subparsers.add_parser(
        "operon", help="Detect operons and sub-operons.")
    operon_basic = operon_parser.add_argument_group(
        'basic arguments')
    operon_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    operon_basic.add_argument(
        "--tss_files", "-t", default=None, nargs='+',
        help="Paths of the TSS gff files.")
    operon_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs etc. Gff files of TSSs, transcripts, terminators, "
        "need to be separately assigned to --tss_files, "
        "--transcript_files, --terminator_files, respectively.")
    operon_basic.add_argument(
        "--transcript_files", "-a", required=True, nargs='+',
        help="Paths of the transcript gff files.")
    operon_basic.add_argument(
        "--terminator_files", "-e", default=None, nargs='+',
        help="Paths of terminator gff files.")
    operon_add = operon_parser.add_argument_group(
        'additional arguments')
    operon_add.add_argument(
        "--tss_tolerance", "-tt", default=5, type=int,
        help="The 5'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated TSSs. Default is 5.")
    operon_add.add_argument(
        "--terminator_tolerance", "-et", default=30, type=int,
        help="The 3'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated terminators. "
        "Default is 30.")
    operon_add.add_argument(
        "--min_length", "-l", default=20, type=int,
        help="The minimum length of operon. Default is 20.")
    operon_parser.set_defaults(func=run_operon)
    
    # Parameters of CircRNA detection
    circrna_parser = subparsers.add_parser(
        "circrna", help="Detect circular RNAs.")
    circrna_basic = circrna_parser.add_argument_group(
        'basic arguments')
    circrna_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    circrna_basic.add_argument(
        "--read_files", "-rp", default=None, nargs='+',
        help="Paths of read fasta or fastq files. ANNOgesic will map the reads via "
        "segemehl (with -S). Required format: $SET_NAME:$READ1,$READ2,... "
        "If multiple data sets need to be assigned, please separated them "
        "by spaces. The read files compressed by bz2 or gz files can be accepted as well. "
        "For using BAM files, please check --bam_files.")
    circrna_basic.add_argument(
        "--bam_files", "-b", default=None, nargs='+',
        help="Path of input BAM files. Required format: "
        "$SET_NAME:$BAM1,$BAM2,... . BAM files need to be generated using the "
        "mapper segemehl with the parameter \"-S\". If multiple data sets "
        "need to be assigned, please separated them by spaces. For using "
        "fasta files of reads, please check --read_files.")
    circrna_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files. ")
    circrna_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs, "
        "tRNAs, rRNAs, etc.")
    circrna_add = circrna_parser.add_argument_group(
        'additional arguments')
    circrna_add.add_argument(
        "--segemehl_path", default="segemehl.x",
        help="Path of segemehl.x in segemehl package.")
    circrna_add.add_argument(
        "--testrealign_path", default="testrealign.x",
        help="Path of testrealign.x in segemehl package.")
    circrna_add.add_argument(
        "--samtools_path", default="samtools",
        help="Path of samtools.")
    circrna_add.add_argument(
        "--parallels", "-p", default=10, type=int,
        help="The number of parallels runs for mapping if "
        "--read_files is assigned. Default is 10.")
    circrna_add.add_argument(
        "--support_reads", "-s", default=10, type=int,
        help="The minimum reads for supporting a circular RNA. Default is 10.")
    circrna_add.add_argument(
        "--start_ratio", "-sr", default=0.5, type=float, 
        help="The minimum ratio -- (reads of circRNA / all reads) "
        "at the starting points of candidates. Default is 0.5.")
    circrna_add.add_argument(
        "--end_ratio", "-er", default=0.5, type=float, 
        help="Similar to --start_ratio. "
        "This value is for the end points of candidates. Default is 0.5.")
    circrna_add.add_argument(
        "--ignore_hypothetical_protein", "-ih", default=False,
        help="For ignoring hypothetical protein in "
        "the genome annotation file. Default is False.")
    circrna_parser.set_defaults(func=run_circrna)
    # Parameters of GO term
    goterm_parser = subparsers.add_parser(
        "go_term", help="Extract GO terms from Uniprot.")
    goterm_basic = goterm_parser.add_argument_group(
        'basic arguments')
    goterm_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    goterm_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containg CDSs.")
    goterm_basic.add_argument(
        "--transcript_files", "-a", default=None, nargs='+',
        help="Paths of the transcript gff files which can "
        "be used to retrieve GO terms based on expressed CDS and all CDS.")
    goterm_basic.add_argument(
        "--uniprot_id", "-u", required=True,
        help="Path of the UniProt ID mapping database.")
    goterm_basic.add_argument(
        "--go_obo", "-go", required=True,
        help="Path of go.obo.")
    goterm_basic.add_argument(
        "--goslim_obo", "-gs", required=True,
        help="Path of goslim.obo.")
    goterm_parser.set_defaults(func=run_goterm)
    
    # Parameters of sRNA target prediction
    starget_parser = subparsers.add_parser(
        "srna_target", help="Detect sRNA-mRNA interactions.")
    starget_basic = starget_parser.add_argument_group(
        'basic arguments')
    starget_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    starget_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing CDSs.")
    starget_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    starget_basic.add_argument(
        "--srna_files", "-s", required=True, nargs='+',
        help="Paths of the sRNA gff files.")
    starget_basic.add_argument(
        "--query_srnas", "-q", default=["all"], nargs='+',
        help="The query sRNA. The input format is "
        "$GENOME:$START:$END:$STRAND. If multiple sRNAs need to be assigned, "
        "please use spaces to separate them. For example, "
        "NC_007795.1:200:534:+ NC_007795.1:6767:6900:-. If all sRNAs of the "
        "reference genome need to be used, please assign \"all\". "
        "Default is all.")
    starget_basic.add_argument(
        "--program", "-p", nargs='+', required=True,
        choices=['RNAplex', 'RNAup', 'IntaRNA'],
        help="The program for detecting sRNA-mRNA interaction. "
        "Please choose \"RNAplex\", \"RNAup\" or \"IntaRNA\". "
        "If multiple programs need to be executed, please use space to "
        "separate them.")
    starget_basic.add_argument(
        "--top", "-t", default=50, type=int,
        help="The ranking number of targets which will be included "
        "to final output. The ranking is based on the binding energy. "
        "Default is 50.")
    starget_add = starget_parser.add_argument_group(
        'additional arguments')
    starget_add.add_argument(
        "--rnaplfold_path", default="RNAplfold",
        help="Path of RNAplfold of the Vienna package.")
    starget_add.add_argument(
        "--rnaplex_path", default="RNAplex",
        help="Path of RNAplex of the Vienna package.")
    starget_add.add_argument(
        "--rnaup_path", default="RNAup",
        help="Path of RNAup of the Vienna package.")
    starget_add.add_argument(
        "--intarna_path", default="IntaRNA",
        help="Path of IntaRNA.")
    starget_add.add_argument(
        "--interaction_length", "-i", default=30, type=int,
        help="Maximum length of an interaction. Default is 30.")
    starget_add.add_argument(
        "--window_size_target_rnaplex", "-wt", default=240, type=int,
        help="The average of the pair probabilities over "
        "windows for mRNA target. It is only applied for \"RNAplex\". "
        "Default is 240.")
    starget_add.add_argument(
        "--span_target_rnaplex", "-st", default=160, type=int,
        help="The maximum allowed separation of a base pair to "
        "span for mRNA target. It is only applied for \"RNAplex\". "
        "Default is 160.")
    starget_add.add_argument(
        "--window_size_srna_rnaplfold", "-ws", default=30, type=int,
        help="Similar to --window_size_target, but for sRNA. "
        "Default is 30.")
    starget_add.add_argument(
        "--span_srna_rnaplfold", "-ss", default=30, type=int,
        help="Similar to --span_target, but for sRNA. Default is 30.")
    starget_add.add_argument(
        "--unstructured_region_rnaplex_target", "-ut", default=30, type=int,
        help="Calculate the mean probability of the unpaired region for "
        "mRNA target. It only works for \"RNAplex\". Default is 30.")
    starget_add.add_argument(
        "--unstructured_region_rnaplex_srna", "-us", default=30, type=int,
        help="Similar to --unstructured_region_rnaplex_target, "
        "but for sRNA. Default is 30.")
    starget_add.add_argument(
        "--unstructured_region_rnaup", "-uu", default=40, type=int,
        help="Compute the mean probability of unpaired region. "
        "It only works for \"RNAup\". Default is 40.")
    starget_add.add_argument(
        "--energy_threshold_rnaplex", "-e", default=-8, type=float,
        help="The minimum energy for a duplex. It only works for "
        "\"RNAplex\". Default is -8.")
    starget_add.add_argument(
        "--duplex_distance_rnaplex", "-d", default=20, type=int,
        help="Distance between target 3'ends of two consecutive duplexes. "
        "It works for \"RNAplex\". Default is 20.")
    starget_add.add_argument(
        "--parallels_rnaplex", "-pp", default=5, type=int,
        help="The number of parallel jobs for running RNAplex. "
        "Default is 5.")
    starget_add.add_argument(
        "--parallels_rnaup", "-pu", default=20, type=int,
        help="The number of parallel jobs for running RNAup. "
        "Default is 20.")
    starget_add.add_argument(
        "--parallels_intarna", "-pi", default=10, type=int,
        help="The number of parallel jobs for running IntaRNA. "
        "Default is 10.")
    starget_add.add_argument(
        "--continue_rnaup", "-cr", default=False, action="store_true",
        help="For running RNAup based on the "
        "previous intermediate results if the previous process stopped. "
        "Default is False.")
    starget_add.add_argument(
        "--slide_window_size_srna_intarna", "-sw", default=150, type=int,
        help="The silding window size of sRNA sequences. 0 will use the "
        "full sequence to execute IntaRNA. Default is 150.")
    starget_add.add_argument(
        "--max_loop_length_srna_intarna", "-ls", default=100, type=int,
        help="The maximal loop length of sRNA. If the value is assigned "
        "by 0, --slide_window_size_srna_intarna will be used for the "
        "maximal loop length of sRNA. Default is 100.")
    starget_add.add_argument(
        "--slide_window_size_target_intarna", "-tw", default=150, type=int,
        help="The silding window size of target sequences. 0 will use the "
        "full sequence to execute IntaRNA. Default is 150.")
    starget_add.add_argument(
        "--max_loop_length_target_intarna", "-lt", default=100, type=int,
        help="The maximal loop length of target. If the value is assigned "
        "by 0, --slide_window_size_target_intarna will be used for the "
        "maximal loop length of target. Default is 100.")
    starget_add.add_argument(
        "--mode_intarna", "-mi", default=None, choices=['H', 'E', 'M'],
        help="The prediction mode of IntaRNA. 'H' is heuristic, 'M' is "
        "exact with long computational time. 'E' is exact with long "
        "computational time and high memory.")
    starget_add.add_argument(
        "--potential_target_start", "-ps", default=200, type=int,
        help="Distance for the extraction of upstream nucleotides of "
        "--target_feature. Default is 200.")
    starget_add.add_argument(
        "--potential_target_end", "-pe", default=150, type=int,
        help="Distance for the extraction of downstream nucleotides of "
        "--target_feature. Default is 150.")
    starget_add.add_argument(
        "--target_feature", "-tf", default=["CDS"], nargs='+',
        help="The feature name of potential targets. "
        "If multiple features need to be assigned, please use spaces to "
        "separate them. For example, CDS exon. Default is CDS.")
    starget_parser.set_defaults(func=sRNA_target)
    
    # Parameters of SNP transcript
    snp_parser = subparsers.add_parser(
        "snp", help="Detect SNP/mutation and generate fasta file if mutations "
        "were found.")
    snp_basic = snp_parser.add_argument_group(
        'basic arguments')
    snp_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    snp_basic.add_argument(
        "--bam_type", "-t", required=True,
        choices=["related_genome", "reference_genome"],
        help="If the BAM files are produced by mapping to a related genome, "
        "please assign \"related_genome\". "
        "the mutations between the related genome and the refernce genome can be "
        "detected for generating sequence of the query "
        "genome. If the BAM files are produced by mapping to the "
        "reference genome, please assign \"reference_genome\". "
        "The mutations of reference genome can be detected.")
    snp_basic.add_argument(
        "--program", "-p", required=True, nargs='+',
        choices=["with_BAQ", "without_BAQ", "extend_BAQ"],
        help="The program for detecting SNPs: "
        "\"with_BAQ\", \"without_BAQ\", \"extend_BAQ\". "
        "Multi-programs can be executed at the same time (separated "
        "by spaces). For example, with_BAQ without_BAQ extend_BAQ.")
    snp_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    snp_basic.add_argument(
        "--bam_files", "-b", required=True, nargs='+',
        help="Path of input BAM files. Required format: "
        "$SET_NAME:$BAM1,$BAM2,... .  If multiple sets need to be assigned, "
        "please use spaces to separate them.")
    snp_add = snp_parser.add_argument_group(
        'additional arguments')
    snp_add.add_argument(
        "--samtools_path", default="samtools",
        help="Path of samtools.")
    snp_add.add_argument(
        "--bcftools_path", default="bcftools",
        help="Path of bcftools.")
    snp_add.add_argument(
        "--quality", "-q", default=40, type=int,
        help="The minimum quality score of a mutation. Default is 40.")
    snp_add.add_argument(
        "--read_depth_range", "-d", default="n_10,none",
        help="The range of read depth for a mutation. The format is "
        "$MIN,$MAX. It can be assigned by different types: "
        "1. real number (\"r\"), 2. times of the number of bam files "
        "(counted from --bam_files) (\"n\") "
        "or 3. times of the average read depth (\"a\"). For example, "
        "n_10,a_3 is assigned. If the average read depth is 70 and 4 "
        "bam files are provided, n_10 will be 40 "
        "and a_3 will be 140 (average read depth * 3). "
        "Based on the same example, if this value is r_10,a_3, the minimum "
        "read depth will become exact 10 reads. If \"none\" is assigned, "
        "read depth will not be considered. Default is n_10,none.")
    snp_add.add_argument(
        "--ploidy", "-pl", default="haploid",
        choices=["haploid", "diploid"],
        help="The query genome is haploid or diploid. Default is haploid.")
    snp_add.add_argument(
        "--rg_tag", "-R", default=False, action="store_true",
        help="For one BAM file which includes multiple samples "
        "(opposite of --ignore-RG in samtools). Default is False.")
    snp_add.add_argument(
        "--caller", "-c", default="m", choices=["c", "m"],
        help="The types of caller - consensus-caller or multiallelic-caller. "
        "For details, please check documentation of bcftools. \"c\" represents "
        "consensus-caller. \"m\" represents multiallelic-caller. "
        "Default is m.")
    snp_add.add_argument(
        "--dp4_cutoff", "-D", default="n_10,0.8",
        help="The cutoff of DP4. The format is "
        "$MIN_SNP_READS_NUMBER:$MIN_SNP_READS_RATIO. $MIN_SNP_READS_NUMBER "
        "can be assigned by three types: 1. fix number (\"r\"), 2. times of "
        "the number of bam files (counted from --bam_files) (\"n\") or 3. "
        "times of average read depth (\"a\"). For example, n_10,0.8. If the "
        "BAM files are 4, it means the minimum mapped reads of a SNP is 40 "
        "(4 * 10), and the minimum ratio of mapped read of a SNP "
        "(mapped reads of a SNP / total reads) is 0.8. Default is n_10,0.8.")
    snp_add.add_argument(
        "--indel_fraction", "-if", default="n_10,0.8",
        help="The minimum IDV and IMF which supports for insertion "
        "of deletion. The minimum IDV can be assigned by different types: "
        "1. fix number (\"r\"), 2. times of the number of bam files "
        "(assigned by --bam_files) (\"n\") or "
        "3. times of the average read depth (\"a\"). The input format is "
        "$MIN_IDF:$MIN_IMF. For example, The value is n_10,0.8 and 4 BAM "
        "files are assigned. The minimum IDV is 40, and the minimum IMF "
        "is 0.8. Default is n_10,0.8.")
    snp_add.add_argument(
        "--filter_tag_info", "-ft",
        default=["RPB_b0.1", "MQSB_b0.1", "MQB_b0.1", "BQB_b0.1"], nargs='+',
        help="For using more filters to improve the detection. "
        "Please assign 1. the name of a tag, 2. bigger (\"b\") or smaller "
        "(\"s\") and 3. the value of the filter. For example, "
        "\"RPB_b0.1,MQ0F_s0\" means that RPB should be bigger than 0.1 "
        "and MQ0F should be smaller than 0. "
        "Default is RPB_b0.1,MQSB_b0.1,MQB_b0.1,BQB_b0.1.")
    snp_parser.set_defaults(func=SNP)
    
    # Parameters of protein-protein interaction network
    ppi_parser = subparsers.add_parser(
        "ppi_network", help="Detect protein-protein interactions suported by "
        "literature.")
    ppi_basic = ppi_parser.add_argument_group(
        'basic arguments')
    ppi_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    ppi_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the genome annotation gff files containing genes and "
        "CDSs with proper locus_tag items in the attributes. ")
    ppi_basic.add_argument(
        "--species_string", "-d", required=True,
        help="Path of the species table of STRING "
        "(species.$VERSION.txt).")
    ppi_basic.add_argument(
        "--query_strains", "-s", nargs='+', required=True,
        help="The name of the input file of the query genomes. "
        "Required format: "
        "$GFF_FILE:$STRAIN_IN_GFF:$STRAIN_IN_STRING:$STRAIN_FOR_PUBMED. "
        "$GFF_FILE is the name of the gff file, $STRAIN_IN_GFF is the name/ID "
        "of the strain in the gff file, $STRAIN_IN_STRING is the strain name "
        "in species table of STRING (species.$VERSION.txt), and "
        "$STRAIN_FOR_PUBMED is the strain name for searching in Pubmed. "
        "If the strain is not available in STRING database, it can be relaced "
        "by a related strain. For example, Staphylococcus_aureus_HG003.gff:"
        "Staphylococcus_aureus_HG003:\"Staphylococcus aureus NCTC 8325\":"
        "\"Staphylococcus aureus\" (Staphylococcus aureus NCTC 8325 is a "
        "related strain of HG003 since HG003 is not available in STRING). "
        "If the assigned name is with spaces, "
        "please use double quotes. For assigning multiple strains, please "
        "separated them by spaces.")
    ppi_basic.add_argument(
        "--query", "-q", default=["all"], nargs='+',
        help="The query proteins. Required format: "
        "$GENOME_NAME_OF_GFF:$START_POINT:$END_POINT:$STRAND. For assigning "
        "multiple proteins, please use spaces to separate them. "
        "For example, Staphylococcus_aureus_HG003:345:456:+ "
        "Staphylococcus_aureus_HG003:2000:3211:-. For computing all proteins "
        "in gff files, just type \"all\". Default is all.")
    ppi_basic.add_argument(
        "--without_strain_pubmed", "-n", default=False, action="store_true",
        help="For retrieving the literature from Pubmed only based on protein "
        "name without assigning strains. Default is False.")
    ppi_add = ppi_parser.add_argument_group(
        'additional arguments')
    ppi_add.add_argument(
        "--score", "-ps", default=0.0, type=float,
        help="The minimum PIE score for searching literature. "
        "The value is from -1 (worst) to 1 (best). Default is 0.")
    ppi_add.add_argument(
        "--node_size", "-ns", default=4000, type=int,
        help="The size of nodes in figure, default is 4000.")
    ppi_parser.set_defaults(func=PPI)
    
    # Parameters of subcellular localization
    sub_local_parser = subparsers.add_parser(
        "localization", help="Predict subcellular localization "
        "of proteins.")
    sub_local_basic = sub_local_parser.add_argument_group(
        'basic arguments')
    sub_local_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    sub_local_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of genome annotation gff files containing CDSs.")
    sub_local_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of genome fasta files.")
    sub_local_basic.add_argument(
        "--transcript_files", "-a", default=None, nargs='+',
        help="Paths of the transcript gff files for detecting the "
        "subcellular localization based on expressed CDS and all CDS.")
    sub_local_basic.add_argument(
        "--bacteria_type", "-b", required=True,
        choices=["positive", "negative"],
        help="The type of bacteria (Gram-positive or Gram-negative). "
        "Please assign 'positive' or 'negative'.")
    sub_local_add = sub_local_parser.add_argument_group(
        'additional arguments')
    sub_local_add.add_argument(
        "--psortb_path", default="psort",
        help="Path of Psortb.")
    sub_local_add.add_argument(
        "--difference_multi", "-d", default=0.5,
        help="For the protein which have multiple predicted locations, "
        "if the difference of psorb scores is smaller than "
        "this value, all locations will be printed out. Default is 0.5. "
        "The maximum value is 10.")
    sub_local_parser.set_defaults(func=Sub_Local)
    
    # Parameters of riboswitch and RNA thermometer
    ribos_parser = subparsers.add_parser(
        "riboswitch_thermometer", help="Predict riboswitches and "
        "RNA thermometers.")
    ribos_basic = ribos_parser.add_argument_group(
        'basic arguments')
    ribos_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    ribos_basic.add_argument(
        "--program", "-p", default="both",
        choices=["riboswitch", "thermometer", "both"],
        help="Please choose the feature for the detection. "
        "The options can be \"riboswitch\", \"thermometer\", \"both\". "
        "Default is both.")
    ribos_basic.add_argument(
        "--riboswitch_id_file", "-ri",
        help="Path of the file which contains the information of riboswitches "
        "in Rfam. Required format of the file: "
        "$RFAM_ID{tab}$RIBOSWITCH_NAME{tab}$DESCRIPTION. "
        "Please check an example in "
        "https://github.com/Sung-Huan/ANNOgesic/blob/master/"
        "database/Rfam_riboswitch_ID.csv")
    ribos_basic.add_argument(
        "--rna_thermometer_id_file", "-ti",
        help="Same format as for -riboswitch_id_file, but for RNA "
        "thermometers. Please check an example in "
        "https://github.com/Sung-Huan/ANNOgesic/blob/master/"
        "database/Rfam_RNA_thermometer_ID.csv")
    ribos_basic.add_argument(
        "--annotation_files", "-g", required=True, nargs='+',
        help="Paths of the annotation gff files.")
    ribos_basic.add_argument(
        "--tss_files", "-t", default=None, nargs='+',
        help="Paths of the TSS gff files.")
    ribos_basic.add_argument(
        "--transcript_files", "-a", required=True, nargs='+',
        help="Paths of the transcript gff files.")
    ribos_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    ribos_basic.add_argument(
        "--rfam_path", "-R", required=True,
        help="Path of the Rfam CM database.")
    ribos_add = ribos_parser.add_argument_group(
        'additional arguments')
    ribos_add.add_argument(
        "--cmscan_path", "-cs", default="cmscan",
        help="Path of cmscan in Infernal package.")
    ribos_add.add_argument(
        "--cmpress_path", "-cp", default="cmpress",
        help="Path of cmpress in Infernal package.")
    ribos_add.add_argument(
        "--utr_length", "-u", default=300, type=int,
        help="The UTR length. Default is 300.")
    ribos_add.add_argument(
        "--cutoff", "-cf", default="e_0.001",
        help="The cutoff of the infernal search. The cutoff can be assigned "
        "by e value (assigned by 'e') or score (assigned by 's'). "
        "For example, 'e_0.001' represents using e value as a cutoff and the "
        "maximum value is 0.001. 's_8' represents using score as a cutoff and "
        "the minimum score is 8. Default is e_0.001.")
    ribos_add.add_argument(
        "--output_all", "-o", default=False, action="store_true",
        help="One query sequence may fit multiple riboswitches or "
        "RNA thermometers. It can print multiple riboswitches or "
        "RNA thermometers. Otherwise, only the highest confident one "
        "will be printed. Default is False.")
    ribos_add.add_argument(
        "--tolerance", "-to", default=10, type=int,
        help="The 5'ends and 3'ends of potential riboswitches or RNA "
        "thermometers will be extended by this value (nucleotides) for "
        "extracting the sequences to search in Rfam. Default is 10.")
    ribos_add.add_argument(
        "--tolerance_rbs", "-tr", default=2, type=int,
        help="The number of nucleotides of ribosome binding site "
        "allow to be different with AGGAGG. Default is 2.")
    ribos_parser.set_defaults(func=Ribos)
    
    # Parameters of CRISPR finder
    crispr_parser = subparsers.add_parser(
        "crispr", help="Predict CRISPR related RNAs.")
    crispr_basic = crispr_parser.add_argument_group(
        'basic arguments')
    crispr_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    crispr_basic.add_argument(
        "--fasta_files", "-f", required=True, nargs='+',
        help="Paths of the genome fasta files.")
    crispr_basic.add_argument(
        "--annotation_files", "-g", default=None, nargs='+',
        help="Paths of the genome gff files containing CDSs for comparing "
        "CRISPRs and the genome annotation to remove the false positives. "
        "Default is None.")
    crispr_add = crispr_parser.add_argument_group(
        'additional arguments')
    crispr_add.add_argument(
        "--crt_path", default="/usr/local/bin/CRT.jar",
        help="Path of CRT.jar. Default is /usr/local/bin/CRT.jar")
    crispr_add.add_argument(
        "--window_size", "-w", default=8, type=int,
        help="Length of the window size for searching CRISPR (range: 6-9). "
        "Default is 8.")
    crispr_add.add_argument(
        "--min_number_repeats", "-mn", default=3, type=int,
        help="Minimum number of repeats that a CRISPR must contain. "
        "Default is 3.")
    crispr_add.add_argument(
        "--min_length_repeat", "-ml", default=23, type=int,
        help="Minimum length of CRISPR repeats. Default is 23.")
    crispr_add.add_argument(
        "--Max_length_repeat", "-Ml", default=47, type=int,
        help="Maximum length of CRISPR repeats. Default is 47.")
    crispr_add.add_argument(
        "--min_length_spacer", "-ms", default=26, type=int,
        help="Minimum length of CRISPR spacers. Default is 26.")
    crispr_add.add_argument(
        "--Max_length_spacer", "-Ms", default=50, type=int,
        help="Maximum length of CRISPR spacers. Default is 50.")
    crispr_add.add_argument(
        "--ignore_hypothetical_protein", "-in", default=False,
        action="store_true",
        help="To ignore hypothetical proteins. "
        "Default is inactive.")
    crispr_parser.set_defaults(func=Crispr)
    
    # Parameters of merge features
    merge_parser = subparsers.add_parser(
        "merge_features", help="Merge all features to one gff file.")
    merge_basic = merge_parser.add_argument_group(
        'basic arguments')
    merge_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    merge_basic.add_argument(
        "--output_prefix", "-op", required=True,
        help="The prefix name of output gff file. "
        "The file name will be $OUTPUT_PREFIX_merge_features.gff.")
    merge_basic.add_argument(
        "--transcript_file", "-a", default=None,
        help="Path of transcript gff file. The parent transcripts "
        "(\"Parent\" in gff attributes) of all features will be generated.")
    merge_basic.add_argument(
        "--other_features_files", "-of", default=None, nargs='+',
        help="Paths of the gff files (besides transcript gff file). "
        "For assigning multiple gff files, please use spaces to separate them.")
    merge_basic.add_argument(
        "--source_for_overlapping", "-s", default=None, nargs='+',
        help="If the locations of features are overlapped, the module will "
        "only keep the feature positions which provided from "
        "--source_for_overlapping. For example, if --source_for_overlapping is "
        "'Ref' 'ANNOgesic', the overlapping features from these two sources "
        "will be both kept. However, if --source_for_overlapping is 'Ref', "
        "only the overlapping features from 'RefSeq will be kept.' The value "
        "of --source_for_overlapping should be the same as the second column "
        "of the gff file. if all the sources of the overlapping features "
        "can not be found in --source_for_overlapping, all the overlapping "
        "features will be kept. Default is keeping all overlap features.")
    merge_add = merge_parser.add_argument_group(
        'additional arguments')
    merge_add.add_argument(
        "--terminator_tolerance", "-et", default=30, type=int,
        help="The 3'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated terminators. "
        "Default is 30.")
    merge_add.add_argument(
        "--tss_tolerance", "-tt", default=5, type=int,
        help="The 5'ends of transcripts will be extended or withdrew by this "
        "value (nucleotides) for searching the associated TSSs. Default is 5.")
    merge_parser.set_defaults(func=Merge)

    # Parameters of generate screenshots
    screen_parser = subparsers.add_parser(
        "screenshot", help="Generate screenshots for selected features using IGV.")
    screen_basic = screen_parser.add_argument_group(
        'basic arguments')
    screen_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    screen_basic.add_argument(
        "--fasta_file", "-f", required=True,
        help="Path of the genome fasta file.")
    screen_basic.add_argument(
        "--main_gff", "-mg", required=True,
        help="Screenshots will be generated based on the positions of genomic "
        "features in this gff file.")
    screen_basic.add_argument(
        "--side_gffs", "-sg", default=None, nargs='+',
        help="The gff files of further genomic features (besides --main_gff). "
        "For assigning multiple files, please use spaces to separated them.")
    screen_basic.add_argument(
        "--tex_notex_libs", "-tl", default=None, nargs='+',
        help="TEX+/- wig files. The format is: "
        "wig_file_path:TEX+/-(tex or notex):condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files.")
    screen_basic.add_argument(
        "--frag_libs", "-fl", default=None, nargs='+',
        help="Wig files of RNA-Seq of fragmented transcripts. The format is: "
        "wig_file_path:frag:condition_id(integer):"
        "replicate_id(alphabet):strand(+ or -). If multiple wig files need "
        "to be assigned, please use spaces to separate the wig files. "
        "For example, my_lib_frag_forward.wig:frag:1:a:+ "
        "my_lib_frag_reverse.wig:frag:1:a:-.")
    screen_basic.add_argument(
        "--output_folder", "-o", required=True,
        help="Path of the output folder. A sub-folder "
        "\"screenshots\" in the --output_folder will be created to store "
        "the results.")
    screen_add = screen_parser.add_argument_group(
        'additional arguments')
    screen_add.add_argument(
        "--height", "-he", default=1500, type=int,
        help="Height of the screenshot. Default is 1500.")
    screen_add.add_argument(
        "--present", "-p", default="expand",
        choices=["expand", "collapse", "squish"],
        help="The presentation types (expand, collapse, or squish) of the "
        "features in the screenshot. Default is expand.")
    screen_parser.set_defaults(func=Screen)

    # Parameter of generating color png
    color_parser = subparsers.add_parser(
        "colorize_screenshot_tracks", help="Add color information to "
        "screenshots (e.g. useful for dRNA-Seq based TSS and PS "
        "detection. It only works after running \"screenshot\" "
        "(after running batch script).")
    color_basic = color_parser.add_argument_group(
        'basic arguments')
    color_basic.add_argument(
        "--project_path", "-pj", required=True,
        help="Path of the project folder.")
    color_basic.add_argument(
        "--screenshot_folder", "-f", required=True,
        help="The folder containing \"screenshots\" which are generated "
        "from the subcommand \"screenshot\".")
    color_basic.add_argument(
        "--track_number", "-t", type=int, required=True,
        help="The number of tracks.")
    color_add = color_parser.add_argument_group(
        'additional arguments')
    color_add.add_argument(
        "--imagemagick_covert_path", "-m", default="convert",
        help="Path of \"convert\" in the ImageMagick package.")
    color_parser.set_defaults(func=color_png)

    args = parser.parse_args()
    
    if args.version is True:
        print("ANNOgesic version " + __version__)
    elif "func" in dir(args):
        controller = Controller(args)
        args.func(controller)
    else:
        parser.print_help()



def create_project(controller):
    controller.create_project(__version__)

def get_input(controller):
    controller.get_input()

def get_target_fasta(controller):
    controller.get_target_fasta()

def run_RATT(controller):
    controller.ratt()

def run_expression(controller):
    controller.expression()

def multiparser(controller):
    controller.multiparser()

def run_TSSpredator(controller):
    controller.tsspredator()

def optimize_TSSpredator(controller):
    controller.optimize()

def color_png(controller):
    controller.color()

def run_Terminator(controller):
    controller.terminator()

def run_Transcript_Assembly(controller):
    controller.transcript()

def run_UTR_detection(controller):
    controller.utr_detection()

def run_sRNA_detection(controller):
    controller.srna_detection()

def run_sORF_detection(controller):
    controller.sorf_detection()

def run_MEME(controller):
    controller.meme()

def run_operon(controller):
    controller.operon()

def run_circrna(controller):
    controller.circrna()

def run_goterm(controller):
    controller.goterm()

def sRNA_target(controller):
    controller.srna_target()

def SNP(controller):
    controller.snp()

def PPI(controller):
    controller.ppi()

def Sub_Local(controller):
    controller.sublocal()

def Ribos(controller):
    controller.ribos()

def Crispr(controller):
    controller.crispr()

def Merge(controller):
    controller.merge()

def Screen(controller):
    controller.screen()
if __name__ == "__main__":
    main()
