diff --git a/CanonicalNucleotideAtoms.py b/CanonicalNucleotideAtoms.py index ef96f434eabbdf042b880b1a6626ceb2a275c183..a608eb0b7732da7603fdbb52479d381d4743009c 100644 --- a/CanonicalNucleotideAtoms.py +++ b/CanonicalNucleotideAtoms.py @@ -17,8 +17,50 @@ def stringToIntTuples(string, tupleLen, offset): return ret +class Nucleotide(Group): + def __init__(self, children, orientation, sequence, resname, factory): + Group.__init__(self, children=children, orientation = orientation) + self.atoms_by_name = {a.name:a for a in self.children} + self._factory = factory + self.sequence = sequence + self.resname = resname + self._first_atomic_index = 0 + + def _get_atomic_index(self, name): + return self._factory._name_to_idx[name] + self._first_atomic_index + + def get_intrahelical_above(self): + """ Returns next nucleotide along 5'-to-3' in same strand-side of helix """ + return self.get_intrahelical_nucleotide(1) + + def get_intrahelical_below(self): + """ Returns last nucleotide along 5'-to-3' in same strand-side of helix """ + return self.get_intrahelical_nucleotide(-1) + + def get_intrahelical_nucleotide(self,offset=1): + """ Returns nucleotide in same strand-side of helix with some offset along 5'-to-3' + + TODO: Does not yet cross intrahelical connections + """ + + strand_piece = self.parent + is_fwd = strand_piece.is_fwd + lo,hi = sorted([strand_piece.start,strand_piece.end]) + + if is_fwd: + nt_idx = strand_piece.children.index(self)+lo + offset + else: + nt_idx = hi-strand_piece.children.index(self) - offset + + seg = strand_piece.segment + if nt_idx < 0 or nt_idx >= seg.num_nt: + return None + else: + return seg._get_atomic_nucleotide(nt_idx, is_fwd) + class CanonicalNucleotideFactory(Group): - DefaultOrientation = rotationAboutAxis([0,0,1], 180) + # DefaultOrientation = rotationAboutAxis([0,0,1], 180) + DefaultOrientation = rotationAboutAxis([1,0,0], 180) def __init__(self, prefix, seq): self.sequence = seq # TODO: used? self.resname = resnames[seq] @@ -50,20 +92,23 @@ class CanonicalNucleotideFactory(Group): atoms = self.children + # minIdx = np.min(self.indices) - 1 self.bonds = np.loadtxt("%s.bonds.txt" % prefix, dtype='i4') # - minIdx self.angles = np.loadtxt("%s.angles.txt" % prefix, dtype='i4') # - minIdx self.dihedrals = np.loadtxt("%s.dihedrals.txt" % prefix, dtype='i4')# - minIdx self.impropers = np.loadtxt("%s.impropers.txt" % prefix, dtype='i4')# - minIdx + self._name_to_idx = { atoms[i].name:i for i in range(len(atoms)) } + def generate(self): - nt = Group(orientation = CanonicalNucleotideFactory.DefaultOrientation) - nt.__dict__['sequence'] = self.sequence - nt.__dict__['resname'] = self.resname - atoms = nt.children = [copy.copy(c) for c in self.children] - for a in atoms: - a.parent = nt - nt.atoms_by_name = {a.name:a for a in nt.children} + atoms = [copy.copy(c) for c in self.children] + nt = Nucleotide( + children = atoms, + orientation = CanonicalNucleotideFactory.DefaultOrientation, + sequence = self.sequence, + resname = self.resname, + factory = self) b = self.bonds a = self.angles diff --git a/segmentmodel.py b/segmentmodel.py index 9e34112dd8bb816c4273ba81027e1296606f4d90..ad712c16214079f252f48cff84d8968e6559eb5f 100644 --- a/segmentmodel.py +++ b/segmentmodel.py @@ -482,10 +482,8 @@ class Segment(ConnectableElement, Group): # if self.ntAt5prime is None and self.ntAt3prime is None: key = key+"singlet" key = seq - if not is_fwd: - nt_dict = canonicalNtFwd - else: - nt_dict = canonicalNtRev + nt_dict = canonicalNtFwd if is_fwd else canonicalNtRev + atoms = nt_dict[ key ].generate() # TODO: clone? atoms.orientation = orientation.dot(atoms.orientation)