Slightly updated routine for circular DNA in case beads overlap

mrdna/segmentmodel.py

@@ -359,6 +359,30 @@ given context, provide the position
             # nt_pos = self.get_nt_position(seg)
             # return seg.nt_pos_to_contour(nt_pos)
 
+    def combine(self, other):
+        assert(other in self.intrahelical_neighbors)
+        assert(self in other.intrahelical_neighbors)
+
+        self.intrahelical_neighbors.remove(other)
+        other.intrahelical_neighbors.remove(self)
+        for b in other.intrahelical_neighbors:
+            b.intrahelical_neighbors.remove(other)
+            b.intrahelical_neighbors.append(self)
+            self.intrahelical_neighbors.append(b)
+        for l in other.locations:
+            self.locations.append(l)
+            l.particle = self
+        
+        ## Remove bead
+        other.parent.children.remove(other)
+        if other in other.parent.beads:
+            other.parent.beads.remove(other)
+        if 'orientation_bead' in other.__dict__:
+            other.parent.children.remove(other.orientation_bead)
+
+        for b in list(other.parent.bonds):
+            if other in b[:2]: other.parent.bonds.remove(b)
+
     def update_position(self, contour_position):
         self.contour_position = contour_position
         self.position = self.parent.contour_to_position(contour_position)
@@ -2285,6 +2309,32 @@ class SegmentModel(ArbdModel):
                 for b in (b1,b2): assert( b is not None )
                 b1.make_intrahelical_neighbor(b2)
 
+        def _combine_zero_sep_beads():
+            skip_keys = []
+            for b1,b2 in [k for k,v in special_seps.items() if v == 0]:
+                if (b1,b2) in skip_keys: continue
+                del special_seps[(b1,b2)]
+                del special_seps[(b2,b1)]
+                skip_keys.append((b2,b1))
+
+                for b in b2.intrahelical_neighbors:
+                    if b is b1: continue
+                    if b is not b1 and b.parent is b2.parent:
+                        sep = np.abs(b2.contour_position - b.contour_position) * b.parent.num_nt
+                        special_seps[(b,b1)] = sep
+                        special_seps[(b1,b)] = sep
+                b1.combine(b2)
+
+                for k in list(special_seps.keys()):
+                    if b2 in k:
+                        if b2 == k[0]:
+                            newkey = (b1,k[1])
+                        else:
+                            newkey = (k[0],b2)
+                        special_seps[newkey] = special_seps[k]
+                        del special_seps[k]
+        _combine_zero_sep_beads()
+
         """ Reassign bead types """
         if self.DEBUG: print("Assigning bead types")
         beadtype_s = dict()
@@ -3267,8 +3317,8 @@ class SegmentModel(ArbdModel):
                     ntsJ1,ntsJ2 = [segJ.contour_to_nt_pos(B.address) for B in (B1,B2)]
                     ntsI = ntsI2-ntsI1+1
                     ntsJ = ntsJ2-ntsJ1+1
-                    assert( np.isclose( ntsI, int(round(ntsI)) ) )
-                    assert( np.isclose( ntsJ, int(round(ntsJ)) ) )
+                    # assert( np.isclose( ntsI, int(round(ntsI)) ) )
+                    # assert( np.isclose( ntsJ, int(round(ntsJ)) ) )
                     ntsI,ntsJ = [int(round(i)) for i in (ntsI,ntsJ)]
 
                     ## Find if dsDNA "segments" are pointing in same direction