How to use the mappy.revcomp function in mappy
To help you get started, we’ve selected a few mappy examples, based on popular ways it is used in public projects.
Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
haotianteng / Chiron / chiron / chiron_label.py View on Github 
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:
print("No Reference sequence found in %s"%(input_f))
return raw_signal,raw_seq,ref_seq,decap_eventfor i, h in enumerate(seqmaps):
if i > 0:
flag = int(FSECONDARY)
elif not h.is_primary:
flag = int(FSUPPLEMENTARY)
else:
flag = 0
leftclip = '{}S'.format(h.q_st) if h.q_st > 0 else ''
rightclip = '{}S'.format(len(seq) - h.q_en) if h.q_en < len(seq) else ''
if h.strand > 0:
seq_f = seq
qual_f = qual
else:
seq_f = mappy.revcomp(seq)
qual_f = qual[::-1]
leftclip, rightclip = rightclip, leftclip
flag |= FREVERSE
fullcigar = leftclip + h.cigar_str + rightclip
yield (name, flag, h.ctg, h.r_st + 1, h.mapq, fullcigar, '*',
0, 0, seq_f, qual_f, 'NM:i:{}'.format(h.NM))
