Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
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))
trans_start: Start position of the transcription(required in RNA mode).
alignment: If requrie alignment.
"""
with h5py.File(input_f,'r') as input_fh:
raw_entry = list(input_fh['/Raw/Reads'].values())[0]
raw_signal = raw_entry['Signal'].value
raw_seq = input_fh[BASECALL_ENTRY+'/BaseCalled_template/Fastq'].value
if mode !=0:
assert trans_start is not None
raw_signal,raw_seq,decap_event = _decap(input_fh,trans_start,raw_signal,raw_seq)
else:
decap_event = input_fh[BASECALL_ENTRY+'/BaseCalled_template/Events'].value
align = None
ref_seq = None
if alignment:
ref = mappy.Aligner(ref_f,preset = "map-ont",best_n = 5)
aligns = ref.map(raw_seq.split(b'\n')[1])
maxmapq = -np.inf
for aln in aligns:
if aln.mapq > maxmapq:
maxmapq = aln.mapq
align = aln
if align is None:
print("FAIL MAPPING "+input_f)
else:
if align.strand == -1:
ref_seq = mappy.revcomp(ref.seq(align.ctg,start = align.r_st,end = align.r_en))
else:
ref_seq = ref.seq(align.ctg,start = align.r_st,end = align.r_en)
if (mode == 1) or (mode == -1):
raw_signal = raw_signal[::-1]
if ref_seq is None and alignment:
import os
import io
import re
import sys
import queue
import traceback
import threading
# pip allows tombo install without correct version of mappy, so check here
try:
import mappy
if sys.version_info[0] > 2:
mappy.Aligner(os.path.devnull).seq('')
else:
mappy.Aligner(os.path.devnull).seq(b'')
except AttributeError:
th.error_message_and_exit('Tombo requires mappy version >= 2.10.')
# Future warning from cython in h5py
import warnings
warnings.simplefilter(action='ignore', category=FutureWarning)
import h5py
import numpy as np
np.seterr(all='raise')
import multiprocessing as mp
from tqdm import tqdm
from tqdm._utils import _term_move_up
from time import sleep
print(" -m INT mininum chaining score")
print(" -k INT k-mer length")
print(" -w INT minimizer window length")
print(" -r INT band width")
sys.exit(1)
preset, min_cnt, min_sc, k, w, bw = None, None, None, None, None, None
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)
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): # traverse hits
print('{}\t{}\t{}'.format(name, len(seq), h))
def create_reference_index(reference_panel):
## Parse the panel of reference sequences used to find the best match
panel_aligner = mp.Aligner(reference_panel, best_n = 1)
if not panel_aligner:
raise Exception("ERROR: failed to load/build index file '{}'".format(reference_panel))
reference_panel_names = []
for name, seq, qual in mp.fastx_read(reference_panel, read_comment=False):
reference_panel_names.append(name)
return (panel_aligner, reference_panel_names)
def __init__(self, indexfile, output, output_layout):
self.aligner = mappy.Aligner(indexfile)
if not self.aligner:
raise Exception('Could not open minimap2 index {}.'.format(indexfile))
self.writers = self.open_writers(indexfile, output, output_layout)