Author: bugman
Date: Fri Mar 15 11:21:12 2013
New Revision: 18828
URL: http://svn.gna.org/viewcvs/relax?rev=18828&view=rev
Log:
Fully implemented the structure.create_rotor_pdb user function.
For this, the generic_fns.structure.geometric.create_rotor_propellers()
function was created.
Modified:
trunk/generic_fns/structure/geometric.py
trunk/user_functions/structure.py
Modified: trunk/generic_fns/structure/geometric.py
URL:
http://svn.gna.org/viewcvs/relax/trunk/generic_fns/structure/geometric.py?rev=18828&r1=18827&r2=18828&view=diff
==============================================================================
--- trunk/generic_fns/structure/geometric.py (original)
+++ trunk/generic_fns/structure/geometric.py Fri Mar 15 11:21:12 2013
@@ -21,7 +21,7 @@
# Python module imports.
from math import cos, pi, sin
-from numpy import arccos, array, dot, eye, float64, transpose, zeros
+from numpy import arccos, array, cross, dot, eye, float64, transpose, zeros
from numpy.linalg import norm
from os import getcwd
from string import ascii_uppercase
@@ -34,7 +34,7 @@
from generic_fns.structure.internal import Internal
from generic_fns.structure.mass import centre_of_mass
from lib.geometry.lines import closest_point_ax
-from maths_fns.rotation_matrix import two_vect_to_R
+from maths_fns.rotation_matrix import axis_angle_to_R, two_vect_to_R
from relax_errors import RelaxError, RelaxNoPdbError, RelaxNoSequenceError,
RelaxNoTensorError, RelaxNoVectorsError
from relax_io import get_file_path, open_write_file
from relax_warnings import RelaxWarning
@@ -585,13 +585,15 @@
status.observers.result_file.notify()
-def create_rotor_pdb(file=None, dir=None, axis=None, axis_pt=True,
centre=None, span=2e-9, blade_length=5e-10, force=False, staggered=False):
+def create_rotor_pdb(file=None, dir=None, rotor_angle=None, axis=None,
axis_pt=True, centre=None, span=2e-9, blade_length=5e-10, force=False,
staggered=False):
"""Create a PDB representation of a rotor motional model.
@keyword file: The name of the PDB file to create.
@type file: str
@keyword dir: The name of the directory to place the PDB file
into.
@type dir: str
+ @keyword rotor_angle: The angle of the rotor motion in degrees.
+ @type rotor_angle: float
@keyword axis: The vector defining the rotor axis.
@type axis: numpy rank-1, 3D array
@keyword axis_pt: A point lying anywhere on the rotor axis. This
is used to define the position of the axis in 3D space.
@@ -608,10 +610,13 @@
@type staggered: bool
"""
- # Convert the arguments.
+ # Convert the arguments to numpy arrays, radians and Angstrom.
axis = array(axis, float64)
axis_pt = array(axis_pt, float64)
centre = array(centre, float64)
+ rotor_angle = rotor_angle / 360.0 * 2.0 * pi
+ span = span * 1e10
+ blade_length = blade_length * 1e10
# Normalise.
axis_norm = axis / norm(axis)
@@ -629,21 +634,24 @@
mol = structure.get_molecule('rotor')
# The central point.
- point = closest_point_ax(line_pt=axis_pt, axis=axis, point=centre)
- mol.atom_add(pdb_record='HETATM', atom_num=1, atom_name='CTR',
res_name='AX', res_num=1, pos=point, element='PT')
+ mid_point = closest_point_ax(line_pt=axis_pt, axis=axis, point=centre)
+ mol.atom_add(pdb_record='HETATM', atom_num=1, atom_name='CTR',
res_name='AX', res_num=1, pos=mid_point, element='PT')
# Centre of the propellers.
- prop1 = point + axis_norm * span * 1e10
- mol.atom_add(pdb_record='HETATM', atom_num=2, atom_name='PRP',
res_name='PR1', res_num=2, pos=prop1, element='SI')
- mol.atom_connect(index1=0, index2=1)
+ prop1 = mid_point + axis_norm * span
+ prop1_index = 1
+ mol.atom_add(pdb_record='HETATM', atom_num=2, atom_name='PRP',
res_name='PRC', res_num=2, pos=prop1, element='O')
+ mol.atom_connect(index1=0, index2=prop1_index)
# Centre of the propellers.
- prop2 = point - axis_norm * span * 1e10
- mol.atom_add(pdb_record='HETATM', atom_num=3, atom_name='PRP',
res_name='PR2', res_num=3, pos=prop2, element='SI')
- mol.atom_connect(index1=0, index2=2)
-
- # Create the PDB file.
- ######################
+ prop2 = mid_point - axis_norm * span
+ prop2_index = 2
+ mol.atom_add(pdb_record='HETATM', atom_num=3, atom_name='PRP',
res_name='PRC', res_num=3, pos=prop2, element='O')
+ mol.atom_connect(index1=0, index2=prop2_index)
+
+ # Create the rotor propellers.
+ create_rotor_propellers(mol=mol, rotor_angle=rotor_angle, centre=prop1,
axis=axis, blade_length=blade_length, staggered=staggered)
+ create_rotor_propellers(mol=mol, rotor_angle=rotor_angle, centre=prop2,
axis=-axis, blade_length=blade_length, staggered=staggered)
# Print out.
print("\nGenerating the PDB file.")
@@ -664,6 +672,100 @@
dir = getcwd()
cdp.result_files.append(['diff_tensor_pdb', 'Diffusion tensor PDB',
get_file_path(file, dir)])
status.observers.result_file.notify()
+
+
+def create_rotor_propellers(mol=None, rotor_angle=None, centre=None,
axis=None, blade_length=5.0, staggered=False):
+ """Create a PDB representation of a rotor motional model.
+
+ @keyword mol: The internal structural object molecule
container to add the atoms to.
+ @type mol: MolContainer instance
+ @keyword rotor_angle: The angle of the rotor motion in radians.
+ @type rotor_angle: float
+ @keyword centre: The central point of the propeller.
+ @type centre: numpy rank-1, 3D array
+ @keyword axis: The vector defining the rotor axis.
+ @type axis: numpy rank-1, 3D array
+ @keyword blade_length: The length of the representative rotor blades in
Angstrom.
+ @type blade_length: float
+ @keyword staggered: A flag which if True will cause the rotor blades
to be staggered. This is used to avoid blade overlap.
+ @type staggered: bool
+ """
+
+ # Init.
+ step_angle = 2.0 / 360.0 * 2.0 * pi
+ R = zeros((3, 3), float64)
+ res_num = mol.last_residue() + 1
+
+ # Blade vectors.
+ blades = zeros((4, 3), float64)
+ if abs(dot(axis, array([0, 0, 1], float64))) == 1.0: # Avoid failures
in artificial situations.
+ blades[0] = cross(axis, array([1, 0, 0], float64))
+ else:
+ blades[0] = cross(axis, array([0, 0, 1], float64))
+ blades[0] = blades[0] / norm(blades[0])
+ blades[1] = cross(axis, blades[0])
+ blades[1] = blades[1] / norm(blades[1])
+ blades[2] = -blades[0]
+ blades[3] = -blades[1]
+ print axis
+ print blades
+
+ # Create the 4 blades.
+ for i in range(len(blades)):
+ # Staggering.
+ if staggered and i % 2:
+ blade_origin = centre - axis * 2
+
+ # Non-staggered.
+ else:
+ blade_origin = centre
+
+ # Add an atom for the blage origin.
+ blade_origin_index = mol.atom_add(pdb_record='HETATM',
atom_name='BLO', res_name='PRB', res_num=res_num, pos=blade_origin,
element='O')
+
+ # The centre edge point of the blade.
+ mid_point = blade_origin + blades[i] * blade_length
+ mid_pt_index = mol.atom_add(pdb_record='HETATM', atom_name='BLD',
res_name='PRB', res_num=res_num, pos=mid_point, element='N')
+ mol.atom_connect(index1=mid_pt_index, index2=blade_origin_index)
+
+ # Build the blade.
+ angle = 0.0
+ pos_last_index = mid_pt_index
+ neg_last_index = mid_pt_index
+ while True:
+ # Increase the angle.
+ angle += step_angle
+
+ # The edge rotation.
+ if angle > rotor_angle:
+ axis_angle_to_R(axis, rotor_angle, R)
+
+ # The normal rotation matrix.
+ else:
+ axis_angle_to_R(axis, angle, R)
+
+ # The positive edge.
+ pos_point = dot(R, mid_point - blade_origin) + blade_origin
+ pos_index = mol.atom_add(pdb_record='HETATM', atom_name='BLD',
res_name='PRB', res_num=res_num, pos=pos_point, element='N')
+ mol.atom_connect(index1=pos_index, index2=pos_last_index)
+ mol.atom_connect(index1=pos_index, index2=blade_origin_index)
+
+ # The negative edge.
+ neg_point = dot(transpose(R), mid_point - blade_origin) +
blade_origin
+ neg_index = mol.atom_add(pdb_record='HETATM', atom_name='BLD',
res_name='PRB', res_num=res_num, pos=neg_point, element='N')
+ mol.atom_connect(index1=neg_index, index2=neg_last_index)
+ mol.atom_connect(index1=neg_index, index2=blade_origin_index)
+
+ # Update the indices.
+ pos_last_index = pos_index
+ neg_last_index = neg_index
+
+ # Finish.
+ if angle > rotor_angle:
+ break
+
+ # Increment the residue number.
+ res_num += 1
def create_vector_dist(length=None, symmetry=True, file=None, dir=None,
force=False):
Modified: trunk/user_functions/structure.py
URL:
http://svn.gna.org/viewcvs/relax/trunk/user_functions/structure.py?rev=18828&r1=18827&r2=18828&view=diff
==============================================================================
--- trunk/user_functions/structure.py (original)
+++ trunk/user_functions/structure.py Fri Mar 15 11:21:12 2013
@@ -255,6 +255,13 @@
can_be_none = True
)
uf.add_keyarg(
+ name = "rotor_angle",
+ default = 0.0,
+ py_type = "float",
+ desc_short = "rotor angle",
+ desc = "The angle of the rotor motion in degrees."
+)
+uf.add_keyarg(
name = "axis",
py_type = "float_array",
dim = 3,
@@ -306,6 +313,10 @@
# Description.
uf.desc.append(Desc_container())
uf.desc[-1].add_paragraph("This creates a PDB file representation of a rotor
motional model. The model axis is defined by a vector and a single point on
the axis. The centre of the representation will be taken as the point on the
rotor axis closest to the given centre position. The size of the
representation is defined by the span, which is the distance from the central
point to the rotors, and the length of the blades.")
+# Prompt examples.
+uf.desc.append(Desc_container("Prompt examples"))
+uf.desc[-1].add_paragraph("The following is a synthetic example:")
+uf.desc[-1].add_prompt("relax> structure.create_rotor_pdb(file='rotor.pdb',
rotor_angle=20.0, axis=[0., 0., 1.], axis_pt=[1., 1., 0.], centre=[0., 0.,
2.], span=2e-9, blade_length=1e-9)")
uf.backend = generic_fns.structure.geometric.create_rotor_pdb
uf.menu_text = "create_&rotor_pdb"
uf.gui_icon = "oxygen.actions.list-add-relax-blue"
r18828 - in /trunk: generic_fns/structure/geometric.py user_functions/structure.py