Author: bugman
Date: Tue Oct 25 14:20:22 2011
New Revision: 14912
URL: http://svn.gna.org/viewcvs/relax?rev=14912&view=rev
Log:
Implemented the back end of the structure.displacement user function.
This will still require some debugging, and the data objects will need to be
stored in the XML save
files.
Modified:
1.3/generic_fns/structure/api_base.py
1.3/generic_fns/structure/internal.py
1.3/generic_fns/structure/main.py
Modified: 1.3/generic_fns/structure/api_base.py
URL:
http://svn.gna.org/viewcvs/relax/1.3/generic_fns/structure/api_base.py?rev=14912&r1=14911&r2=14912&view=diff
==============================================================================
--- 1.3/generic_fns/structure/api_base.py (original)
+++ 1.3/generic_fns/structure/api_base.py Tue Oct 25 14:20:22 2011
@@ -29,6 +29,9 @@
"""
# Python module imports.
+from math import pi
+from numpy import array, dot, float64, outer, transpose, zeros
+from numpy.linalg import norm, svd
from os import sep
from re import match
from types import MethodType
@@ -36,6 +39,7 @@
# relax module import.
from data.relax_xml import fill_object_contents, xml_to_object
+from maths_fns.rotation_matrix import R_to_axis_angle
from relax_errors import RelaxError, RelaxFileError,
RelaxFromXMLNotEmptyError, RelaxImplementError
from relax_io import file_root
from relax_warnings import RelaxWarning
@@ -181,6 +185,21 @@
raise RelaxImplementError
+ def calc_displacement(self, model_from=None, model_to=None,
atom_id=None):
+ """Calculate the rotational and translational displacement between
two structural models.
+
+ @keyword model_from: The optional model number for the
starting position of the displacement.
+ @type model_from: int or None
+ @keyword model_to: The optional model number for the ending
position of the displacement.
+ @type model_to: int or None
+ @keyword atom_id: The molecule, residue, and atom
identifier string. This matches the spin ID string format.
+ @type atom_id: str or None
+ """
+
+ # Raise the error.
+ raise RelaxImplementError
+
+
def delete(self):
"""Prototype method stub for deleting all structural data from the
current data pipe."""
@@ -797,6 +816,133 @@
+class Displacements:
+ """A special object for representing rotational and translational
displacements between models."""
+
+ def __init__(self):
+ """Initialise the storage objects."""
+
+ # The displacement structures.
+ self._translation_vector = {}
+ self._translation_distance = {}
+ self._rotation_matrix = {}
+ self._rotation_axis = {}
+ self._rotation_angle = {}
+
+
+ def _calc_centriod(self, coords):
+ """Calculate the centroid of the structure.
+
+ @keyword coord: The atomic coordinates.
+ @type coord: numpy rank-2, Nx3 array
+ @return: The centroid.
+ @rtype: numpy rank-1, 3D array
+ """
+
+ # The sum.
+ centroid = coords.sum(0) / coords.shape[0]
+
+ # Return.
+ return centroid
+
+
+ def _calc_rotation(self, coord_from=None, coord_to=None,
centroid_from=None, centroid_to=None):
+ """Calculate the rotation via SVD.
+
+ @keyword coord_from: The list of atomic coordinates for the
starting structure.
+ @type coord_from: numpy rank-2, Nx3 array
+ @keyword coord_to: The list of atomic coordinates for the
ending structure.
+ @type coord_to: numpy rank-2, Nx3 array
+ @keyword centroid_from: The starting centroid.
+ @type centroid_from: numpy rank-1, 3D array
+ @keyword centroid_to: The ending centroid.
+ @type centroid_to: numpy rank-1, 3D array
+ @return: The rotation matrix.
+ @rtype: numpy rank-2, 3D array
+ """
+
+ # Initialise the H matrix.
+ H = zeros((3, 3), float64)
+
+ # Loop over the atoms.
+ for i in range(coord_from.shape[0]):
+ # The positions shifted to the origin.
+ orig_from = coord_from[i] - centroid_from
+ orig_to = coord_to[i] - centroid_to
+
+ # The outer product.
+ H += outer(orig_from, orig_to)
+
+ # SVD.
+ U, S, V = svd(H)
+
+ # The rotation.
+ R = dot(V, transpose(U))
+
+ # Return the rotation.
+ return R
+
+
+ def _calculate(self, model_from=None, model_to=None, coord_from=None,
coord_to=None):
+ """Calculate the rotational and translational displacements using
the given coordinate sets.
+
+ @keyword model_from: The model number of the starting structure.
+ @type model_from: int
+ @keyword model_to: The model number of the ending structure.
+ @type model_to: int
+ @keyword coord_from: The list of atomic coordinates for the
starting structure.
+ @type coord_from: numpy rank-2, Nx3 array
+ @keyword coord_to: The list of atomic coordinates for the
ending structure.
+ @type coord_to: numpy rank-2, Nx3 array
+ """
+
+ # Calculate the centroids.
+ centroid_from = self._calc_centriod(coord_from)
+ centroid_to = self._calc_centriod(coord_to)
+
+ # The translation.
+ trans_vect = centroid_to - centroid_from
+ trans_dist = norm(trans_vect)
+
+ # Calculate the rotation.
+ R = self._calc_rotation(coord_from=coord_from, coord_to=coord_to,
centroid_from=centroid_from, centroid_to=centroid_to)
+ axis, angle = R_to_axis_angle(R)
+
+ # Print out.
+ print("\n\nCalculating the rotational and translational
displacements from model %s to model %s.\n" % (model_from, model_to))
+ print("Start centroid: [%20.15f, %20.15f, %20.15f]" %
(centroid_from[0], centroid_from[1], centroid_from[2]))
+ print("End centroid: [%20.15f, %20.15f, %20.15f]" %
(centroid_to[0], centroid_to[1], centroid_to[2]))
+ print("Translation vector: [%20.15f, %20.15f, %20.15f]" %
(trans_vect[0], trans_vect[1], trans_vect[2]))
+ print("Translation distance: %.15f" % trans_dist)
+ print("Rotation matrix:")
+ print(" [[%20.15f, %20.15f, %20.15f]" % (R[0, 0], R[0, 1], R[0,
2]))
+ print(" [%20.15f, %20.15f, %20.15f]" % (R[1, 0], R[1, 1], R[1,
2]))
+ print(" [%20.15f, %20.15f, %20.15f]]" % (R[2, 0], R[2, 1], R[2,
2]))
+ print("Rotation axis: [%20.15f, %20.15f, %20.15f]" %
(axis[0], axis[1], axis[2]))
+ print("Rotation angle (deg): %.15f" % (angle / 2.0 / pi * 360.0))
+
+ # Initialise structures if necessary.
+ if not self._translation_vector.has_key(model_from):
+ self._translation_vector[model_from] = {}
+ if not self._translation_distance.has_key(model_from):
+ self._translation_distance[model_from] = {}
+ if not self._rotation_matrix.has_key(model_from):
+ self._rotation_matrix[model_from] = {}
+ if not self._rotation_axis.has_key(model_from):
+ self._rotation_axis[model_from] = {}
+ if not self._rotation_angle.has_key(model_from):
+ self._rotation_angle[model_from] = {}
+
+ # Store the data.
+ self._translation_vector[model_from][model_to] = trans_vect
+ self._translation_distance[model_from][model_to] = trans_dist
+ self._rotation_matrix[model_from][model_to] = R
+ self._rotation_axis[model_from][model_to] = axis
+ self._rotation_angle[model_from][model_to] = angle
+
+
+
+
class ModelList(list):
"""List type data container for the different structural models.
Modified: 1.3/generic_fns/structure/internal.py
URL:
http://svn.gna.org/viewcvs/relax/1.3/generic_fns/structure/internal.py?rev=14912&r1=14911&r2=14912&view=diff
==============================================================================
--- 1.3/generic_fns/structure/internal.py (original)
+++ 1.3/generic_fns/structure/internal.py Tue Oct 25 14:20:22 2011
@@ -32,7 +32,7 @@
from warnings import warn
# relax module imports.
-from api_base import Base_struct_API, ModelList
+from api_base import Base_struct_API, ModelList, Displacements
from data.relax_xml import fill_object_contents, xml_to_object
from generic_fns import pipes, relax_re
from generic_fns.mol_res_spin import spin_loop
@@ -584,22 +584,20 @@
return data
- def atom_loop(self, atom_id=None, str_id=None, model_num_flag=False,
mol_name_flag=False, res_num_flag=False, res_name_flag=False,
atom_num_flag=False, atom_name_flag=False, element_flag=False,
pos_flag=False, ave=False):
+ def atom_loop(self, atom_id=None, str_id=None, model_num=None,
model_num_flag=False, mol_name_flag=False, res_num_flag=False,
res_name_flag=False, atom_num_flag=False, atom_name_flag=False,
element_flag=False, pos_flag=False, ave=False):
"""Generator function for looping over all atoms in the internal
relax structural object.
- @keyword atom_id: The molecule, residue, and atom
identifier string. Only atoms
- matching this selection will be yielded.
+ @keyword atom_id: The molecule, residue, and atom
identifier string. Only atoms matching this selection will be yielded.
@type atom_id: str
- @keyword str_id: The structure identifier. This can be
the file name, model
- number, or structure number. If None,
then all structures will
- be looped over.
+ @keyword str_id: The structure identifier. This can be
the file name, model number, or structure number. If None, then all
structures will be looped over.
@type str_id: str, int, or None
+ @keyword model_num: Only loop over a specific model.
+ @type model_num: int or None
@keyword model_num_flag: A flag which if True will cause the
model number to be yielded.
@type model_num_flag: bool
@keyword mol_name_flag: A flag which if True will cause the
molecule name to be yielded.
@type mol_name_flag: bool
- @keyword res_num_flag: A flag which if True will cause the
residue number to be
- yielded.
+ @keyword res_num_flag: A flag which if True will cause the
residue number to be yielded.
@type res_num_flag: bool
@keyword res_name_flag: A flag which if True will cause the
residue name to be yielded.
@type res_name_flag: bool
@@ -609,16 +607,12 @@
@type atom_name_flag: bool
@keyword element_flag: A flag which if True will cause the
element name to be yielded.
@type element_flag: bool
- @keyword pos_flag: A flag which if True will cause the
atomic position to be
- yielded.
+ @keyword pos_flag: A flag which if True will cause the
atomic position to be yielded.
@type pos_flag: bool
- @keyword ave: A flag which if True will result in this
method returning the
- average atom properties across all
loaded structures.
+ @keyword ave: A flag which if True will result in this
method returning the average atom properties across all loaded structures.
@type ave: bool
@return: A tuple of atomic information, as
described in the docstring.
- @rtype: tuple consisting of optional molecule
name (str), residue number
- (int), residue name (str), atom number
(int), atom name(str),
- element name (str), and atomic position
(array of len 3).
+ @rtype: tuple consisting of optional molecule
name (str), residue number (int), residue name (str), atom number (int), atom
name(str), element name (str), and atomic position (array of len 3).
"""
# Check that the structure is loaded.
@@ -629,7 +623,7 @@
sel_obj = Selection(atom_id)
# Model loop.
- for model in self.model_loop():
+ for model in self.model_loop(model_num):
# Loop over the molecules.
for mol_index in range(len(model.mol)):
mol = model.mol[mol_index]
@@ -785,6 +779,56 @@
# Return the data.
return data
+
+
+ def calc_displacement(self, model_from=None, model_to=None,
atom_id=None):
+ """Calculate the rotational and translational displacement between
two structural models.
+
+ @keyword model_from: The optional model number for the
starting position of the displacement.
+ @type model_from: int or None
+ @keyword model_to: The optional model number for the ending
position of the displacement.
+ @type model_to: int or None
+ @keyword atom_id: The molecule, residue, and atom
identifier string. This matches the spin ID string format.
+ @type atom_id: str or None
+ """
+
+ # Convert the model_from and model_to args to lists, is supplied.
+ if model_from != None:
+ model_from = [model_from]
+ if model_to != None:
+ model_to = [model_to]
+
+ # Create a list of all models.
+ models = []
+ for model in self.model_loop():
+ models.append(model.num)
+
+ # Set model_from or model_to to all models if None.
+ if model_from == None:
+ model_from = models
+ if model_to == None:
+ model_to = models
+
+ # Initialise the data structure.
+ if not hasattr(self, 'displacements'):
+ self.displacements = Displacements()
+
+ # Loop over the starting models.
+ for i in range(len(model_from)):
+ # Assemble the atomic coordinates.
+ coord_from = []
+ for pos in self.atom_loop(atom_id=atom_id,
model_num=model_from[i], pos_flag=True):
+ coord_from.append(pos[0])
+
+ # Loop over the ending models.
+ for j in range(len(model_to)):
+ # Assemble the atomic coordinates.
+ coord_to = []
+ for pos in self.atom_loop(atom_id=atom_id,
model_num=model_to[j], pos_flag=True):
+ coord_to.append(pos[0])
+
+ # Send to the base container for the calculations.
+ self.displacements._calculate(model_from=model_from[i],
model_to=model_to[j], coord_from=array(coord_from), coord_to=array(coord_to))
def delete(self):
Modified: 1.3/generic_fns/structure/main.py
URL:
http://svn.gna.org/viewcvs/relax/1.3/generic_fns/structure/main.py?rev=14912&r1=14911&r2=14912&view=diff
==============================================================================
--- 1.3/generic_fns/structure/main.py (original)
+++ 1.3/generic_fns/structure/main.py Tue Oct 25 14:20:22 2011
@@ -59,6 +59,22 @@
del spin.bond_vect
if hasattr(spin, 'xh_vect'):
del spin.xh_vect
+
+
+def displacement(model_from=None, model_to=None, atom_id=None):
+ """Calculate the rotational and translational displacement between two
structural models.
+
+ This will just redirect straight to the API.
+
+ @keyword model_from: The optional model number for the starting
position of the displacement.
+ @type model_from: int or None
+ @keyword model_to: The optional model number for the ending
position of the displacement.
+ @type model_to: int or None
+ @keyword atom_id: The molecule, residue, and atom identifier
string. This matches the spin ID string format.
+ @type atom_id: str or None
+ """
+
+ cdp.structure.calc_displacement(model_from=model_from,
model_to=model_to, atom_id=atom_id)
def get_pos(spin_id=None, str_id=None, ave_pos=False):
r14912 - in /1.3/generic_fns/structure: api_base.py internal.py main.py