Source code for prody.atomic.dihedral

# -*- coding: utf-8 -*-
"""This module defines :class:`Dihedral` for dealing with dihedral information provided
by using :meth:`.AtomGroup.setDihedrals` method."""

from numbers import Integral
import numpy as np

RAD2DEG = 180 / np.pi

__all__ = ['Dihedral']


[docs]class Dihedral(object): """A pointer class for dihedrald atoms. Following built-in functions are customized for this class: * :func:`len` returns dihedral length, i.e. :meth:`getSize` * :func:`iter` yields :class:`~.Atom` instances""" __slots__ = ['_ag', '_acsi', '_indices'] def __init__(self, ag, indices, acsi=None): self._ag = ag self._indices = np.array(indices) if acsi is None: self._acsi = ag.getACSIndex() else: self._acsi = acsi def __repr__(self): one, two, three, four = self._indices names = self._ag._getNames() return '<Dihedral: {0}({1})--{2}({3})--{4}({5}--{6}({7})) from {8}>'.format( names[one], one, names[two], two, names[three], three, names[four], four, str(self._ag)) def __str__(self): one, two, three, four = self._indices names = self._ag._getNames() return '{0}({1})--{2}({3})--{4}({5})--{6}({7})'.format( names[one], one, names[two], two, names[three], three, names[four], four) def __eq__(self, other): return (isinstance(other, Dihedral) and other.getAtomGroup() is self._ag and (np.all(other.getIndices() == self._indices) or np.all(other.getIndices() == list(reversed(self._indices))))) def __ne__(self, other): return not self.__eq__(other) def __size__(self): return self.getSize() def __iter__(self): for index in self._indices: yield self._ag[index]
[docs] def getAtomGroup(self): """Returns atom group.""" return self._ag
[docs] def getAtoms(self): """Returns dihedrald atoms.""" return (self._ag[self._indices[0]], self._ag[self._indices[1]], self._ag[self._indices[2]])
[docs] def getIndices(self): """Returns indices of dihedrald atoms.""" return self._indices.copy()
[docs] def getSize(self, radian=False): """Returns dihedral size.""" a1, a2, a3 = self.getVectors() v1 = np.cross(a1, a2) v1 = v1 / (v1 * v1).sum(-1)**0.5 v2 = np.cross(a2, a3) v2 = v2 / (v2 * v2).sum(-1)**0.5 porm = np.sign((v1 * a3).sum(-1)) rad = np.arccos((v1*v2).sum(-1) / ((v1**2).sum(-1) * (v2**2).sum(-1))**0.5) if not porm == 0: rad = rad * porm if radian: return rad else: return rad * RAD2DEG
[docs] def getVectors(self): """Returns bond vectors that originate from the central atom.""" one, two, three, four = self._indices acsi = self.getACSIndex() vector1 = self._ag._coords[acsi, two] - self._ag._coords[acsi, one] vector2 = self._ag._coords[acsi, three] - self._ag._coords[acsi, two] vector3 = self._ag._coords[acsi, four] - self._ag._coords[acsi, three] return vector1, vector2, vector3
[docs] def getACSIndex(self): """Returns index of the coordinate set.""" acsi = self._acsi if acsi >= self._ag._n_csets: raise ValueError('{0} has fewer coordsets than assumed by {1}' .format(str(self._ag), str(self))) return acsi
[docs] def setACSIndex(self, index): """Set the coordinate set at *index* active.""" if self._ag._coords is None: raise AttributeError('coordinates are not set') if not isinstance(index, Integral): raise TypeError('index must be an integer') n_csets = self._ag._n_csets if n_csets <= index or n_csets < abs(index): raise IndexError('coordinate set index is out of range') if index < 0: index += n_csets self._acsi = index
def evalDihedrals(dihedrals, n_atoms): """Returns an array mapping atoms to their dihedrald neighbors and an array that stores number of dihedrals made by each atom.""" numdihedrals = np.bincount(dihedrals.reshape((dihedrals.shape[0] * 4))) dmap = np.zeros((n_atoms, numdihedrals.max(), 3), int) dmap.fill(-1) index = np.zeros(n_atoms, int) for dihedral in dihedrals: a, b, c, d = dihedral dmap[a, index[a]] = [b, c, d] dmap[b, index[b]] = [a, c, d] dmap[c, index[c]] = [a, b, d] dmap[d, index[d]] = [a, b, c] index[dihedral] += 1 return dmap, numdihedrals def trimDihedrals(dihedrals, indices): """Returns dihedrals between atoms at given indices.""" iset = set(indices) dihedrals = [dihedral for dihedral in dihedrals if dihedral[0] in iset and dihedral[1] in iset and dihedral[2] in iset and dihedral[3] in iset] if dihedrals: newindices = np.zeros(indices.max()+1, int) newindices[indices] = np.arange(len(indices)) return newindices[np.array(dihedrals)]