lib.structure.represent.rotor

1 ############################################################################### 2 # # 3 # Copyright (C) 2003-2014 Edward d'Auvergne # 4 # # 5 # This file is part of the program relax (http://www.nmr-relax.com). # 6 # # 7 # This program is free software: you can redistribute it and/or modify # 8 # it under the terms of the GNU General Public License as published by # 9 # the Free Software Foundation, either version 3 of the License, or # 10 # (at your option) any later version. # 11 # # 12 # This program is distributed in the hope that it will be useful, # 13 # but WITHOUT ANY WARRANTY; without even the implied warranty of # 14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # 15 # GNU General Public License for more details. # 16 # # 17 # You should have received a copy of the GNU General Public License # 18 # along with this program. If not, see <http://www.gnu.org/licenses/>. # 19 # # 20 ############################################################################### 21 22 # Module docstring. 23 """Representations of a mechanical rotor.""" 24 25 # Python module imports. 26 from math import pi 27 from numpy import array, cross, dot, float64, transpose, zeros 28 from numpy.linalg import norm 29 30 # relax module imports. 31 from lib.geometry.lines import closest_point_ax 32 from lib.geometry.rotations import axis_angle_to_R 33 34

35 -def rotor(structure=None, rotor_angle=None, axis=None, axis_pt=True, label=None, centre=None, span=2e-9, blade_length=5e-10, model_num=None, staggered=False, half_rotor=False):

36 """Create a PDB representation of a rotor motional model. 37 38 @keyword structure: The internal structural object instance to add the rotor to as a molecule. 39 @type structure: lib.structure.internal.object.Internal instance 40 @keyword rotor_angle: The angle of the rotor motion in radian. 41 @type rotor_angle: float 42 @keyword axis: The vector defining the rotor axis. 43 @type axis: numpy rank-1, 3D array 44 @keyword axis_pt: A point lying anywhere on the rotor axis. This is used to define the position of the axis in 3D space. 45 @type axis_pt: numpy rank-1, 3D array 46 @keyword label: The optional label for the rotor axis. If supplied, this cannot be longer than 4 characters due to the PDB format restriction. 47 @type label: str 48 @keyword centre: The central point of the representation. If this point is not on the rotor axis, then the closest point on the axis will be used for the centre. 49 @type centre: numpy rank-1, 3D array 50 @keyword span: The distance from the central point to the rotor blades (meters). 51 @type span: float 52 @keyword blade_length: The length of the representative rotor blades. 53 @type blade_length: float 54 @keyword model_num: The structural model number to add the rotor to. If not supplied, the same rotor structure will be added to all models. 55 @type model_num: int or None 56 @keyword staggered: A flag which if True will cause the rotor blades to be staggered. This is used to avoid blade overlap. 57 @type staggered: bool 58 @keyword half_rotor: A flag which if True will cause only the positive half of the rotor to be created. 59 @type half_rotor: bool 60 """ 61 62 # Convert the arguments to numpy arrays, radians and Angstrom. 63 axis = array(axis, float64) 64 axis_pt = array(axis_pt, float64) 65 centre = array(centre, float64) 66 span = span * 1e10 67 blade_length = blade_length * 1e10 68 69 # Normalise. 70 axis_norm = axis / norm(axis) 71 72 # Add a structure (handling up to 3 rotors). 73 if structure.has_molecule(model_num=model_num, name='rotor') and structure.has_molecule(model_num=model_num, name='rotor2'): 74 structure.add_molecule(model_num=model_num, name='rotor3') 75 mol_name = 'rotor3' 76 elif structure.has_molecule(model_num=model_num, name='rotor'): 77 structure.add_molecule(model_num=model_num, name='rotor2') 78 mol_name = 'rotor2' 79 else: 80 structure.add_molecule(model_num=model_num, name='rotor') 81 mol_name = 'rotor' 82 83 # Loop over the models. 84 for model in structure.model_loop(model_num): 85 # Alias the single molecule from the single model. 86 mol = structure.get_molecule(mol_name, model=model.num) 87 88 # The central point. 89 mid_point = closest_point_ax(line_pt=axis_pt, axis=axis, point=centre) 90 pos_index = mol.atom_add(pdb_record='HETATM', atom_name='CTR', res_name='RTX', res_num=1, pos=mid_point, element='PT') 91 92 # Centre of the propellers, positive half. 93 prop1 = mid_point + axis_norm * span 94 prop1_index = pos_index + 1 95 mol.atom_add(pdb_record='HETATM', atom_name='PRP', res_name='RTX', res_num=2, pos=prop1, element='O') 96 mol.atom_connect(index1=pos_index, index2=prop1_index) 97 98 # Centre of the propellers, negative half. 99 if not half_rotor: 100 prop2 = mid_point - axis_norm * span 101 prop2_index = pos_index + 2 102 mol.atom_add(pdb_record='HETATM', atom_name='PRP', res_name='RTX', res_num=3, pos=prop2, element='O') 103 mol.atom_connect(index1=pos_index, index2=prop2_index) 104 105 # Create the rotor propellers. 106 rotor_propellers(mol=mol, rotor_angle=rotor_angle, centre=prop1, axis=axis, blade_length=blade_length, staggered=staggered) 107 if not half_rotor: 108 rotor_propellers(mol=mol, rotor_angle=rotor_angle, centre=prop2, axis=-axis, blade_length=blade_length, staggered=staggered) 109 110 # Add atoms for the labels. 111 res_num = mol.res_num[-1]+1 112 label_pos1 = mid_point + axis_norm * (span + 2.0) 113 mol.atom_add(pdb_record='HETATM', atom_name=label, res_name='RTL', res_num=res_num, pos=label_pos1, element='H') 114 if not half_rotor: 115 label_pos2 = mid_point - axis_norm * (span + 2.0) 116 mol.atom_add(pdb_record='HETATM', atom_name=label, res_name='RTL', res_num=res_num, pos=label_pos2, element='H')

117 118

119 -def rotor_propellers(mol=None, rotor_angle=None, centre=None, axis=None, blade_length=5.0, step_angle=2.0, staggered=False):

120 """Create a PDB representation of a rotor motional model. 121 122 This function will create a fixed number of atoms, placing the propeller blade steps outside of the rotor angle on the edge. This is to allow for model support whereby the rotor angle between models can be different but the atomic count and connectivity must be the same in all models. 123 124 125 @keyword mol: The internal structural object molecule container to add the atoms to. 126 @type mol: MolContainer instance 127 @keyword rotor_angle: The angle of the rotor motion in radians. 128 @type rotor_angle: float 129 @keyword centre: The central point of the propeller. 130 @type centre: numpy rank-1, 3D array 131 @keyword axis: The vector defining the rotor axis. 132 @type axis: numpy rank-1, 3D array 133 @keyword blade_length: The length of the representative rotor blades in Angstrom. 134 @type blade_length: float 135 @keyword step_angle: The angle between propeller blades, in degrees. 136 @type step_angle: float 137 @keyword staggered: A flag which if True will cause the rotor blades to be staggered. This is used to avoid blade overlap. 138 @type staggered: bool 139 """ 140 141 # Init. 142 step_angle = step_angle / 360.0 * 2.0 * pi 143 step_num = int(pi / step_angle) 144 R = zeros((3, 3), float64) 145 res_num = mol.last_residue() + 1 146 147 # Blade vectors. 148 blades = zeros((4, 3), float64) 149 if abs(dot(axis, array([0, 0, 1], float64))) == 1.0: # Avoid failures in artificial situations. 150 blades[0] = cross(axis, array([1, 0, 0], float64)) 151 else: 152 blades[0] = cross(axis, array([0, 0, 1], float64)) 153 blades[0] = blades[0] / norm(blades[0]) 154 blades[1] = cross(axis, blades[0]) 155 blades[1] = blades[1] / norm(blades[1]) 156 blades[2] = -blades[0] 157 blades[3] = -blades[1] 158 159 # Create the 4 blades. 160 for i in range(len(blades)): 161 # Staggering. 162 if staggered and i % 2: 163 blade_origin = centre - axis * 2 164 165 # Non-staggered. 166 else: 167 blade_origin = centre 168 169 # Add an atom for the blage origin. 170 blade_origin_index = mol.atom_add(pdb_record='HETATM', atom_name='BLO', res_name='RTB', res_num=res_num, pos=blade_origin, element='O') 171 172 # The centre edge point of the blade. 173 mid_point = blade_origin + blades[i] * blade_length 174 mid_pt_index = mol.atom_add(pdb_record='HETATM', atom_name='BLD', res_name='RTB', res_num=res_num, pos=mid_point, element='N') 175 mol.atom_connect(index1=mid_pt_index, index2=blade_origin_index) 176 177 # Build the blade. 178 angle = 0.0 179 pos_last_index = mid_pt_index 180 neg_last_index = mid_pt_index 181 for j in range(step_num): 182 # Increase the angle. 183 angle += step_angle 184 185 # The edge rotation (place all points outside of the rotor angle on the edge). 186 if angle > rotor_angle: 187 axis_angle_to_R(axis, rotor_angle, R) 188 189 # The normal rotation matrix. 190 else: 191 axis_angle_to_R(axis, angle, R) 192 193 # The positive edge. 194 pos_point = dot(R, mid_point - blade_origin) + blade_origin 195 pos_index = mol.atom_add(pdb_record='HETATM', atom_name='BLD', res_name='RTB', res_num=res_num, pos=pos_point, element='N') 196 mol.atom_connect(index1=pos_index, index2=pos_last_index) 197 mol.atom_connect(index1=pos_index, index2=blade_origin_index) 198 199 # The negative edge. 200 neg_point = dot(transpose(R), mid_point - blade_origin) + blade_origin 201 neg_index = mol.atom_add(pdb_record='HETATM', atom_name='BLD', res_name='RTB', res_num=res_num, pos=neg_point, element='N') 202 mol.atom_connect(index1=neg_index, index2=neg_last_index) 203 mol.atom_connect(index1=neg_index, index2=blade_origin_index) 204 205 # Update the indices. 206 pos_last_index = pos_index 207 neg_last_index = neg_index 208 209 # Increment the residue number. 210 res_num += 1

211

Source Code for Module lib.structure.represent.rotor