generic_fns.angles

1 ############################################################################### 2 # # 3 # Copyright (C) 2003-2012 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 """Module for the manipulation of angular information.""" 24 25 # Python module imports. 26 from math import acos, pi, sin 27 from numpy import dot 28 from warnings import warn 29 30 # relax module imports. 31 from generic_fns import pipes 32 from generic_fns.interatomic import interatomic_loop 33 from generic_fns.mol_res_spin import exists_mol_res_spin_data, generate_spin_id, spin_loop 34 from relax_errors import RelaxError, RelaxNoPdbError, RelaxNoSequenceError, RelaxNoTensorError 35 from relax_warnings import RelaxWarning 36 37

38 -def angle_diff_frame():

39 """Function for calculating the angle defining the XH vector in the diffusion frame.""" 40 41 # Test if the current data pipe exists. 42 pipes.test() 43 44 # Test if the PDB file has been loaded. 45 if not hasattr(cdp, 'structure'): 46 raise RelaxNoPdbError 47 48 # Test if sequence data is loaded. 49 if not exists_mol_res_spin_data(): 50 raise RelaxNoSequenceError 51 52 # Test if the diffusion tensor data is loaded. 53 if not hasattr(cdp, 'diff_tensor'): 54 raise RelaxNoTensorError('diffusion') 55 56 # Sphere. 57 if cdp.diff_tensor.type == 'sphere': 58 return 59 60 # Spheroid. 61 elif cdp.diff_tensor.type == 'spheroid': 62 spheroid_frame() 63 64 # Ellipsoid. 65 elif cdp.diff_tensor.type == 'ellipsoid': 66 raise RelaxError("No coded yet.")

67 68

69 -def ellipsoid_frame():

70 """Calculate the spherical angles of the bond vector in the ellipsoid frame.""" 71 72 # Get the unit vectors Dx, Dy, and Dz of the diffusion tensor axes. 73 Dx, Dy, Dz = diffusion_tensor.unit_axes() 74 75 # Spin loop. 76 for spin, mol_name, res_num, res_name in spin_loop(full_info=True): 77 # Test if the vector exists. 78 if not hasattr(spin, 'xh_vect'): 79 # Get the spin id string. 80 spin_id = generate_spin_id(mol_name=mol_name, res_num=res_num, res_name=res_name, spin_num=spin.num, spin_name=spin.name) 81 82 # Throw a warning. 83 warn(RelaxWarning("No angles could be calculated for the spin " + repr(spin_id) + ".")) 84 85 # Skip the spin. 86 continue 87 88 # dz and dx direction cosines. 89 dz = dot(Dz, spin.xh_vect) 90 dx = dot(Dx, spin.xh_vect) 91 92 # Calculate the polar angle theta. 93 spin.theta = acos(dz) 94 95 # Calculate the azimuthal angle phi. 96 spin.phi = acos(dx / sin(spin.theta))

97 98

99 -def fold_spherical_angles(theta, phi, theta_lower=0, theta_upper=2*pi, theta_window=2*pi, phi_lower=0, phi_upper=2*pi, phi_window=2*pi):

100 """Fold the spherical angles taking symmetry into account. 101 102 The angles will be folded between:: 103 104 0 <= theta <= pi, 105 0 <= phi <= 2*pi, 106 107 @param theta: The azimuthal angle. 108 @type theta: float 109 @param phi: The polar angle. 110 @type phi: float 111 @param theta_lower: The theta angle lower bound (defaults to 0). 112 @type theta_lower: float 113 @param theta_upper: The theta angle upper bound (defaults to 2*pi). 114 @type theta_upper: float 115 @param theta_window: The size of the theta angle window where symmetry exists (defaults to 2*pi). 116 @type theta_window: float 117 @param phi_lower: The phi angle lower bound (defaults to 0). 118 @type phi_lower: float 119 @param phi_upper: The phi angle upper bound (defaults to 2*pi). 120 @type phi_upper: float 121 @param phi_window: The size of the phi angle window where symmetry exists (defaults to 2*pi). 122 @type phi_window: float 123 @return: The folded angles, theta and phi. 124 @rtype: float 125 """ 126 127 # Check the bounds and window. 128 if theta_window - (theta_upper - theta_lower) > 1e-7: 129 raise RelaxError("The theta angle lower and upper bounds [%s, %s] do not match the window size of %s." % (theta_lower, theta_upper, theta_window)) 130 if phi_window - (phi_upper - phi_lower) > 1e-7: 131 raise RelaxError("The phi angle lower and upper bounds [%s, %s] do not match the window size of %s." % (phi_lower, phi_upper, phi_window)) 132 133 # First wrap the angles. 134 theta = wrap_angles(theta, theta_lower, theta_upper, theta_window) 135 phi = wrap_angles(phi, phi_lower, phi_upper, phi_window) 136 137 # Then remove the symmetry to the lower half of phi. 138 if phi >= phi_upper - phi_window/2.0: 139 theta = pi - theta 140 phi = phi - pi 141 142 # Wrap again if necessary. 143 theta = wrap_angles(theta, theta_lower, theta_upper, theta_window) 144 phi = wrap_angles(phi, phi_lower, phi_upper, phi_window) 145 146 # Return the folded angles. 147 return theta, phi

148 149

150 -def spheroid_frame():

151 """Function for calculating the angle alpha of the XH vector within the spheroid frame.""" 152 153 # Loop over the interatomic info. 154 for interatom in interatomic_loop(): 155 # Test if the vector exists. 156 if not hasattr(interatom, 'vector'): 157 # Throw a warning. 158 warn(RelaxWarning("No angles could be calculated for the spin pair '%s' and '%s'." % (interatom.spin_id1, interatom.spin_id2))) 159 160 # Skip the container. 161 continue 162 163 # Calculate alpha. 164 interatom.alpha = acos(dot(cdp.diff_tensor.Dpar_unit, interatom.vector))

165 166

167 -def wrap_angles(angle, lower, upper, window=2*pi):

168 """Convert the given angle to be between the lower and upper values. 169 170 @param angle: The starting angle. 171 @type angle: float 172 @param lower: The lower bound. 173 @type lower: float 174 @param upper: The upper bound. 175 @type upper: float 176 @param window: The size of the window where symmetry exists (defaults to 2pi). 177 @type window: float 178 @return: The wrapped angle. 179 @rtype: float 180 """ 181 182 # Check the bounds and window. 183 if window - (upper - lower) > 1e-7: 184 raise RelaxError("The lower and upper bounds [%s, %s] do not match the window size of %s." % (lower, upper, window)) 185 186 # Keep wrapping until the angle is within the limits. 187 while True: 188 # The angle is too big. 189 if angle > upper: 190 angle = angle - window 191 192 # The angle is too small. 193 elif angle < lower: 194 angle = angle + window 195 196 # Inside the window, so stop wrapping. 197 else: 198 break 199 200 # Return the wrapped angle. 201 return angle

202

Source Code for Module generic_fns.angles