lib.structure.represent.rotor

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

36 -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):

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

118 119

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

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

212

Source Code for Module lib.structure.represent.rotor