Skip to content

Commit 44d3697

Browse files
committed
fix incorrect trim point and trim overlap size too
1 parent af86243 commit 44d3697

File tree

1 file changed

+11
-4
lines changed

1 file changed

+11
-4
lines changed

src/scripts/fix_haplogaps.py

+11-4
Original file line numberDiff line numberDiff line change
@@ -161,7 +161,7 @@ def parse_current_alns(tip_support, alns):
161161
assert ">" + node not in node_cuts
162162
assert "<" + node not in node_cuts
163163
assert ">" + node not in node_trims
164-
print("S\t" + node + "_trim" + "\t" + node_seqs[node][node_trims[">" + node]:])
164+
print("S\t" + node + "_trim" + "\t" + node_seqs[node][node_trims["<" + node]:])
165165
if (">" + node) in node_cuts:
166166
assert ">" + node not in node_trims
167167
assert "<" + node not in node_trims
@@ -183,10 +183,17 @@ def parse_current_alns(tip_support, alns):
183183
for edge2 in edges[edge]:
184184
fromnode = edge
185185
tonode = edge2
186+
overlap=edge_overlaps[canon(edge, edge2)]
186187
if ">" + fromnode[1:] in node_trims or "<" + fromnode[1:] in node_trims:
187-
fromnode = fromnode + "_trim"
188+
if ">" + fromnode[1:] in node_trims and ">" + fromnode[1:] in canon(edge, edge2): overlap -= node_trims[">" + fromnode[1:]]
189+
if "<" + fromnode[1:] in node_trims and "<" + fromnode[1:] in canon(edge, edge2): overlap -= node_trims["<" + fromnode[1:]]
190+
191+
fromnode = fromnode + "_trim"
188192
if ">" + tonode[1:] in node_trims or "<" + tonode[1:] in node_trims:
189-
tonode = tonode + "_trim"
193+
if ">" + tonode[1:] in node_trims and ">" + tonode[1:] in canon(edge, edge2): overlap -= node_trims[">" + tonode[1:]]
194+
if "<" + tonode[1:] in node_trims and "<" + tonode[1:] in canon(edge, edge2): overlap -= node_trims["<" + tonode[1:]]
195+
196+
tonode = tonode + "_trim"
190197
if fromnode[0] == ">" and (edge in node_cuts or revnode(edge) in node_cuts):
191198
fromnode = fromnode + "_end"
192199
elif fromnode[0] == "<" and (edge in node_cuts or revnode(edge) in node_cuts):
@@ -195,7 +202,7 @@ def parse_current_alns(tip_support, alns):
195202
tonode = tonode + "_beg"
196203
elif tonode[0] == "<" and (edge2 in node_cuts or revnode(edge2) in node_cuts):
197204
tonode = tonode + "_end"
198-
print("L\t" + fromnode[1:] + "\t" + ("+" if fromnode[0] == ">" else "-") + "\t" + tonode[1:] + "\t" + ("+" if tonode[0] == ">" else "-") + "\t" + str(edge_overlaps[canon(edge, edge2)]) + "M")
205+
print("L\t" + fromnode[1:] + "\t" + ("+" if fromnode[0] == ">" else "-") + "\t" + tonode[1:] + "\t" + ("+" if tonode[0] == ">" else "-") + "\t" + str(overlap) + "M")
199206

200207
next_fake_node_id = 0
201208
for edge in extra_edges:

0 commit comments

Comments
 (0)