Merge pull request #52 from iqbal-lab-org/nucmer_block_bug_fix

Nucmer block bug fix

pling/align_snakemake/integerise_plasmids.py

@@ -77,7 +77,7 @@ def get_coverage(tree):
     return coverage
 
 def integerise_plasmids(plasmid_1: Path, plasmid_2: Path, prefix: str, plasmid_1_name, plasmid_2_name, identity_threshold=80, length_threshold=200):
-    subprocess.check_call(f"perl -w $(which dnadiff) {plasmid_1} {plasmid_2} -p {prefix} 2>/dev/null", shell=True)
+    subprocess.check_call(f"nucmer --diagdiff 20 --breaklen 500  --maxmatch -p {prefix} {plasmid_1} {plasmid_2} && delta-filter -1 {prefix}.delta > {prefix}.1delta", shell=True)
     show_coords_output = subprocess.check_output(f"show-coords -TrcldH -I {identity_threshold} {prefix}.1delta", shell=True).strip().split(b'\n')  # TODO: what about this threshold?
     show_snps_output = subprocess.check_output(f"show-snps -TrH {prefix}.1delta", shell=True).strip().split(b'\n')
 

pling/align_snakemake/matches.py

@@ -246,6 +246,8 @@ def resolve_overlaps(self, overlap_threshold):
                     self.sort(False)
                     i = self.list.index(current_match)
             i = i+1
+            if i>100000:
+                raise Exception("The resolve_overlaps function in matches.py has been stuck in the while loop for 100000 iterations! There is likely a bug, please raise an issue.")
         self.sort(True)
         i=0
         finished = False
@@ -254,7 +256,6 @@ def resolve_overlaps(self, overlap_threshold):
             try:
                 overlap = self[i].rend-self[i+1].rstart
                 not_null = self[i].rstart!=self[i].rend and self[i+1].rstart!=self[i+1].rend
-                not_the_same = self[i].rstart!=self[i+1].rstart and self[i].rend!=self[i+1].rend #ignore multicopy
                 containment = self[i+1].rend<self[i].rend
             except IndexError:
                 finished = True
@@ -273,6 +274,8 @@ def resolve_overlaps(self, overlap_threshold):
                     self.sort(True)
                     i = self.list.index(current_match)
             i = i+1
+            if i>100000:
+                raise Exception("The resolve_overlaps function in matches.py has been stuck in the while loop for 100000 iterations! There is likely a bug, please raise an issue.")
 
 testing = False
 if testing == True: