def kma(inputfile_1, out_path, databases, db_path_kma, min_cov=0.6,
            threshold=0.9, kma_path="cge/kma/kma", sample_name="",
            inputfile_2=None, kma_mrs=None, kma_gapopen=None,
            kma_gapextend=None, kma_penalty=None, kma_reward=None, kma_pm=None,
            kma_fpm=None, kma_memmode=False, kma_nanopore=False, debug=False,
            kma_add_args=None):
        """
           I expect that there will only be one hit pr gene, but if there are
           more, I assume that the sequence of the hits are the same in the res
           file and the aln file.
        """
        threshold = threshold * 100
        min_cov = min_cov * 100

        kma_results = dict()
        kma_results["excluded"] = dict()

        if(sample_name):
           sample_name = "_" + sample_name

        # Initiate output dicts.
        gene_align_sbjct = {}
        gene_align_query = {}
        gene_align_homo = {}

        for db in databases:
            kma_db = db_path_kma + "/" + db
            kma_outfile = out_path + "/kma_" + db + sample_name
            kma_cmd = ("%s -t_db %s -o %s -e 1.0" % (kma_path,
                       kma_db, kma_outfile))
            if(inputfile_2 is not None):
                kma_cmd += " -ipe " + inputfile_1 + " " + inputfile_2
            else:
                kma_cmd += " -i " + inputfile_1
            if(kma_mrs is not None):
                kma_cmd += " -mrs " + str(kma_mrs)
            if(kma_gapopen is not None):
                kma_cmd += " -gapopen " + str(kma_gapopen)
            if(kma_gapextend is not None):
                kma_cmd += " -gapextend " + str(kma_gapextend)
            if(kma_penalty is not None):
                kma_cmd += " -penalty " + str(kma_penalty)
            if(kma_reward is not None):
                kma_cmd += " -reward " + str(kma_reward)
            if(kma_pm is not None):
                kma_cmd += " -pm " + kma_pm
            if(kma_fpm is not None):
                kma_cmd += " -fpm " + kma_fpm
            if (kma_memmode):
                kma_cmd += " -mem_mode "
            if (kma_nanopore):
                kma_cmd += " -bcNano "
                kma_cmd += " -mp 20 "
            if (kma_add_args is not None):
                kma_cmd += " " + kma_add_args + " "

            # kma output files
            align_filename = kma_outfile + ".aln"
            res_filename = kma_outfile + ".res"

            # If .res file exists then skip mapping
            if os.path.isfile(res_filename) and os.access(res_filename, os.R_OK):
                print("Found " + res_filename + " skipping DB.")
            else:
                # Call KMA
                if(debug):
                    print("KMA cmd: " + kma_cmd)

                process = subprocess.Popen(kma_cmd, shell=True,
                                           stdout=subprocess.PIPE,
                                           stderr=subprocess.PIPE)
                out, err = process.communicate()

            kma_results[db] = 'No hit found'

            # Open res file
            try:
                res_file = open(res_filename, "r")
                header = res_file.readline()
            except IOError as error:
                sys.exit("Error: KMA did not run as expected.\n" +
                         "KMA finished with the following response:" +
                         "\n{}\n{}".format(out.decode("utf-8"),
                                           err.decode("utf-8")))

            for line in res_file:

                if kma_results[db] == 'No hit found':
                    kma_results[db] = dict()
                    # kma_results[db]["excluded"] = dict()
                    # continue

                data = [data.strip() for data in line.split("\t")]
                gene = data[0]

                sbjct_len = int(data[3])
                sbjct_ident = float(data[4])
                coverage = float(data[5])
                depth = float(data[-3])
                q_value = float(data[-2])
                p_value = float(data[-1])


                if gene not in kma_results[db]:
                    hit = gene
                else:
                    hit = gene + "_" + str(len(kma_results[db][gene]) + 1)

                exclude_reasons = []

                if(coverage < min_cov or sbjct_ident < threshold):
                    exclude_reasons.append(coverage)
                    exclude_reasons.append(sbjct_ident)

                if(exclude_reasons):
                    # kma_results[db]["excluded"][hit] = exclude_reasons
                    kma_results["excluded"][hit] = exclude_reasons

                kma_results[db][hit] = dict()
                kma_results[db][hit]['sbjct_length'] = sbjct_len
                kma_results[db][hit]["perc_coverage"] = coverage
                kma_results[db][hit]["sbjct_string"] = []
                kma_results[db][hit]["query_string"] = []
                kma_results[db][hit]["homo_string"] = []
                kma_results[db][hit]["sbjct_header"] = gene
                kma_results[db][hit]["perc_ident"] = sbjct_ident
                kma_results[db][hit]["query_start"] = "NA"
                kma_results[db][hit]["query_end"] = "NA"
                kma_results[db][hit]["contig_name"] = "NA"
                kma_results[db][hit]["HSP_length"] = ""
                kma_results[db][hit]["cal_score"] = q_value
                kma_results[db][hit]["depth"] = depth
                kma_results[db][hit]["p_value"] = p_value
            res_file.close()


            if kma_results[db] == 'No hit found':
                continue

            # Open align file
            with open(align_filename, "r") as align_file:
                hit_no = dict()
                gene = ""
                # Parse through alignments
                for line in align_file:
                    # Skip empty lines
                    if(not line.strip()):
                        continue

                    # Check when a new gene alignment start
                    if line.startswith("#"):
                        gene = line[1:].strip()

                        if gene not in hit_no:
                            hit_no[gene] = str(1)
                        else:
                            hit_no[gene] += str(int(hit_no[gene]) + 1)

                    else:
                        # Check if gene one of the user specified genes
                        if hit_no[gene] == '1':
                            hit = gene
                        else:
                            hit = gene + "_" + hit_no[gene]

                        if hit in kma_results[db]:
                            line_data = line.split("\t")[-1].strip()
                            if line.startswith("template"):
                                kma_results[db][hit]["sbjct_string"] += (
                                    [line_data])
                            elif line.startswith("query"):
                                kma_results[db][hit]["query_string"] += (
                                    [line_data])
                            else:
                                kma_results[db][hit]["homo_string"] += (
                                    [line_data])
                        else:
                            print(hit + " not in results: ", kma_results)

            # concatinate all sequences lists and find subject start
            # and subject end

            gene_align_sbjct[db] = {}
            gene_align_query[db] = {}
            gene_align_homo[db] = {}

            for hit in kma_results[db]:
                # if(hit == "excluded"):
                # continue
                align_sbjct = "".join(kma_results[db][hit]['sbjct_string'])
                align_query = "".join(kma_results[db][hit]['query_string'])
                align_homo = "".join(kma_results[db][hit]['homo_string'])

                # Extract only aligned sequences
                start = re.search("^-*(\w+)", align_query).start(1)
                end = re.search("\w+(-*)$", align_query).start(1)

                kma_results[db][hit]['sbjct_string'] = align_sbjct[start:end]
                kma_results[db][hit]['query_string'] = align_query[start:end]
                kma_results[db][hit]['homo_string'] = align_homo[start:end]

                # Save align start and stop positions relative to
                # subject sequence
                kma_results[db][hit]['sbjct_start'] = start + 1
                kma_results[db][hit]["sbjct_end"] = end + 1
                kma_results[db][hit]["HSP_length"] = end - start

                # Count gaps in the alignment
                kma_results[db][hit]["gaps"] = (
                    kma_results[db][hit]['sbjct_string'].count("-") +
                    kma_results[db][hit]['query_string'].count("-"))

                # Save sequences covering the entire subject sequence
                # in seperate variables
                gene_align_sbjct[db][hit] = align_sbjct
                gene_align_query[db][hit] = align_query
                gene_align_homo[db][hit] = align_homo

        return FinderResult(kma_results, gene_align_sbjct, gene_align_query,
                            gene_align_homo)
I expect that there will only be one hit pr gene, but if there are
           more, I assume that the sequence of the hits are the same in the res
           file and the aln file.