minimap2/python/minimap2.py at v2.18 · Triadai/minimap2

Repository navigation

#!/usr/bin/env python

import sys

import getopt

import mappy as mp

def main(argv):

opts, args = getopt.getopt(argv[1:], "x:n:m:k:w:r:c")

if len(args) < 2:

print("Usage: minimap2.py [options] <ref.fa>|<ref.mmi> <query.fq>")

print("Options:")

print(" -x STR preset: sr, map-pb, map-ont, asm5, asm10 or splice")

print(" -n INT mininum number of minimizers")

print(" -m INT mininum chaining score")

print(" -k INT k-mer length")

print(" -w INT minimizer window length")

print(" -r INT band width")

print(" -c output the cs tag")

sys.exit(1)

preset = min_cnt = min_sc = k = w = bw = None

out_cs = False

for opt, arg in opts:

if opt == '-x': preset = arg

elif opt == '-n': min_cnt = int(arg)

elif opt == '-m': min_chain_score = int(arg)

elif opt == '-r': bw = int(arg)

elif opt == '-k': k = int(arg)

elif opt == '-w': w = int(arg)

elif opt == '-c': out_cs = True

a = mp.Aligner(args[0], preset=preset, min_cnt=min_cnt, min_chain_score=min_sc, k=k, w=w, bw=bw)

if not a: raise Exception("ERROR: failed to load/build index file '{}'".format(args[0]))

for name, seq, qual in mp.fastx_read(args[1]): # read one sequence

for h in a.map(seq, cs=out_cs): # traverse hits

print('{}\t{}\t{}'.format(name, len(seq), h))

if __name__ == "__main__":

main(sys.argv)