Package maths_fns :: Module frame_order
[hide private]
[frames] | no frames]

Source Code for Module maths_fns.frame_order

  1  ############################################################################### 
  2  #                                                                             # 
  3  # Copyright (C) 2009-2011 Edward d'Auvergne                                   # 
  4  #                                                                             # 
  5  # This file is part of the program relax.                                     # 
  6  #                                                                             # 
  7  # relax 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 2 of the License, or           # 
 10  # (at your option) any later version.                                         # 
 11  #                                                                             # 
 12  # relax 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 relax; if not, write to the Free Software                        # 
 19  # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA   # 
 20  #                                                                             # 
 21  ############################################################################### 
 22   
 23  # Module docstring. 
 24  """Module containing the target functions of the Frame Order theories.""" 
 25   
 26  # Python module imports. 
 27  from copy import deepcopy 
 28  from numpy import array, dot, float64, ones, transpose, zeros 
 29   
 30  # relax module imports. 
 31  from generic_fns.frame_order import print_frame_order_2nd_degree 
 32  from maths_fns.alignment_tensor import to_5D, to_tensor 
 33  from maths_fns.chi2 import chi2 
 34  from maths_fns.frame_order_matrix_ops import compile_2nd_matrix_free_rotor, compile_2nd_matrix_iso_cone, compile_2nd_matrix_iso_cone_free_rotor, compile_2nd_matrix_iso_cone_torsionless, compile_2nd_matrix_pseudo_ellipse, compile_2nd_matrix_pseudo_ellipse_free_rotor, compile_2nd_matrix_pseudo_ellipse_torsionless, compile_2nd_matrix_rotor, reduce_alignment_tensor 
 35  from maths_fns.rotation_matrix import euler_to_R_zyz as euler_to_R 
 36  from relax_errors import RelaxError 
 37   
 38   
39 -class Frame_order:
40 """Class containing the target function of the optimisation of Frame Order matrix components.""" 41
42 - def __init__(self, model=None, init_params=None, full_tensors=None, red_tensors=None, red_errors=None, full_in_ref_frame=None, frame_order_2nd=None):
43 """Set up the target functions for the Frame Order theories. 44 45 @keyword model: The name of the Frame Order model. 46 @type model: str 47 @keyword init_params: The initial parameter values. 48 @type init_params: numpy float64 array 49 @keyword full_tensors: An array of the {Sxx, Syy, Sxy, Sxz, Syz} values for all full alignment tensors. The format is [Sxx1, Syy1, Sxy1, Sxz1, Syz1, Sxx2, Syy2, Sxy2, Sxz2, Syz2, ..., Sxxn, Syyn, Sxyn, Sxzn, Syzn]. 50 @type full_tensors: numpy nx5D, rank-1 float64 array 51 @keyword red_tensors: An array of the {Sxx, Syy, Sxy, Sxz, Syz} values for all reduced tensors. The array format is the same as for full_tensors. 52 @type red_tensors: numpy nx5D, rank-1 float64 array 53 @keyword red_errors: An array of the {Sxx, Syy, Sxy, Sxz, Syz} errors for all reduced tensors. The array format is the same as for full_tensors. 54 @type red_errors: numpy nx5D, rank-1 float64 array 55 @keyword full_in_ref_frame: An array of flags specifying if the tensor in the reference frame is the full or reduced tensor. 56 @type full_in_ref_frame: numpy rank-1 array 57 @keyword frame_order_2nd: The numerical values of the 2nd degree Frame Order matrix. If supplied, the target functions will optimise directly to these values. 58 @type frame_order_2nd: None or numpy 9D, rank-2 array 59 """ 60 61 # Model test. 62 if not model: 63 raise RelaxError("The type of Frame Order model must be specified.") 64 65 # Store the initial parameter (as a copy). 66 self.params = deepcopy(init_params) 67 68 # Store the agrs. 69 self.model = model 70 self.full_tensors = full_tensors 71 self.red_tensors = red_tensors 72 self.red_errors = red_errors 73 self.full_in_ref_frame = full_in_ref_frame 74 75 # Mix up. 76 if full_tensors != None and frame_order_2nd != None: 77 raise RelaxError("Tensors and Frame Order matrices cannot be supplied together.") 78 if full_tensors == None and frame_order_2nd == None: 79 raise RelaxError("The arguments have been incorrectly supplied, cannot determine the optimisation mode.") 80 81 # Tensor setup. 82 if full_tensors != None: 83 self.__init_tensors() 84 85 # Optimisation to the 2nd degree Frame Order matrix components directly. 86 if model == 'iso cone, free rotor' and frame_order_2nd != None: 87 self.__init_iso_cone_elements(frame_order_2nd) 88 89 # The target function aliases. 90 if model == 'pseudo-ellipse': 91 self.func = self.func_pseudo_ellipse 92 elif model == 'pseudo-ellipse, torsionless': 93 self.func = self.func_pseudo_ellipse_torsionless 94 elif model == 'pseudo-ellipse, free rotor': 95 self.func = self.func_pseudo_ellipse_free_rotor 96 elif model == 'iso cone': 97 self.func = self.func_iso_cone 98 elif model == 'iso cone, torsionless': 99 self.func = self.func_iso_cone_torsionless 100 elif model == 'iso cone, free rotor': 101 self.func = self.func_iso_cone_free_rotor 102 elif model == 'line': 103 self.func = self.func_line 104 elif model == 'line, torsionless': 105 self.func = self.func_line_torsionless 106 elif model == 'line, free rotor': 107 self.func = self.func_line_free_rotor 108 elif model == 'rotor': 109 self.func = self.func_rotor 110 elif model == 'rigid': 111 self.func = self.func_rigid 112 elif model == 'free rotor': 113 self.func = self.func_free_rotor
114 115
116 - def __init_tensors(self):
117 """Set up isotropic cone optimisation against the alignment tensor data.""" 118 119 # Some checks. 120 if self.full_tensors == None or not len(self.full_tensors): 121 raise RelaxError("The full_tensors argument " + repr(self.full_tensors) + " must be supplied.") 122 if self.red_tensors == None or not len(self.red_tensors): 123 raise RelaxError("The red_tensors argument " + repr(self.red_tensors) + " must be supplied.") 124 if self.red_errors == None or not len(self.red_errors): 125 raise RelaxError("The red_errors argument " + repr(self.red_errors) + " must be supplied.") 126 if self.full_in_ref_frame == None or not len(self.full_in_ref_frame): 127 raise RelaxError("The full_in_ref_frame argument " + repr(self.full_in_ref_frame) + " must be supplied.") 128 129 # Tensor set up. 130 self.num_tensors = len(self.full_tensors) / 5 131 self.red_tensors_bc = zeros(self.num_tensors*5, float64) 132 133 # The rotation to the Frame Order eigenframe. 134 self.rot = zeros((3, 3), float64) 135 self.tensor_3D = zeros((3, 3), float64) 136 137 # The cone axis storage and molecular frame z-axis. 138 self.cone_axis = zeros(3, float64) 139 self.z_axis = array([0, 0, 1], float64) 140 141 # Initialise the Frame Order matrices. 142 self.frame_order_2nd = zeros((9, 9), float64)
143 144
145 - def __init_iso_cone_elements(self, frame_order_2nd):
146 """Set up isotropic cone optimisation against the 2nd degree Frame Order matrix elements. 147 148 @keyword frame_order_2nd: The numerical values of the 2nd degree Frame Order matrix. If supplied, the target functions will optimise directly to these values. 149 @type frame_order_2nd: numpy 9D, rank-2 array 150 """ 151 152 # Store the real matrix components. 153 self.data = frame_order_2nd 154 155 # The errors. 156 self.errors = ones((9, 9), float64) 157 158 # The cone axis storage and molecular frame z-axis. 159 self.cone_axis = zeros(3, float64) 160 self.z_axis = array([0, 0, 1], float64) 161 162 # The rotation. 163 self.rot = zeros((3, 3), float64) 164 165 # Initialise the Frame Order matrices. 166 self.frame_order_2nd = zeros((9, 9), float64) 167 168 # Alias the target function. 169 self.func = self.func_iso_cone_elements
170 171
172 - def func_free_rotor(self, params):
173 """Target function for free rotor model optimisation. 174 175 This function optimises against alignment tensors. The Euler angles for the tensor rotation are the first 3 parameters optimised in this model, followed by the polar and azimuthal angles of the cone axis. 176 177 @param params: The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma, theta, phi}. 178 @type params: list of float 179 @return: The chi-squared or SSE value. 180 @rtype: float 181 """ 182 183 # Unpack the parameters. 184 ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi = params 185 186 # Generate the 2nd degree Frame Order super matrix. 187 frame_order_2nd = compile_2nd_matrix_free_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi) 188 189 # Reduce and rotate the tensors. 190 self.reduce_and_rot(0.0, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 191 192 # Return the chi-squared value. 193 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
194 195
196 - def func_iso_cone_elements(self, params):
197 """Target function for isotropic cone model optimisation using the Frame Order matrix. 198 199 This function optimises by directly matching the elements of the 2nd degree Frame Order 200 super matrix. The cone axis spherical angles theta and phi and the cone angle theta are the 201 3 parameters optimised in this model. 202 203 @param params: The vector of parameter values {theta, phi, theta_cone} where the first two are the polar and azimuthal angles of the cone axis theta_cone is the isotropic cone angle. 204 @type params: list of float 205 @return: The chi-squared or SSE value. 206 @rtype: float 207 """ 208 209 # Break up the parameters. 210 theta, phi, theta_cone = params 211 212 # Generate the 2nd degree Frame Order super matrix. 213 self.frame_order_2nd = compile_2nd_matrix_iso_cone_free_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, theta, phi, theta_cone) 214 215 # Make the Frame Order matrix contiguous. 216 self.frame_order_2nd = self.frame_order_2nd.copy() 217 218 # Reshape the numpy arrays for use in the chi2() function. 219 self.data.shape = (81,) 220 self.frame_order_2nd.shape = (81,) 221 self.errors.shape = (81,) 222 223 # Get the chi-squared value. 224 val = chi2(self.data, self.frame_order_2nd, self.errors) 225 226 # Reshape the arrays back to normal. 227 self.data.shape = (9, 9) 228 self.frame_order_2nd.shape = (9, 9) 229 self.errors.shape = (9, 9) 230 231 # Return the chi2 value. 232 return val
233 234
235 - def func_iso_cone(self, params):
236 """Target function for isotropic cone model optimisation. 237 238 This function optimises against alignment tensors. 239 240 @param params: The vector of parameter values {beta, gamma, theta, phi, s1} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and s1 is the isotropic cone order parameter. 241 @type params: list of float 242 @return: The chi-squared or SSE value. 243 @rtype: float 244 """ 245 246 # Unpack the parameters. 247 ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta, sigma_max = params 248 249 # Generate the 2nd degree Frame Order super matrix. 250 frame_order_2nd = compile_2nd_matrix_iso_cone(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta, sigma_max) 251 252 # Reduce and rotate the tensors. 253 self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 254 255 # Return the chi-squared value. 256 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
257 258
259 - def func_iso_cone_free_rotor(self, params):
260 """Target function for free rotor isotropic cone model optimisation. 261 262 This function optimises against alignment tensors. 263 264 @param params: The vector of parameter values {beta, gamma, theta, phi, s1} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and s1 is the isotropic cone order parameter. 265 @type params: list of float 266 @return: The chi-squared or SSE value. 267 @rtype: float 268 """ 269 270 # Unpack the parameters. 271 ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi, cone_s1 = params 272 273 # Generate the 2nd degree Frame Order super matrix. 274 frame_order_2nd = compile_2nd_matrix_iso_cone_free_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi, cone_s1) 275 276 # Reduce and rotate the tensors. 277 self.reduce_and_rot(0.0, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 278 279 # Return the chi-squared value. 280 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
281 282
283 - def func_iso_cone_torsionless(self, params):
284 """Target function for torsionless isotropic cone model optimisation. 285 286 This function optimises against alignment tensors. 287 288 @param params: The vector of parameter values {beta, gamma, theta, phi, cone_theta} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and cone_theta is cone opening angle. 289 @type params: list of float 290 @return: The chi-squared or SSE value. 291 @rtype: float 292 """ 293 294 # Unpack the parameters. 295 ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi, cone_theta = params 296 297 # Generate the 2nd degree Frame Order super matrix. 298 frame_order_2nd = compile_2nd_matrix_iso_cone_torsionless(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi, cone_theta) 299 300 # Reduce and rotate the tensors. 301 self.reduce_and_rot(0.0, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 302 303 # Return the chi-squared value. 304 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
305 306
307 - def func_pseudo_ellipse(self, params):
308 """Target function for pseudo-elliptic cone model optimisation. 309 310 @param params: The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 3 are the pseudo-elliptic cone geometric parameters. 311 @type params: list of float 312 @return: The chi-squared or SSE value. 313 @rtype: float 314 """ 315 316 # Unpack the parameters. 317 ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max = params 318 319 # Generate the 2nd degree Frame Order super matrix. 320 frame_order_2nd = compile_2nd_matrix_pseudo_ellipse(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max) 321 322 # Reduce and rotate the tensors. 323 self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 324 325 # Return the chi-squared value. 326 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
327 328
329 - def func_pseudo_ellipse_free_rotor(self, params):
330 """Target function for free_rotor pseudo-elliptic cone model optimisation. 331 332 @param params: The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 2 are the free_rotor pseudo-elliptic cone geometric parameters. 333 @type params: list of float 334 @return: The chi-squared or SSE value. 335 @rtype: float 336 """ 337 338 # Unpack the parameters. 339 ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y = params 340 341 # Generate the 2nd degree Frame Order super matrix. 342 frame_order_2nd = compile_2nd_matrix_pseudo_ellipse_free_rotor(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y) 343 344 # Reduce and rotate the tensors. 345 self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 346 347 # Return the chi-squared value. 348 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
349 350
351 - def func_pseudo_ellipse_torsionless(self, params):
352 """Target function for torsionless pseudo-elliptic cone model optimisation. 353 354 @param params: The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 2 are the torsionless pseudo-elliptic cone geometric parameters. 355 @type params: list of float 356 @return: The chi-squared or SSE value. 357 @rtype: float 358 """ 359 360 # Unpack the parameters. 361 ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y = params 362 363 # Generate the 2nd degree Frame Order super matrix. 364 frame_order_2nd = compile_2nd_matrix_pseudo_ellipse_torsionless(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y) 365 366 # Reduce and rotate the tensors. 367 self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 368 369 # Return the chi-squared value. 370 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
371 372
373 - def func_rigid(self, params):
374 """Target function for rigid model optimisation. 375 376 This function optimises against alignment tensors. The Euler angles for the tensor rotation are the 3 parameters optimised in this model. 377 378 @param params: The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma}. 379 @type params: list of float 380 @return: The chi-squared or SSE value. 381 @rtype: float 382 """ 383 384 # Unpack the parameters. 385 ave_pos_alpha, ave_pos_beta, ave_pos_gamma = params 386 387 # Reduce and rotate the tensors. 388 self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma) 389 390 # Return the chi-squared value. 391 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
392 393
394 - def func_rotor(self, params):
395 """Target function for rotor model optimisation. 396 397 This function optimises against alignment tensors. The Euler angles for the tensor rotation are the first 3 parameters optimised in this model, followed by the polar and azimuthal angles of the cone axis and the torsion angle restriction. 398 399 @param params: The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma, theta, phi, sigma_max}. 400 @type params: list of float 401 @return: The chi-squared or SSE value. 402 @rtype: float 403 """ 404 405 # Unpack the parameters. 406 ave_pos_alpha, ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi, sigma_max = params 407 408 # Generate the 2nd degree Frame Order super matrix. 409 frame_order_2nd = compile_2nd_matrix_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi, sigma_max) 410 411 # Reduce and rotate the tensors. 412 self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 413 414 # Return the chi-squared value. 415 return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors)
416 417
418 - def reduce_and_rot(self, ave_pos_alpha=None, ave_pos_beta=None, ave_pos_gamma=None, daeg=None):
419 """Reduce and rotate the alignments tensors using the frame order matrix and Euler angles. 420 421 @keyword ave_pos_alpha: The alpha Euler angle describing the average domain position, the tensor rotation. 422 @type ave_pos_alpha: float 423 @keyword ave_pos_beta: The beta Euler angle describing the average domain position, the tensor rotation. 424 @type ave_pos_beta: float 425 @keyword ave_pos_gamma: The gamma Euler angle describing the average domain position, the tensor rotation. 426 @type ave_pos_gamma: float 427 @keyword daeg: The 2nd degree frame order matrix. 428 @type daeg: rank-2, 9D array 429 """ 430 431 # Alignment tensor rotation. 432 euler_to_R(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, self.rot) 433 434 # Back calculate the rotated tensors. 435 for i in range(self.num_tensors): 436 # Tensor indices. 437 index1 = i*5 438 index2 = i*5+5 439 440 # Reduction. 441 if daeg != None: 442 # Reduce the tensor. 443 reduce_alignment_tensor(daeg, self.full_tensors[index1:index2], self.red_tensors_bc[index1:index2]) 444 445 # Convert the reduced tensor to 3D, rank-2 form. 446 to_tensor(self.tensor_3D, self.red_tensors_bc[index1:index2]) 447 448 # No reduction: 449 else: 450 # Convert the original tensor to 3D, rank-2 form. 451 to_tensor(self.tensor_3D, self.full_tensors[index1:index2]) 452 453 # Rotate the tensor (normal R.X.RT rotation). 454 if self.full_in_ref_frame[i]: 455 tensor_3D = dot(self.rot, dot(self.tensor_3D, transpose(self.rot))) 456 457 # Rotate the tensor (inverse RT.X.R rotation). 458 else: 459 tensor_3D = dot(transpose(self.rot), dot(self.tensor_3D, self.rot)) 460 461 # Convert the tensor back to 5D, rank-1 form, as the back-calculated reduced tensor. 462 to_5D(self.red_tensors_bc[index1:index2], tensor_3D)
463