Skip to content

Commit

Permalink
Fixed seek/tell bug for dysgu-call
Browse files Browse the repository at this point in the history
  • Loading branch information
kcleal committed Mar 25, 2024
1 parent 3faafbe commit 503345e
Show file tree
Hide file tree
Showing 4 changed files with 18 additions and 19 deletions.
11 changes: 6 additions & 5 deletions dysgu/call_component.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -930,7 +930,7 @@ cdef tuple break_ops(positions, precise, int limit, float median_pos):
ops.append((i + v + 1, -1, 1))
else:
ops.append((i + v, -1, 1))
ops = sorted(ops, reverse=limit < 0)
ops.sort(reverse=limit < 0)
cdef int cum_sum = 0
cdef int max_sum = 0
cdef int max_sum_i = 0
Expand Down Expand Up @@ -1386,10 +1386,10 @@ cdef call_from_reads(u_reads_info, v_reads_info, int insert_size, int insert_std
svtypes_counts[2].append(i)
elif i.svtype == "TRA":
svtypes_counts[3].append(i)
svtypes_res = sorted(svtypes_counts, key=sort_by_length, reverse=True)
svtypes_counts.sort(key=sort_by_length, reverse=True)
results = []
cdef EventResult_t er
for sub_informative in svtypes_res:
for sub_informative in svtypes_counts:
if len(sub_informative) >= min_support:
svtype = sub_informative[0].svtype
if svtype == "INV" or svtype == "TRA":
Expand Down Expand Up @@ -1617,7 +1617,8 @@ cdef get_reads(infile, nodes_info, buffered_reads, n2n, bint add_to_buffer, site
aligns.append((n, buffered_reads[int_node]))
continue
fpos.append((n, int_node))
for node, int_node in sorted(fpos, key=fpos_srt):
fpos.sort(key=fpos_srt)
for node, int_node in fpos:
infile.seek(node.tell)
try:
a = next(infile)
Expand All @@ -1643,7 +1644,7 @@ cdef get_reads(infile, nodes_info, buffered_reads, n2n, bint add_to_buffer, site
if add_to_buffer:
buffered_reads[int_node] = a
break
# else:
# else:
# if has_index:
# nn = n2n[int_node]
# hash_names[(nn.hash_name, nn.flag, nn.pos, nn.chrom)] = nn
Expand Down
4 changes: 2 additions & 2 deletions dysgu/cluster.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -395,8 +395,8 @@ def merge_events(potential, max_dist, tree, paired_end=False, try_rev=False, pic
node_to_event = {i.event_id: i for i in potential}
cdef int k
for grp in components:
c = [node_to_event[n] for n in grp]
best = sorted(c, key=srt_func, reverse=True)
best = [node_to_event[n] for n in grp]
best.sort(key=srt_func, reverse=True)
w0 = best[0]
if not pick_best:
weight = w0.pe + w0.supp + w0.spanning
Expand Down
20 changes: 9 additions & 11 deletions dysgu/coverage.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -161,10 +161,10 @@ cdef class GenomeScanner:
else:
f_iter = self.input_bam.fetch(until_eof=True)
for aln in f_iter:
# if aln.flag & 1284 or aln.mapq < mq_thresh or aln.cigartuples is None:
# continue

cigar_l = aln._delegate.core.n_cigar
if cigar_l == 0 or aln.flag & 1284:
tell = 0 if self.no_tell else self.input_bam.tell()
continue

if aln.rname != self.current_tid:
Expand All @@ -182,7 +182,6 @@ cdef class GenomeScanner:
good_read = True
elif not aln.flag & 2:
good_read = True
# cigar_l = aln._delegate.core.n_cigar
cigar_p = bam_get_cigar(aln._delegate)
for i in range(cigar_l):
cigar_value = cigar_p[i]
Expand All @@ -199,12 +198,11 @@ cdef class GenomeScanner:
elif opp == 0 or opp == 7 or opp == 8:
self.cpp_cov_track.add(pos + index_start, pos + index_start + length)
index_start += length
if aln.mapq < mq_thresh:
continue
if not good_read:
continue
if not self.cpp_cov_track.cov_val_good(self.current_tid, aln.rname, pos):

if aln.mapq < mq_thresh or not good_read or not self.cpp_cov_track.cov_val_good(self.current_tid, aln.rname, pos):
tell = 0 if self.no_tell else self.input_bam.tell()
continue

self._add_to_bin_buffer(aln, tell)
tell = 0 if self.no_tell else self.input_bam.tell()
while len(self.staged_reads) > 0:
Expand Down Expand Up @@ -460,9 +458,9 @@ cdef class GenomeScanner:
cdef float current_coverage
if self.current_chrom != rname:
self.current_chrom = rname
in_roi = False
if self.overlap_regions:
in_roi = intersecter(self.overlap_regions, a.rname, apos, apos+1)
# in_roi = False
# if self.overlap_regions:
# in_roi = intersecter(self.overlap_regions, a.rname, apos, apos+1)
if rname == self.current_chrom and bin_pos == self.current_pos:
self.current_bin.append((a, tell))
else:
Expand Down
2 changes: 1 addition & 1 deletion dysgu/graph.pyx
Original file line number Diff line number Diff line change
Expand Up @@ -704,7 +704,7 @@ class AlignmentsSA:
query_end = start_temp + query_end - query_start
query_start = start_temp
query_aligns.append(AlnBlock(query_start, query_end, ref_start, ref_end, gettid(sa[0]), int(sa[4]), sa[2], False))
query_aligns = sorted(query_aligns)
query_aligns.sort()
cdef int index = 0
for i, item in enumerate(query_aligns):
if item.this:
Expand Down

0 comments on commit 503345e

Please sign in to comment.