Added _set_cadnano2_sequence method to SegmentModel

mrdna/readers/cadnano_segments.py

@@ -362,15 +362,20 @@ class cadnano_part(SegmentModel):
                 else:
                     raise Exception("Segment could not be found")
 
-                def callback(segment, bp_to_zidx=bp_to_zidx):
+                seg._cadnano_helix = hid
+                seg._cadnano_start = zid1
+                seg._cadnano_end   = zid2
+                seg._cadnano_bp_to_zidx = bp_to_zidx
+
+                def callback(segment):
                     for b in segment.beads:
                         bp = int(round(b.get_nt_position(segment)))
                         if bp < 0: bp = 0
                         if bp >= segment.num_nt: bp = segment.num_nt-1
                         try:
-                            b.beta = bp_to_zidx[bp]
+                            b.beta = segment._cadnano_bp_to_zidx[bp]
                             if 'orientation_bead' in b.__dict__:
-                                b.orientation_bead.beta = bp_to_zidx[bp]
+                                b.orientation_bead.beta = segment._cadnano_bp_to_zidx[bp]
                         except:
                             pass
                 seg._generate_bead_callbacks.append(callback)

mrdna/segmentmodel.py

@@ -17,6 +17,9 @@ from .model.CanonicalNucleotideAtoms import enmTemplateHC, enmTemplateSQ, enmCor
 
 from .model.spring_from_lp import k_angle as angle_spring_from_lp
 
+import csv
+
+
 # import pdb
 """
 TODO:
@@ -2791,6 +2794,35 @@ class SegmentModel(ArbdModel):
         for s in self.segments:
             s.randomize_unset_sequence()
 
+    def _set_cadnano2_sequence(self, cadnano_sequence_csv_file, replace_unset=None):
+        try:
+            self.segments[0]._cadnano_helix
+        except:
+            raise Exception("Cannot set sequence of a non-cadnano model using a cadnano sequence file")
+
+        starts = dict()
+        for strand in self.strands[1:]:
+            if strand.is_circular: continue
+            s = strand.strand_segments[0]
+            hid = s.segment._cadnano_helix
+            idx = s.segment._cadnano_bp_to_zidx[int(s.start)]
+            starts["{:d}[{:d}]".format(hid,idx)] = strand
+
+        base_set = tuple('ATCG')
+        with open(cadnano_sequence_csv_file) as fh:
+            reader = csv.reader(fh)
+            for values in reader:
+                if values[0] == "Start": continue
+                start,end,seq = values[:3]
+                strand = starts[start]
+                if replace_unset in base_set:
+                    seq = [replace_unset if s == '?' else s for s in seq]
+                else:
+                    seq = [random.choice(base_set) if s == '?' else s for s in seq]
+                if len([s for s in seq if s not in base_set]) > 0:
+                    print(seq)
+                strand.set_sequence(seq)
+
     def _generate_strands(self):
         ## clear strands
         try: