Package generic_fns :: Module align_tensor
[hide private]
[frames] | no frames]

Source Code for Module generic_fns.align_tensor

   1  ############################################################################### 
   2  #                                                                             # 
   3  # Copyright (C) 2003-2012 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 functions for the handling of alignment tensors.""" 
  25   
  26  # Python module imports. 
  27  from copy import deepcopy 
  28  from math import pi, sqrt 
  29  from numpy import arccos, dot, float64, linalg, zeros 
  30  from numpy.linalg import norm 
  31  from re import search 
  32  import sys 
  33   
  34  # relax module imports. 
  35  from angles import wrap_angles 
  36  from data.align_tensor import AlignTensorList 
  37  from generic_fns import pipes 
  38  from physical_constants import g1H, h_bar, kB, mu0, return_gyromagnetic_ratio 
  39  from relax_errors import RelaxError, RelaxNoTensorError, RelaxStrError, RelaxTensorError, RelaxUnknownParamCombError, RelaxUnknownParamError 
  40   
  41   
42 -def align_data_exists(tensor, pipe=None):
43 """Function for determining if alignment data exists in the current data pipe. 44 45 @param tensor: The alignment tensor identification string. 46 @type tensor: str 47 @param pipe: The data pipe to search for data in. 48 @type pipe: str 49 @return: The answer to the question. 50 @rtype: bool 51 """ 52 53 # The data pipe to check. 54 if pipe == None: 55 pipe = pipes.cdp_name() 56 57 # Get the data pipe. 58 pipe = pipes.get_pipe(pipe) 59 60 # Test if an alignment tensor corresponding to the arg 'tensor' exists. 61 if hasattr(pipe, 'align_tensors'): 62 for data in pipe.align_tensors: 63 if data.name == tensor: 64 return True 65 else: 66 return False
67 68
69 -def all_tensors_fixed():
70 """Determine if all alignment tensors are fixed. 71 72 @return: True if all tensors are fixed, False otherwise. 73 @rtype: bool 74 """ 75 76 # Loop over the tensors. 77 for i in range(len(cdp.align_tensors)): 78 # Not fixed, so return False. 79 if not cdp.align_tensors[i].fixed: 80 return False 81 82 # All tensors are fixed. 83 return True
84 85
86 -def calc_chi_tensor(A, B0, T):
87 """Convert the alignment tensor into the magnetic susceptibility (chi) tensor. 88 89 A can be either the full tensor (3D or 5D), a component Aij of the tensor, Aa, or Ar, anything that can be multiplied by the constants to convert from one to the other. 90 91 92 @param A: The alignment tensor or alignment tensor component. 93 @type A: numpy array or float 94 @param B0: The magnetic field strength in Hz. 95 @type B0: float 96 @param T: The temperature in Kalvin. 97 @type T: float 98 @return: A multiplied by the PCS constant. 99 @rtype: numpy array or float 100 """ 101 102 # B0 in Tesla. 103 B0 = 2.0 * pi * B0 / g1H 104 105 # The conversion factor. 106 conv = 15.0 * mu0 * kB * T / B0**2 107 108 # Return the converted value. 109 return conv * A
110 111
112 -def copy(tensor_from=None, pipe_from=None, tensor_to=None, pipe_to=None):
113 """Function for copying alignment tensor data from one data pipe to another. 114 115 @param tensor_from: The identification string of the alignment tensor to copy the data from. 116 @type tensor_from: str 117 @param pipe_from: The data pipe to copy the alignment tensor data from. This defaults to the 118 current data pipe. 119 @type pipe_from: str 120 @param tensor_to: The identification string of the alignment tensor to copy the data to. 121 @type tensor_to: str 122 @param pipe_to: The data pipe to copy the alignment tensor data to. This defaults to the 123 current data pipe. 124 @type pipe_to: str 125 """ 126 127 # Defaults. 128 if tensor_from == tensor_to and pipe_from == None and pipe_to == None: 129 raise RelaxError("The pipe_from and pipe_to arguments cannot both be set to None when the tensor names are the same.") 130 elif pipe_from == None: 131 pipe_from = pipes.cdp_name() 132 elif pipe_to == None: 133 pipe_to = pipes.cdp_name() 134 135 # Test if the pipe_from and pipe_to data pipes exist. 136 pipes.test(pipe_from) 137 pipes.test(pipe_to) 138 139 # Get the data pipes. 140 dp_from = pipes.get_pipe(pipe_from) 141 dp_to = pipes.get_pipe(pipe_to) 142 143 # Test if pipe_from contains alignment tensor data. 144 if not align_data_exists(tensor_from, pipe_from): 145 raise RelaxNoTensorError('alignment') 146 147 # Test if pipe_to contains alignment tensor data. 148 if align_data_exists(tensor_to, pipe_to): 149 raise RelaxTensorError('alignment') 150 151 # Create the align_tensors dictionary if it doesn't yet exist. 152 if not hasattr(dp_to, 'align_tensors'): 153 dp_to.align_tensors = AlignTensorList() 154 155 # Find the tensor index. 156 index_from = get_tensor_index(tensor_from, pipe_from) 157 index_to = get_tensor_index(tensor_to, pipe_to) 158 159 # Copy the data. 160 if index_to == None: 161 dp_to.align_tensors.append(deepcopy(dp_from.align_tensors[index_from])) 162 else: 163 dp_to.align_tensors[index_to] = deepcopy(dp_from.align_tensors[index_from])
164 165
166 -def data_names():
167 """Function for returning a list of names of data structures associated with the sequence.""" 168 169 names = [ 'align_params' ] 170 171 return names
172 173
174 -def default_value(param):
175 """Return the default values for the alignment tensor parameters. 176 177 @param param: The name of the parameter. 178 @type param: str 179 @return: The default value, which for all parameters is set to zero. 180 @rtype: float 181 """ 182 183 # Return 0.0. 184 return 0.0
185 186 # User function documentation. 187 __default_value_prompt_doc__ = """ 188 Alignment tensor parameter default values 189 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 190 191 ________________________________________________________________________ 192 | | | | 193 | Data type | Object name | Value | 194 |________________________|____________________|________________________| 195 | | | | 196 | Axx | 'Axx' | 0.0 | 197 | | | | 198 | Ayy | 'Ayy' | 0.0 | 199 | | | | 200 | Azz | 'Azz' | 0.0 | 201 | | | | 202 | Axxyy | 'Axxyy' | 0.0 | 203 | | | | 204 | Axy | 'Axy' | 0.0 | 205 | | | | 206 | Axz | 'Axz' | 0.0 | 207 | | | | 208 | Ayz | 'Ayz' | 0.0 | 209 | | | | 210 | alpha | 'alpha' | 0.0 | 211 | | | | 212 | beta | 'beta' | 0.0 | 213 | | | | 214 | gamma | 'gamma' | 0.0 | 215 |________________________|____________________|________________________| 216 217 """ 218 219
220 -def delete(tensor=None):
221 """Function for deleting alignment tensor data. 222 223 @param tensor: The alignment tensor identification string. 224 @type tensor: str or None 225 """ 226 227 # Test if the current data pipe exists. 228 pipes.test() 229 230 # Test if alignment tensor data exists. 231 if tensor and not align_data_exists(tensor): 232 raise RelaxNoTensorError('alignment') 233 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'): 234 raise RelaxNoTensorError('alignment') 235 236 # The tensor list. 237 if tensor: 238 tensors = [tensor] 239 else: 240 tensors = cdp.align_ids 241 242 # Loop over the tensors. 243 for tensor in tensors: 244 # Print out. 245 print("Removing the '%s' tensor." % tensor) 246 247 # Find the tensor index. 248 index = get_tensor_index(tensor) 249 250 # Delete the alignment data. 251 cdp.align_tensors.pop(index) 252 253 # Delete the alignment tensor list if empty. 254 if not len(cdp.align_tensors): 255 del(cdp.align_tensors)
256 257
258 -def display(tensor):
259 """Function for displaying the alignment tensor. 260 261 @param tensor: The alignment tensor identification string. 262 @type tensor: str or None 263 """ 264 265 # Test if the current data pipe exists. 266 pipes.test() 267 268 # Test if alignment tensor data exists. 269 if tensor and not align_data_exists(tensor): 270 raise RelaxNoTensorError('alignment') 271 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'): 272 raise RelaxNoTensorError('alignment') 273 274 # Construct the tensor list. 275 tensor_list = [] 276 if not tensor: 277 for tensor_cont in cdp.align_tensors: 278 tensor_list.append(tensor_cont.name) 279 else: 280 tensor_list.append(tensor) 281 282 # Loop over the tensors. 283 for tensor in tensor_list: 284 # Test if alignment tensor data exists. 285 if not align_data_exists(tensor): 286 raise RelaxNoTensorError('alignment') 287 288 # Pull out the tensor. 289 data = get_tensor_object(tensor) 290 291 # Header. 292 head = "# Tensor: %s #" % tensor 293 print("\n\n\n" + '#' * len(head) + "\n" + head + "\n" + '#' * len(head)) 294 295 296 # The Saupe matrix. 297 ################### 298 299 title = "# Saupe order matrix." 300 print("\n\n" + title + '\n' + '#'*len(title) + '\n') 301 302 # The parameter set {Sxx, Syy, Sxy, Sxz, Syz}. 303 print("# 5D, rank-1 notation {Sxx, Syy, Sxy, Sxz, Syz}:") 304 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Sxx, data.Syy, data.Sxy, data.Sxz, data.Syz))) 305 306 # The parameter set {Szz, Sxx-yy, Sxy, Sxz, Syz}. 307 print("# 5D, rank-1 notation {Szz, Sxx-yy, Sxy, Sxz, Syz} (the Pales default format).") 308 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Szz, data.Sxxyy, data.Sxy, data.Sxz, data.Syz))) 309 310 # 3D form. 311 print("# 3D, rank-2 notation.") 312 print("%s" % (data.S)) 313 314 315 # The alignment tensor. 316 ####################### 317 318 title = "# Alignment tensor." 319 print("\n\n" + title + '\n' + '#'*len(title) + '\n') 320 321 # The parameter set {Axx, Ayy, Axy, Axz, Ayz}. 322 print("# 5D, rank-1 notation {Axx, Ayy, Axy, Axz, Ayz}:") 323 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Axx, data.Ayy, data.Axy, data.Axz, data.Ayz))) 324 325 # The parameter set {Azz, Axx-yy, Axy, Axz, Ayz}. 326 print("# 5D, rank-1 notation {Azz, Axx-yy, Axy, Axz, Ayz} (the Pales default format).") 327 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Azz, data.Axxyy, data.Axy, data.Axz, data.Ayz))) 328 329 # 3D form. 330 print("# 3D, rank-2 notation.") 331 print("%s" % data.A) 332 333 334 # The probability tensor. 335 ######################### 336 337 title = "# Probability tensor." 338 print("\n\n" + title + '\n' + '#'*len(title) + '\n') 339 340 # The parameter set {Pxx, Pyy, Pxy, Pxz, Pyz}. 341 print("# 5D, rank-1 notation {Pxx, Pyy, Pxy, Pxz, Pyz}:") 342 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Pxx, data.Pyy, data.Pxy, data.Pxz, data.Pyz))) 343 344 # The parameter set {Pzz, Pxx-yy, Pxy, Pxz, Pyz}. 345 print("# 5D, rank-1 notation {Pzz, Pxx-yy, Pxy, Pxz, Pyz}.") 346 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Pzz, data.Pxxyy, data.Pxy, data.Pxz, data.Pyz))) 347 348 # 3D form. 349 print("# 3D, rank-2 notation.") 350 print("%s" % data.P) 351 352 353 # The magnetic susceptibility tensor. 354 ##################################### 355 356 title = "# Magnetic susceptibility tensor." 357 print("\n\n" + title + '\n' + '#'*len(title) + '\n') 358 chi_tensor = True 359 360 # The field strength. 361 print("# The magnetic field strength (MHz):") 362 if hasattr(cdp, 'frq') and tensor in cdp.frq: 363 print(("%s\n" % (cdp.frq[tensor] / 1e6))) 364 else: 365 print(("Not set.\n")) 366 chi_tensor = False 367 368 # The temperature. 369 print("# The temperature (K):") 370 if hasattr(cdp, 'temperature') and tensor in cdp.temperature: 371 print(("%s\n" % cdp.temperature[tensor])) 372 else: 373 print(("Not set.\n")) 374 chi_tensor = False 375 376 # No chi tensor. 377 if not chi_tensor: 378 print("# The chi tensor:\nN/A.\n") 379 380 # Calculate the chi tensor. 381 else: 382 # Conversions. 383 chi_xx = calc_chi_tensor(data.Axx, cdp.frq[tensor], cdp.temperature[tensor]) 384 chi_xy = calc_chi_tensor(data.Axy, cdp.frq[tensor], cdp.temperature[tensor]) 385 chi_xz = calc_chi_tensor(data.Axz, cdp.frq[tensor], cdp.temperature[tensor]) 386 chi_yy = calc_chi_tensor(data.Ayy, cdp.frq[tensor], cdp.temperature[tensor]) 387 chi_yz = calc_chi_tensor(data.Ayz, cdp.frq[tensor], cdp.temperature[tensor]) 388 chi_zz = calc_chi_tensor(data.Azz, cdp.frq[tensor], cdp.temperature[tensor]) 389 chi_xxyy = calc_chi_tensor(data.Axxyy, cdp.frq[tensor], cdp.temperature[tensor]) 390 chi = calc_chi_tensor(data.A, cdp.frq[tensor], cdp.temperature[tensor]) 391 392 # The parameter set {chi_xx, chi_yy, chi_xy, chi_xz, chi_yz}. 393 print("# 5D, rank-1 notation {chi_xx, chi_yy, chi_xy, chi_xz, chi_yz}:") 394 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (chi_xx, chi_yy, chi_xy, chi_xz, chi_yz))) 395 396 # The parameter set {chi_zz, chi_xx-yy, chi_xy, chi_xz, chi_yz}. 397 print("# 5D, rank-1 notation {chi_zz, chi_xx-yy, chi_xy, chi_xz, chi_yz}.") 398 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (chi_zz, chi_xxyy, chi_xy, chi_xz, chi_yz))) 399 400 # 3D form. 401 print("# 3D, rank-2 notation.") 402 print("%s" % chi) 403 404 405 # The Eigensystem. 406 ################## 407 408 title = "# Eigensystem." 409 print("\n\n" + title + '\n' + '#'*len(title) + '\n') 410 411 # Eigenvalues. 412 print("# Saupe order matrix eigenvalues {Sxx, Syy, Szz}.") 413 print(("[%25.12e, %25.12e, %25.12e]\n" % (data.S_diag[0, 0], data.S_diag[1, 1], data.S_diag[2, 2]))) 414 print("# Alignment tensor eigenvalues {Axx, Ayy, Azz}.") 415 print(("[%25.12e, %25.12e, %25.12e]\n" % (data.A_diag[0, 0], data.A_diag[1, 1], data.A_diag[2, 2]))) 416 print("# Probability tensor eigenvalues {Pxx, Pyy, Pzz}.") 417 print(("[%25.12e, %25.12e, %25.12e]\n" % (data.P_diag[0, 0], data.P_diag[1, 1], data.P_diag[2, 2]))) 418 if chi_tensor: 419 chi_diag = calc_chi_tensor(data.A_diag, cdp.frq[tensor], cdp.temperature[tensor]) 420 print("# Magnetic susceptibility eigenvalues {chi_xx, chi_yy, chi_zz}.") 421 print(("[%25.12e, %25.12e, %25.12e]\n" % (chi_diag[0, 0], chi_diag[1, 1], chi_diag[2, 2]))) 422 423 # Eigenvectors. 424 print("# Eigenvector x.") 425 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.unit_x[0], data.unit_x[1], data.unit_x[2]))) 426 print("# Eigenvector y.") 427 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.unit_y[0], data.unit_y[1], data.unit_y[2]))) 428 print("# Eigenvector z.") 429 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.unit_z[0], data.unit_z[1], data.unit_z[2]))) 430 431 # Rotation matrix. 432 print("# Rotation matrix.") 433 print("%s\n" % data.rotation) 434 435 # zyz. 436 print("# Euler angles in zyz notation {alpha, beta, gamma}.") 437 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.euler[0], data.euler[1], data.euler[2]))) 438 439 440 # Geometric description. 441 ######################## 442 443 title = "# Geometric description." 444 print("\n\n" + title + '\n' + '#'*len(title) + '\n') 445 446 # The GDO. 447 print("# Generalized degree of order (GDO).") 448 print("GDO = %-25.12e\n" % gdo(data.A)) 449 450 # Anisotropy. 451 print("# Alignment tensor axial component (Aa = 3/2 * Azz, where Aii are the eigenvalues).") 452 print("Aa = %-25.12e\n" % data.Aa) 453 454 # Rhombicity. 455 print("# Rhombic component (Ar = Axx - Ayy, where Aii are the eigenvalues).") 456 print("Ar = %-25.12e\n" % data.Ar) 457 print("# Rhombicity (R = Ar / Aa).") 458 print("R = %-25.12f\n" % data.R) 459 print("# Asymmetry parameter (eta = (Axx - Ayy) / Azz, where Aii are the eigenvalues).") 460 print("eta = %-25.12f\n" % data.eta) 461 462 # Magnetic susceptibility tensor. 463 if chi_tensor: 464 # Chi tensor anisotropy. 465 print("# Magnetic susceptibility axial parameter (chi_ax = chi_zz - (chi_xx + chi_yy)/2, where chi_ii are the eigenvalues).") 466 print("chi_ax = %-25.12e\n" % (chi_diag[2, 2] - (chi_diag[0, 0] + chi_diag[1, 1])/2.0)) 467 468 # Chi tensor rhombicity. 469 print("# Magnetic susceptibility rhombicity parameter (chi_rh = chi_xx - chi_yy, where chi_ii are the eigenvalues).") 470 print("chi_rh = %-25.12e\n" % (chi_diag[0, 0] - chi_diag[1, 1])) 471 472 # Some white space. 473 print("\n\n\n")
474 475
476 -def fix(id=None, fixed=True):
477 """Fix the alignment tensor during optimisation. 478 479 @keyword id: The alignment tensor ID string. If set to None, then all alignment tensors will be fixed. 480 @type id: str or None 481 @keyword fixed: If True, the alignment tensor will be fixed during optimisation. If False, the alignment tensors will be optimised. 482 @type fixed: bool 483 """ 484 485 # Test if the current data pipe exists. 486 pipes.test() 487 488 # Test if alignment tensor data exists. 489 if not hasattr(cdp, 'align_tensors') or not hasattr(cdp, 'align_ids'): 490 raise RelaxNoTensorError('alignment') 491 492 # Loop over the tensors. 493 for i in range(len(cdp.align_tensors)): 494 # ID match. 495 if id and cdp.align_tensors[i].name == id: 496 cdp.align_tensors[i].fixed = fixed 497 498 # Set all tensor flags. 499 if id == None: 500 cdp.align_tensors[i].fixed = fixed
501 502
503 -def fold_angles(sim_index=None):
504 """Wrap the Euler angles and remove the glide reflection and translational symmetries. 505 506 Wrap the angles such that:: 507 508 0 <= alpha <= 2pi, 509 0 <= beta <= pi, 510 0 <= gamma <= 2pi. 511 512 For the simulated values, the angles are wrapped as:: 513 514 alpha - pi <= alpha_sim <= alpha + pi 515 beta - pi/2 <= beta_sim <= beta + pi/2 516 gamma - pi <= gamma_sim <= gamma + pi 517 518 519 @param sim_index: The simulation index. If set to None then the actual values will be folded 520 rather than the simulation values. 521 @type sim_index: int or None 522 """ 523 524 525 # Wrap the angles. 526 ################## 527 528 # Get the current angles. 529 alpha = cdp.align_tensors.alpha 530 beta = cdp.align_tensors.beta 531 gamma = cdp.align_tensors.gamma 532 533 # Simulated values. 534 if sim_index != None: 535 alpha_sim = cdp.align_tensors.alpha_sim[sim_index] 536 beta_sim = cdp.align_tensors.beta_sim[sim_index] 537 gamma_sim = cdp.align_tensors.gamma_sim[sim_index] 538 539 # Normal value. 540 if sim_index == None: 541 cdp.align_tensors.alpha = wrap_angles(alpha, 0.0, 2.0*pi) 542 cdp.align_tensors.beta = wrap_angles(beta, 0.0, 2.0*pi) 543 cdp.align_tensors.gamma = wrap_angles(gamma, 0.0, 2.0*pi) 544 545 # Simulation values. 546 else: 547 cdp.align_tensors.alpha_sim[sim_index] = wrap_angles(alpha_sim, alpha - pi, alpha + pi) 548 cdp.align_tensors.beta_sim[sim_index] = wrap_angles(beta_sim, beta - pi, beta + pi) 549 cdp.align_tensors.gamma_sim[sim_index] = wrap_angles(gamma_sim, gamma - pi, gamma + pi) 550 551 552 # Remove the glide reflection and translational symmetries. 553 ########################################################### 554 555 # Normal value. 556 if sim_index == None: 557 # Fold beta inside 0 and pi. 558 if cdp.align_tensors.beta >= pi: 559 cdp.align_tensors.alpha = pi - cdp.align_tensors.alpha 560 cdp.align_tensors.beta = cdp.align_tensors.beta - pi 561 562 # Simulation values. 563 else: 564 # Fold beta_sim inside beta-pi/2 and beta+pi/2. 565 if cdp.align_tensors.beta_sim[sim_index] >= cdp.align_tensors.beta + pi/2.0: 566 cdp.align_tensors.alpha_sim[sim_index] = pi - cdp.align_tensors.alpha_sim[sim_index] 567 cdp.align_tensors.beta_sim[sim_index] = cdp.align_tensors.beta_sim[sim_index] - pi 568 elif cdp.align_tensors.beta_sim[sim_index] <= cdp.align_tensors.beta - pi/2.0: 569 cdp.align_tensors.alpha_sim[sim_index] = pi - cdp.align_tensors.alpha_sim[sim_index] 570 cdp.align_tensors.beta_sim[sim_index] = cdp.align_tensors.beta_sim[sim_index] + pi
571 572
573 -def gdo(A):
574 """Calculate the generalized degree of order (GDO) for the given alignment tensor. 575 576 @param A: The alignment tensor. 577 @type A: rank-2, 3D numpy array 578 @return: The GDO value. 579 @rtype: float 580 """ 581 582 # The matrix norm. 583 gdo = sqrt(3.0/2.0) * norm(A) 584 585 # Return the GDO. 586 return gdo
587 588
589 -def get_tensor_index(tensor, pipe=None):
590 """Function for returning the index corresponding to the 'tensor' argument. 591 592 @param tensor: The alignment tensor identification string. 593 @type tensor: str 594 @param pipe: The data pipe to search for data in. 595 @type pipe: str 596 @return: The index corresponding to the 'tensor' arg. 597 @rtype: int 598 """ 599 600 # The data pipe to check. 601 if pipe == None: 602 pipe = pipes.cdp_name() 603 604 # Get the data pipe. 605 dp = pipes.get_pipe(pipe) 606 607 # Init. 608 index = None 609 610 # Loop over the tensors. 611 for i in xrange(len(dp.align_tensors)): 612 if dp.align_tensors[i].name == tensor: 613 index = i 614 615 # Return the index. 616 return index
617 618
619 -def get_tensor_object(tensor, pipe=None):
620 """Function for returning the AlignTensorData instance corresponding to the 'tensor' argument. 621 622 @param tensor: The alignment tensor identification string. 623 @type tensor: str 624 @param pipe: The data pipe to search for data in. 625 @type pipe: str 626 @return: The alignment tensor object corresponding to the 'tensor' arg. 627 @rtype: AlignTensorData instance 628 """ 629 630 # The data pipe to check. 631 if pipe == None: 632 pipe = pipes.cdp_name() 633 634 # Init. 635 data = None 636 637 # Loop over the tensors. 638 for i in xrange(len(cdp.align_tensors)): 639 if cdp.align_tensors[i].name == tensor: 640 data = cdp.align_tensors[i] 641 642 # Return the object. 643 return data
644 645
646 -def init(tensor=None, params=None, scale=1.0, angle_units='deg', param_types=0, errors=False):
647 """Function for initialising the alignment tensor. 648 649 @keyword tensor: The alignment tensor identification string. 650 @type tensor: str 651 @keyword params: The alignment tensor parameters. 652 @type params: float 653 @keyword scale: The alignment tensor eigenvalue scaling value. 654 @type scale: float 655 @keyword angle_units: The units for the angle parameters (either 'deg' or 'rad'). 656 @type angle_units: str 657 @keyword param_types: The type of parameters supplied. The flag values correspond to, 0: 658 {Axx, Ayy, Axy, Axz, Ayz}, and 1: {Azz, Axx-yy, Axy, Axz, Ayz}. 659 @type param_types: int 660 @keyword errors: A flag which determines if the alignment tensor data or its errors are 661 being input. 662 @type errors: bool 663 """ 664 665 # Test if the current data pipe exists. 666 pipes.test() 667 668 # Test if alignment tensor data already exists. 669 if errors and not align_data_exists(tensor): 670 raise RelaxNoTensorError('alignment') 671 672 # Check the validity of the angle_units argument. 673 valid_types = ['deg', 'rad'] 674 if not angle_units in valid_types: 675 raise RelaxError("The alignment tensor 'angle_units' argument " + repr(angle_units) + " should be either 'deg' or 'rad'.") 676 677 # Add the tensor ID to the current data pipe. 678 if not hasattr(cdp, 'align_ids'): 679 cdp.align_ids = [] 680 if tensor not in cdp.align_ids: 681 cdp.align_ids.append(tensor) 682 683 # Add the align_tensors object to the data pipe. 684 if not errors: 685 # Initialise the super structure. 686 if not hasattr(cdp, 'align_tensors'): 687 cdp.align_tensors = AlignTensorList() 688 689 # Add the tensor, if it doesn't already exist. 690 if tensor not in cdp.align_tensors.names(): 691 cdp.align_tensors.add_item(tensor) 692 693 # Get the tensor. 694 tensor_obj = get_tensor_object(tensor) 695 696 # {Sxx, Syy, Sxy, Sxz, Syz}. 697 if param_types == 0: 698 # Unpack the tuple. 699 Sxx, Syy, Sxy, Sxz, Syz = params 700 701 # Scaling. 702 Sxx = Sxx * scale 703 Syy = Syy * scale 704 Sxy = Sxy * scale 705 Sxz = Sxz * scale 706 Syz = Syz * scale 707 708 # Set the parameters. 709 set(tensor=tensor_obj, value=[Sxx, Syy, Sxy, Sxz, Syz], param=['Sxx', 'Syy', 'Sxy', 'Sxz', 'Syz'], errors=errors) 710 711 # {Szz, Sxx-yy, Sxy, Sxz, Syz}. 712 elif param_types == 1: 713 # Unpack the tuple. 714 Szz, Sxxyy, Sxy, Sxz, Syz = params 715 716 # Scaling. 717 Szz = Szz * scale 718 Sxxyy = Sxxyy * scale 719 Sxy = Sxy * scale 720 Sxz = Sxz * scale 721 Syz = Syz * scale 722 723 # Set the parameters. 724 set(tensor=tensor_obj, value=[Szz, Sxxyy, Sxy, Sxz, Syz], param=['Szz', 'Sxxyy', 'Sxy', 'Sxz', 'Syz'], errors=errors) 725 726 # {Axx, Ayy, Axy, Axz, Ayz}. 727 elif param_types == 2: 728 # Unpack the tuple. 729 Axx, Ayy, Axy, Axz, Ayz = params 730 731 # Scaling. 732 Axx = Axx * scale 733 Ayy = Ayy * scale 734 Axy = Axy * scale 735 Axz = Axz * scale 736 Ayz = Ayz * scale 737 738 # Set the parameters. 739 set(tensor=tensor_obj, value=[Axx, Ayy, Axy, Axz, Ayz], param=['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], errors=errors) 740 741 # {Azz, Axx-yy, Axy, Axz, Ayz}. 742 elif param_types == 3: 743 # Unpack the tuple. 744 Azz, Axxyy, Axy, Axz, Ayz = params 745 746 # Scaling. 747 Azz = Azz * scale 748 Axxyy = Axxyy * scale 749 Axy = Axy * scale 750 Axz = Axz * scale 751 Ayz = Ayz * scale 752 753 # Set the parameters. 754 set(tensor=tensor_obj, value=[Azz, Axxyy, Axy, Axz, Ayz], param=['Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz'], errors=errors) 755 756 # {Axx, Ayy, Axy, Axz, Ayz}. 757 elif param_types == 4: 758 # Unpack the tuple. 759 Axx, Ayy, Axy, Axz, Ayz = params 760 761 # Get the bond length. 762 r = None 763 for spin in spin_loop(): 764 # First spin. 765 if r == None: 766 r = spin.r 767 768 # Different value. 769 if r != spin.r: 770 raise RelaxError("Not all spins have the same bond length.") 771 772 # Scaling. 773 scale = scale / kappa() * r**3 774 Axx = Axx * scale 775 Ayy = Ayy * scale 776 Axy = Axy * scale 777 Axz = Axz * scale 778 Ayz = Ayz * scale 779 780 # Set the parameters. 781 set(tensor=tensor_obj, value=[Axx, Ayy, Axy, Axz, Ayz], param=['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], errors=errors) 782 783 # {Azz, Axx-yy, Axy, Axz, Ayz}. 784 elif param_types == 5: 785 # Unpack the tuple. 786 Azz, Axxyy, Axy, Axz, Ayz = params 787 788 # Get the bond length. 789 r = None 790 for spin in spin_loop(): 791 # First spin. 792 if r == None: 793 r = spin.r 794 795 # Different value. 796 if r != spin.r: 797 raise RelaxError("Not all spins have the same bond length.") 798 799 # Scaling. 800 scale = scale / kappa() * r**3 801 Azz = Azz * scale 802 Axxyy = Axxyy * scale 803 Axy = Axy * scale 804 Axz = Axz * scale 805 Ayz = Ayz * scale 806 807 # Set the parameters. 808 set(tensor=tensor_obj, value=[Azz, Axxyy, Axy, Axz, Ayz], param=['Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz'], errors=errors) 809 810 # {Pxx, Pyy, Pxy, Pxz, Pyz}. 811 elif param_types == 6: 812 # Unpack the tuple. 813 Pxx, Pyy, Pxy, Pxz, Pyz = params 814 815 # Scaling. 816 Pxx = Pxx * scale 817 Pyy = Pyy * scale 818 Pxy = Pxy * scale 819 Pxz = Pxz * scale 820 Pyz = Pyz * scale 821 822 # Set the parameters. 823 set(tensor=tensor_obj, value=[Pxx, Pyy, Pxy, Pxz, Pyz], param=['Pxx', 'Pyy', 'Pxy', 'Pxz', 'Pyz'], errors=errors) 824 825 # {Pzz, Pxx-yy, Pxy, Pxz, Pyz}. 826 elif param_types == 7: 827 # Unpack the tuple. 828 Pzz, Pxxyy, Pxy, Pxz, Pyz = params 829 830 # Scaling. 831 Pzz = Pzz * scale 832 Pxxyy = Pxxyy * scale 833 Pxy = Pxy * scale 834 Pxz = Pxz * scale 835 Pyz = Pyz * scale 836 837 # Set the parameters. 838 set(tensor=tensor_obj, value=[Pzz, Pxxyy, Pxy, Pxz, Pyz], param=['Pzz', 'Pxxyy', 'Pxy', 'Pxz', 'Pyz'], errors=errors) 839 840 # Unknown parameter combination. 841 else: 842 raise RelaxUnknownParamCombError('param_types', param_types)
843 844
845 -def map_bounds(param):
846 """The function for creating bounds for the mapping function.""" 847 848 # {Axx, Ayy, Azz, Axxyy, Axy, Axz, Ayz}. 849 if param in ['Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz']: 850 return [-50, 50] 851 852 # alpha. 853 elif param == 'alpha': 854 return [0, 2*pi] 855 856 # beta. 857 elif param == 'beta': 858 return [0, pi] 859 860 # gamma. 861 elif param == 'gamma': 862 return [0, 2*pi]
863 864
865 -def kappa(nuc1='15N', nuc2='1H'):
866 """Function for calculating the kappa constant. 867 868 The kappa constant is:: 869 870 kappa = -3/(8pi^2).gI.gS.mu0.h_bar, 871 872 where gI and gS are the gyromagnetic ratios of the I and S spins, mu0 is the permeability of 873 free space, and h_bar is Planck's constant divided by 2pi. 874 875 @param nuc1: The first nucleus type. 876 @type nuc1: str 877 @param nuc2: The first nucleus type. 878 @type nuc2: str 879 @return: The kappa constant value. 880 @rtype: float 881 """ 882 883 # Gyromagnetic ratios. 884 gI = return_gyromagnetic_ratio(nuc1) 885 gS = return_gyromagnetic_ratio(nuc2) 886 887 # Kappa. 888 return -3.0/(8.0*pi**2) * gI * gS * mu0 * h_bar
889 890
891 -def map_labels(index, params, bounds, swap, inc):
892 """Function for creating labels, tick locations, and tick values for an OpenDX map. 893 894 @param index: The index (which isn't used here?!?). 895 @type index: int 896 @param params: The list of parameter names. 897 @type params: list of str 898 @param bounds: The bounds of the map. 899 @type bounds: list of lists (of a float and bin) 900 @param swap: An array for switching axes around. 901 @type swap: list of int 902 @param inc: The number of increments of one dimension in the map. 903 @type inc: list of int 904 """ 905 906 # Initialise. 907 labels = "{" 908 tick_locations = [] 909 tick_values = [] 910 n = len(params) 911 axis_incs = 5 912 loc_inc = inc / axis_incs 913 914 # Increment over the model parameters. 915 for i in xrange(n): 916 # Parameter conversion factors. 917 factor = return_conversion_factor(params[swap[i]]) 918 919 # Parameter units. 920 units = return_units(params[swap[i]]) 921 922 # Labels. 923 if units: 924 labels = labels + "\"" + params[swap[i]] + " (" + units + ")\"" 925 else: 926 labels = labels + "\"" + params[swap[i]] + "\"" 927 928 # Tick values. 929 vals = bounds[swap[i], 0] / factor 930 val_inc = (bounds[swap[i], 1] - bounds[swap[i], 0]) / (axis_incs * factor) 931 932 if i < n - 1: 933 labels = labels + " " 934 else: 935 labels = labels + "}" 936 937 # Tick locations. 938 string = "{" 939 val = 0.0 940 for j in xrange(axis_incs + 1): 941 string = string + " " + repr(val) 942 val = val + loc_inc 943 string = string + " }" 944 tick_locations.append(string) 945 946 # Tick values. 947 string = "{" 948 for j in xrange(axis_incs + 1): 949 string = string + "\"" + "%.2f" % vals + "\" " 950 vals = vals + val_inc 951 string = string + "}" 952 tick_values.append(string) 953 954 return labels, tick_locations, tick_values
955 956
957 -def matrix_angles(basis_set=0, tensors=None):
958 """Function for calculating the 5D angles between the alignment tensors. 959 960 The basis set used for the 5D vector construction changes the angles calculated. 961 962 @param basis_set: The basis set to use for constructing the 5D vectors. If set to 0, the 963 basis set is {Sxx, Syy, Sxy, Sxz, Syz}. If 1, then the basis set is {Szz, 964 Sxxyy, Sxy, Sxz, Syz}. 965 @type basis_set: int 966 @param tensors: An array of tensors to apply SVD to. If None, all tensors will be used. 967 @type tensors: None or array of str 968 """ 969 970 # Test that alignment tensor data exists. 971 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'): 972 raise RelaxNoTensorError('alignment') 973 974 # Count the number of tensors. 975 tensor_num = 0 976 for tensor in cdp.align_tensors: 977 if tensors and tensor.name not in tensors: 978 continue 979 tensor_num = tensor_num + 1 980 981 # Create the matrix which contains the 5D vectors. 982 matrix = zeros((tensor_num, 5), float64) 983 984 # Loop over the tensors. 985 i = 0 986 for tensor in cdp.align_tensors: 987 # Skip tensors. 988 if tensors and tensor.name not in tensors: 989 continue 990 991 # Unitary basis set. 992 if basis_set == 0: 993 # Pack the elements. 994 matrix[i, 0] = tensor.Sxx 995 matrix[i, 1] = tensor.Syy 996 matrix[i, 2] = tensor.Sxy 997 matrix[i, 3] = tensor.Sxz 998 matrix[i, 4] = tensor.Syz 999 1000 # Geometric basis set. 1001 elif basis_set == 1: 1002 # Pack the elements. 1003 matrix[i, 0] = tensor.Szz 1004 matrix[i, 1] = tensor.Sxxyy 1005 matrix[i, 2] = tensor.Sxy 1006 matrix[i, 3] = tensor.Sxz 1007 matrix[i, 4] = tensor.Syz 1008 1009 # Normalisation. 1010 norm = linalg.norm(matrix[i]) 1011 matrix[i] = matrix[i] / norm 1012 1013 # Increment the index. 1014 i = i + 1 1015 1016 # Initialise the matrix of angles. 1017 cdp.align_tensors.angles = zeros((tensor_num, tensor_num), float64) 1018 1019 # Header print out. 1020 sys.stdout.write("\nData pipe: " + repr(pipes.cdp_name()) + "\n") 1021 sys.stdout.write("\n5D angles in deg between the vectors ") 1022 if basis_set == 0: 1023 sys.stdout.write("{Sxx, Syy, Sxy, Sxz, Syz}") 1024 elif basis_set == 1: 1025 sys.stdout.write("{Szz, Sxx-yy, Sxy, Sxz, Syz}") 1026 sys.stdout.write(":\n") 1027 sys.stdout.write("%8s" % '') 1028 for i in xrange(tensor_num): 1029 sys.stdout.write("%8i" % i) 1030 sys.stdout.write("\n") 1031 1032 # First loop. 1033 for i in xrange(tensor_num): 1034 # Print out. 1035 sys.stdout.write("%8i" % i) 1036 1037 # Second loop. 1038 for j in xrange(tensor_num): 1039 # Dot product. 1040 delta = dot(matrix[i], matrix[j]) 1041 1042 # Check. 1043 if delta > 1: 1044 delta = 1 1045 1046 # The angle (in rad). 1047 cdp.align_tensors.angles[i, j] = arccos(delta) 1048 1049 # Print out the angles in degrees. 1050 sys.stdout.write("%8.1f" % (cdp.align_tensors.angles[i, j]*180.0/pi)) 1051 1052 # Print out. 1053 sys.stdout.write("\n")
1054 1055
1056 -def num_tensors(skip_fixed=True):
1057 """Count the number of tensors. 1058 1059 @keyword skip_fixed: If set to True, then only the tensors without the fixed flag will be counted. If set to False, then all tensors will be counted. 1060 @type skip_fixed: bool 1061 @return: The number of tensors (excluding fixed tensors by default). 1062 @rtype: int 1063 """ 1064 1065 # Init. 1066 count = 0 1067 1068 # Loop over the tensors. 1069 for tensor_cont in cdp.align_tensors: 1070 # Skip fixed tensors. 1071 if skip_fixed and tensor_cont.fixed: 1072 continue 1073 1074 # Increment. 1075 count += 1 1076 1077 # Return the count. 1078 return count
1079 1080
1081 -def reduction(full_tensor=None, red_tensor=None):
1082 """Specify which tensor is a reduction of which other tensor. 1083 1084 @param full_tensor: The full alignment tensor. 1085 @type full_tensor: str 1086 @param red_tensor: The reduced alignment tensor. 1087 @type red_tensor: str 1088 """ 1089 1090 # Test if the current data pipe exists. 1091 pipes.test() 1092 1093 # Test if alignment tensor data exists. 1094 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'): 1095 raise RelaxNoTensorError('alignment') 1096 1097 # Tensor information. 1098 match_full = False 1099 match_red = False 1100 i = 0 1101 for tensor_cont in cdp.align_tensors: 1102 # Test the tensor names. 1103 if tensor_cont.name == full_tensor: 1104 match_full = True 1105 index_full = i 1106 if tensor_cont.name == red_tensor: 1107 match_red = True 1108 index_red = i 1109 1110 # Increment. 1111 i = i + 1 1112 1113 # No match. 1114 if not match_full: 1115 raise RelaxNoTensorError('alignment', full_tensor) 1116 if not match_red: 1117 raise RelaxNoTensorError('alignment', red_tensor) 1118 1119 # Store. 1120 if not hasattr(cdp.align_tensors, 'reduction'): 1121 cdp.align_tensors.reduction = [] 1122 cdp.align_tensors.reduction.append([index_full, index_red])
1123 1124
1125 -def return_conversion_factor(param):
1126 """Function for returning the factor of conversion between different parameter units. 1127 1128 @param param: The parameter name. 1129 @type param: str 1130 @return: The conversion factor. 1131 @rtype: float 1132 """ 1133 1134 # Get the object name. 1135 object_name = return_data_name(param) 1136 1137 # {Axx, Ayy, Azz, Axxyy, Axy, Axz, Ayz}. 1138 if object_name in ['Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz']: 1139 return 1.0 1140 1141 # Angles. 1142 elif object_name in ['alpha', 'beta', 'gamma']: 1143 return (2.0*pi) / 360.0 1144 1145 # No conversion factor. 1146 else: 1147 return 1.0
1148 1149
1150 -def return_data_name(name):
1151 """Return the parameter name. 1152 1153 @param name: The name of the parameter to return the name of. 1154 @type name: str 1155 @return: The parameter name. 1156 @rtype: str 1157 """ 1158 1159 # Enforce that the name must be a string. 1160 if not isinstance(name, str): 1161 raise RelaxStrError('name', name) 1162 1163 # Sxx. 1164 if search('^[Ss]xx$', name): 1165 return 'Sxx' 1166 1167 # Syy. 1168 if search('^[Ss]yy$', name): 1169 return 'Syy' 1170 1171 # Szz. 1172 if search('^[Ss]zz$', name): 1173 return 'Szz' 1174 1175 # Sxy. 1176 if search('^[Ss]xy$', name): 1177 return 'Sxy' 1178 1179 # Sxz. 1180 if search('^[Ss]xz$', name): 1181 return 'Sxz' 1182 1183 # Syz. 1184 if search('^[Ss]yz$', name): 1185 return 'Syz' 1186 1187 # Sxx-yy. 1188 if search('^[Ss]xxyy$', name): 1189 return 'Sxxyy' 1190 1191 # Axx. 1192 if search('^[Aa]xx$', name): 1193 return 'Axx' 1194 1195 # Ayy. 1196 if search('^[Aa]yy$', name): 1197 return 'Ayy' 1198 1199 # Azz. 1200 if search('^[Aa]zz$', name): 1201 return 'Azz' 1202 1203 # Axy. 1204 if search('^[Aa]xy$', name): 1205 return 'Axy' 1206 1207 # Axz. 1208 if search('^[Aa]xz$', name): 1209 return 'Axz' 1210 1211 # Ayz. 1212 if search('^[Aa]yz$', name): 1213 return 'Ayz' 1214 1215 # Axx-yy. 1216 if search('^[Aa]xxyy$', name): 1217 return 'Axxyy' 1218 1219 # Pxx. 1220 if search('^[Pp]xx$', name): 1221 return 'Pxx' 1222 1223 # Pyy. 1224 if search('^[Pp]yy$', name): 1225 return 'Pyy' 1226 1227 # Pzz. 1228 if search('^[Pp]zz$', name): 1229 return 'Pzz' 1230 1231 # Pxy. 1232 if search('^[Pp]xy$', name): 1233 return 'Pxy' 1234 1235 # Pxz. 1236 if search('^[Pp]xz$', name): 1237 return 'Pxz' 1238 1239 # Pyz. 1240 if search('^[Pp]yz$', name): 1241 return 'Pyz' 1242 1243 # Pxx-yy. 1244 if search('^[Pp]xxyy$', name): 1245 return 'Pxxyy' 1246 1247 # alpha. 1248 if search('^a$', name) or search('alpha', name): 1249 return 'alpha' 1250 1251 # beta. 1252 if search('^b$', name) or search('beta', name): 1253 return 'beta' 1254 1255 # gamma. 1256 if search('^g$', name) or search('gamma', name): 1257 return 'gamma' 1258 1259 # No parameter? 1260 raise RelaxUnknownParamError(name)
1261 1262 # User function documentation. 1263 __return_data_name_prompt_doc__ = """ 1264 Alignment tensor parameter string matching patterns 1265 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1266 1267 ____________________________________________________________________________________________ 1268 | | | | 1269 | Data type | Object name | Patterns | 1270 |________________________________________________________|______________|__________________| 1271 | | | | 1272 | The xx component of the Saupe order matrix - Sxx | 'Sxx' | '^[Sa]xx$' | 1273 | | | | 1274 | The yy component of the Saupe order matrix - Syy | 'Syy' | '^[Sa]yy$' | 1275 | | | | 1276 | The zz component of the Saupe order matrix - Szz | 'Szz' | '^[Sa]zz$' | 1277 | | | | 1278 | The xy component of the Saupe order matrix - Sxy | 'Sxy' | '^[Sa]xy$' | 1279 | | | | 1280 | The xz component of the Saupe order matrix - Sxz | 'Sxz' | '^[Sa]xz$' | 1281 | | | | 1282 | The yz component of the Saupe order matrix - Syz | 'Syz' | '^[Sa]yz$' | 1283 | | | | 1284 | The xx-yy component of the Saupe order matrix - Sxx-yy | 'Sxxyy' | '^[Sa]xxyy$' | 1285 | | | | 1286 | The xx component of the alignment tensor - Axx | 'Axx' | '^[Aa]xx$' | 1287 | | | | 1288 | The yy component of the alignment tensor - Ayy | 'Ayy' | '^[Aa]yy$' | 1289 | | | | 1290 | The zz component of the alignment tensor - Azz | 'Azz' | '^[Aa]zz$' | 1291 | | | | 1292 | The xy component of the alignment tensor - Axy | 'Axy' | '^[Aa]xy$' | 1293 | | | | 1294 | The xz component of the alignment tensor - Axz | 'Axz' | '^[Aa]xz$' | 1295 | | | | 1296 | The yz component of the alignment tensor - Ayz | 'Ayz' | '^[Aa]yz$' | 1297 | | | | 1298 | The xx-yy component of the alignment tensor - Axx-yy | 'Axxyy' | '^[Aa]xxyy$' | 1299 | | | | 1300 | The xx component of the probability matrix - Pxx | 'Pxx' | '^[Pa]xx$' | 1301 | | | | 1302 | The yy component of the probability matrix - Pyy | 'Pyy' | '^[Pa]yy$' | 1303 | | | | 1304 | The zz component of the probability matrix - Pzz | 'Pzz' | '^[Pa]zz$' | 1305 | | | | 1306 | The xy component of the probability matrix - Pxy | 'Pxy' | '^[Pa]xy$' | 1307 | | | | 1308 | The xz component of the probability matrix - Pxz | 'Pxz' | '^[Pa]xz$' | 1309 | | | | 1310 | The yz component of the probability matrix - Pyz | 'Pyz' | '^[Pa]yz$' | 1311 | | | | 1312 | The xx-yy component of the probability matrix - Pxx-yy | 'Pxxyy' | '^[Pa]xxyy$' | 1313 | | | | 1314 | The first Euler angle of the alignment tensor - alpha | 'alpha' | '^a$' or 'alpha' | 1315 | | | | 1316 | The second Euler angle of the alignment tensor - beta | 'beta' | '^b$' or 'beta' | 1317 | | | | 1318 | The third Euler angle of the alignment tensor - gamma | 'gamma' | '^g$' or 'gamma' | 1319 |________________________________________________________|______________|__________________| 1320 """ 1321 1322
1323 -def return_tensor(index, skip_fixed=True):
1324 """Return the tensor container for the given index, skipping fixed tensors if required. 1325 1326 @param index: The index of the tensor (if skip_fixed is True, then fixed tensors are not included in the index count). 1327 @type index: int 1328 @keyword skip_fixed: A flag which if True will exclude fixed tensors from the indexation. 1329 @type skip_fixed: bool 1330 @return: The tensor corresponding to the index. 1331 @rtype: data.align_tensor.AlignTensorData instance 1332 """ 1333 1334 # Init. 1335 count = 0 1336 1337 # Loop over the tensors. 1338 for tensor_cont in cdp.align_tensors: 1339 # Skip fixed tensors. 1340 if skip_fixed and tensor_cont.fixed: 1341 continue 1342 1343 # Index match, so return the container. 1344 if index == count: 1345 return tensor_cont 1346 1347 # Increment. 1348 count += 1 1349 1350 # Return False if the container was not found. 1351 return False
1352 1353
1354 -def return_units(param):
1355 """Function for returning a string representing the parameters units. 1356 1357 @param param: The parameter name. 1358 @type param: str 1359 @return: The string representation of the units. 1360 @rtype: str 1361 """ 1362 1363 # Get the object name. 1364 object_name = return_data_name(param) 1365 1366 # {Axx, Ayy, Azz, Axxyy, Axy, Axz, Ayz}. 1367 if object_name in ['Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz']: 1368 return 'Hz' 1369 1370 # Angles. 1371 elif object_name in ['alpha', 'beta', 'gamma']: 1372 return 'deg'
1373 1374
1375 -def set(tensor=None, value=None, param=None, errors=False):
1376 """Set the tensor. 1377 1378 @keyword tensor: The alignment tensor object. 1379 @type tensor: AlignTensorData instance 1380 @keyword value: The list of values to set the parameters to. 1381 @type value: list of float 1382 @keyword param: The list of parameter names. 1383 @type param: list of str 1384 @keyword errors: A flag which determines if the alignment tensor data or its errors are being 1385 input. 1386 @type errors: bool 1387 """ 1388 1389 # Initialise. 1390 geo_params = [] 1391 geo_values = [] 1392 orient_params = [] 1393 orient_values = [] 1394 1395 # Loop over the parameters. 1396 for i in xrange(len(param)): 1397 # Get the object name. 1398 param[i] = return_data_name(param[i]) 1399 1400 # Unknown parameter. 1401 if param[i] == None: 1402 raise RelaxUnknownParamError("alignment tensor", 'None') 1403 1404 # Default value. 1405 if value[i] == None: 1406 value[i] = default_value(object_names[i]) 1407 1408 # Geometric parameter. 1409 if param[i] in ['Sxx', 'Syy', 'Szz', 'Sxxyy', 'Sxy', 'Sxz', 'Syz', 'Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz', 'Pxx', 'Pyy', 'Pzz', 'Pxxyy', 'Pxy', 'Pxz', 'Pyz']: 1410 geo_params.append(param[i]) 1411 geo_values.append(value[i]) 1412 1413 # Orientational parameter. 1414 if param[i] in ['alpha', 'beta', 'gamma']: 1415 orient_params.append(param[i]) 1416 orient_values.append(value[i]) 1417 1418 1419 # Geometric parameters. 1420 ####################### 1421 1422 # A single geometric parameter. 1423 if len(geo_params) == 1: 1424 # Saupe order matrix. 1425 ##################### 1426 1427 # The single parameter Sxx. 1428 if geo_params[0] == 'Sxx': 1429 if errors: 1430 tensor.Sxx_err = geo_values[0] 1431 else: 1432 tensor.Sxx = geo_values[0] 1433 1434 # The single parameter Syy. 1435 elif geo_params[0] == 'Syy': 1436 if errors: 1437 tensor.Syy_err = geo_values[0] 1438 else: 1439 tensor.Syy = geo_values[0] 1440 1441 # The single parameter Sxy. 1442 elif geo_params[0] == 'Sxy': 1443 if errors: 1444 tensor.Sxy_err = geo_values[0] 1445 else: 1446 tensor.Sxy = geo_values[0] 1447 1448 # The single parameter Sxz. 1449 elif geo_params[0] == 'Sxz': 1450 if errors: 1451 tensor.Sxz_err = geo_values[0] 1452 else: 1453 tensor.Sxz = geo_values[0] 1454 1455 # The single parameter Syz. 1456 elif geo_params[0] == 'Syz': 1457 if errors: 1458 tensor.Syz_err = geo_values[0] 1459 else: 1460 tensor.Syz = geo_values[0] 1461 1462 1463 # Alignment tensor. 1464 ################### 1465 1466 # The single parameter Axx. 1467 elif geo_params[0] == 'Axx': 1468 if errors: 1469 tensor.Sxx_err = 3.0/2.0 * geo_values[0] 1470 else: 1471 tensor.Sxx = 3.0/2.0 * geo_values[0] 1472 1473 # The single parameter Ayy. 1474 elif geo_params[0] == 'Ayy': 1475 if errors: 1476 tensor.Syy_err = 3.0/2.0 * geo_values[0] 1477 else: 1478 tensor.Syy = 3.0/2.0 * geo_values[0] 1479 1480 # The single parameter Axy. 1481 elif geo_params[0] == 'Axy': 1482 if errors: 1483 tensor.Sxy_err = 3.0/2.0 * geo_values[0] 1484 else: 1485 tensor.Sxy = 3.0/2.0 * geo_values[0] 1486 1487 # The single parameter Axz. 1488 elif geo_params[0] == 'Axz': 1489 if errors: 1490 tensor.Sxz_err = 3.0/2.0 * geo_values[0] 1491 else: 1492 tensor.Sxz = 3.0/2.0 * geo_values[0] 1493 1494 # The single parameter Ayz. 1495 elif geo_params[0] == 'Ayz': 1496 if errors: 1497 tensor.Syz_err = 3.0/2.0 * geo_values[0] 1498 else: 1499 tensor.Syz = 3.0/2.0 * geo_values[0] 1500 1501 1502 # Probability tensor. 1503 ##################### 1504 1505 # The single parameter Pxx. 1506 elif geo_params[0] == 'Pxx': 1507 if errors: 1508 tensor.Sxx_err = 3.0/2.0 * (geo_values[0] - 1.0/3.0) 1509 else: 1510 tensor.Sxx = 3.0/2.0 * (geo_values[0] - 1.0/3.0) 1511 1512 # The single parameter Pyy. 1513 elif geo_params[0] == 'Pyy': 1514 if errors: 1515 tensor.Syy_err = 3.0/2.0 * (geo_values[0] - 1.0/3.0) 1516 else: 1517 tensor.Syy = 3.0/2.0 * (geo_values[0] - 1.0/3.0) 1518 1519 # The single parameter Pxy. 1520 elif geo_params[0] == 'Pxy': 1521 if errors: 1522 tensor.Sxy_err = 3.0/2.0 * geo_values[0] 1523 else: 1524 tensor.Sxy = 3.0/2.0 * geo_values[0] 1525 1526 # The single parameter Pxz. 1527 elif geo_params[0] == 'Pxz': 1528 if errors: 1529 tensor.Sxz_err = 3.0/2.0 * geo_values[0] 1530 else: 1531 tensor.Sxz = 3.0/2.0 * geo_values[0] 1532 1533 # The single parameter Pyz. 1534 elif geo_params[0] == 'Pyz': 1535 if errors: 1536 tensor.Syz_err = 3.0/2.0 * geo_values[0] 1537 else: 1538 tensor.Syz = 3.0/2.0 * geo_values[0] 1539 1540 # Cannot set the single parameter. 1541 else: 1542 raise RelaxError("The geometric alignment parameter " + repr(geo_params[0]) + " cannot be set.") 1543 1544 # 5 geometric parameters. 1545 elif len(geo_params) == 5: 1546 # The geometric parameter set {Sxx, Syy, Sxy, Sxz, Syz}. 1547 if geo_params.count('Sxx') == 1 and geo_params.count('Syy') == 1 and geo_params.count('Sxy') == 1 and geo_params.count('Sxz') == 1 and geo_params.count('Syz') == 1: 1548 # The parameters. 1549 Sxx = geo_values[geo_params.index('Sxx')] 1550 Syy = geo_values[geo_params.index('Syy')] 1551 Sxy = geo_values[geo_params.index('Sxy')] 1552 Sxz = geo_values[geo_params.index('Sxz')] 1553 Syz = geo_values[geo_params.index('Syz')] 1554 1555 # Set the internal parameter values. 1556 if errors: 1557 tensor.Axx_err = 2.0/3.0 * Sxx 1558 tensor.Ayy_err = 2.0/3.0 * Syy 1559 tensor.Axy_err = 2.0/3.0 * Sxy 1560 tensor.Axz_err = 2.0/3.0 * Sxz 1561 tensor.Ayz_err = 2.0/3.0 * Syz 1562 else: 1563 tensor.Axx = 2.0/3.0 * Sxx 1564 tensor.Ayy = 2.0/3.0 * Syy 1565 tensor.Axy = 2.0/3.0 * Sxy 1566 tensor.Axz = 2.0/3.0 * Sxz 1567 tensor.Ayz = 2.0/3.0 * Syz 1568 1569 # The geometric parameter set {Szz, Sxxyy, Sxy, Sxz, Syz}. 1570 elif geo_params.count('Szz') == 1 and geo_params.count('Sxxyy') == 1 and geo_params.count('Sxy') == 1 and geo_params.count('Sxz') == 1 and geo_params.count('Syz') == 1: 1571 # The parameters. 1572 Szz = geo_values[geo_params.index('Szz')] 1573 Sxxyy = geo_values[geo_params.index('Sxxyy')] 1574 Sxy = geo_values[geo_params.index('Sxy')] 1575 Sxz = geo_values[geo_params.index('Sxz')] 1576 Syz = geo_values[geo_params.index('Syz')] 1577 1578 # Set the internal parameter values. 1579 if errors: 1580 tensor.Axx_err = 2.0/3.0 * -0.5*(Szz-Sxxyy) 1581 tensor.Ayy_err = 2.0/3.0 * -0.5*(Szz+Sxxyy) 1582 tensor.Axy_err = 2.0/3.0 * Sxy 1583 tensor.Axz_err = 2.0/3.0 * Sxz 1584 tensor.Ayz_err = 2.0/3.0 * Syz 1585 else: 1586 tensor.Axx = 2.0/3.0 * -0.5*(Szz-Sxxyy) 1587 tensor.Ayy = 2.0/3.0 * -0.5*(Szz+Sxxyy) 1588 tensor.Axy = 2.0/3.0 * Sxy 1589 tensor.Axz = 2.0/3.0 * Sxz 1590 tensor.Ayz = 2.0/3.0 * Syz 1591 1592 # The geometric parameter set {Axx, Ayy, Axy, Axz, Ayz}. 1593 elif geo_params.count('Axx') == 1 and geo_params.count('Ayy') == 1 and geo_params.count('Axy') == 1 and geo_params.count('Axz') == 1 and geo_params.count('Ayz') == 1: 1594 # The parameters. 1595 Axx = geo_values[geo_params.index('Axx')] 1596 Ayy = geo_values[geo_params.index('Ayy')] 1597 Axy = geo_values[geo_params.index('Axy')] 1598 Axz = geo_values[geo_params.index('Axz')] 1599 Ayz = geo_values[geo_params.index('Ayz')] 1600 1601 # Set the internal parameter values. 1602 if errors: 1603 tensor.Axx_err = Axx 1604 tensor.Ayy_err = Ayy 1605 tensor.Axy_err = Axy 1606 tensor.Axz_err = Axz 1607 tensor.Ayz_err = Ayz 1608 else: 1609 tensor.Axx = Axx 1610 tensor.Ayy = Ayy 1611 tensor.Axy = Axy 1612 tensor.Axz = Axz 1613 tensor.Ayz = Ayz 1614 1615 # The geometric parameter set {Azz, Axxyy, Axy, Axz, Ayz}. 1616 elif geo_params.count('Azz') == 1 and geo_params.count('Axxyy') == 1 and geo_params.count('Axy') == 1 and geo_params.count('Axz') == 1 and geo_params.count('Ayz') == 1: 1617 # The parameters. 1618 Azz = geo_values[geo_params.index('Azz')] 1619 Axxyy = geo_values[geo_params.index('Axxyy')] 1620 Axy = geo_values[geo_params.index('Axy')] 1621 Axz = geo_values[geo_params.index('Axz')] 1622 Ayz = geo_values[geo_params.index('Ayz')] 1623 1624 # Set the internal parameter values. 1625 if errors: 1626 tensor.Axx_err = -0.5*(Azz-Axxyy) 1627 tensor.Ayy_err = -0.5*(Azz+Axxyy) 1628 tensor.Axy_err = Axy 1629 tensor.Axz_err = Axz 1630 tensor.Ayz_err = Ayz 1631 else: 1632 tensor.Axx = -0.5*(Azz-Axxyy) 1633 tensor.Ayy = -0.5*(Azz+Axxyy) 1634 tensor.Axy = Axy 1635 tensor.Axz = Axz 1636 tensor.Ayz = Ayz 1637 1638 # The geometric parameter set {Pxx, Pyy, Pxy, Pxz, Pyz}. 1639 elif geo_params.count('Pxx') == 1 and geo_params.count('Pyy') == 1 and geo_params.count('Pxy') == 1 and geo_params.count('Pxz') == 1 and geo_params.count('Pyz') == 1: 1640 # The parameters. 1641 Pxx = geo_values[geo_params.index('Pxx')] 1642 Pyy = geo_values[geo_params.index('Pyy')] 1643 Pxy = geo_values[geo_params.index('Pxy')] 1644 Pxz = geo_values[geo_params.index('Pxz')] 1645 Pyz = geo_values[geo_params.index('Pyz')] 1646 1647 # Set the internal parameter values. 1648 if errors: 1649 tensor.Axx_err = Pxx - 1.0/3.0 1650 tensor.Ayy_err = Pyy - 1.0/3.0 1651 tensor.Axy_err = Pxy 1652 tensor.Axz_err = Pxz 1653 tensor.Ayz_err = Pyz 1654 else: 1655 tensor.Axx = Pxx - 1.0/3.0 1656 tensor.Ayy = Pyy - 1.0/3.0 1657 tensor.Axy = Pxy 1658 tensor.Axz = Pxz 1659 tensor.Ayz = Pyz 1660 1661 # The geometric parameter set {Pzz, Pxxyy, Pxy, Pxz, Pyz}. 1662 elif geo_params.count('Pzz') == 1 and geo_params.count('Pxxyy') == 1 and geo_params.count('Pxy') == 1 and geo_params.count('Pxz') == 1 and geo_params.count('Pyz') == 1: 1663 # The parameters. 1664 Pzz = geo_values[geo_params.index('Pzz')] 1665 Pxxyy = geo_values[geo_params.index('Pxxyy')] 1666 Pxy = geo_values[geo_params.index('Pxy')] 1667 Pxz = geo_values[geo_params.index('Pxz')] 1668 Pyz = geo_values[geo_params.index('Pyz')] 1669 1670 # Set the internal parameter values. 1671 if errors: 1672 tensor.Axx_err = -0.5*(Pzz-Pxxyy) - 1.0/3.0 1673 tensor.Ayy_err = -0.5*(Pzz+Pxxyy) - 1.0/3.0 1674 tensor.Axy_err = Pxy 1675 tensor.Axz_err = Pxz 1676 tensor.Ayz_err = Pyz 1677 else: 1678 tensor.Axx = -0.5*(Pzz-Pxxyy) - 1.0/3.0 1679 tensor.Ayy = -0.5*(Pzz+Pxxyy) - 1.0/3.0 1680 tensor.Axy = Pxy 1681 tensor.Axz = Pxz 1682 tensor.Ayz = Pyz 1683 1684 # Unknown parameter combination. 1685 else: 1686 raise RelaxUnknownParamCombError('geometric parameter set', geo_params) 1687 1688 1689 # Unknown geometric parameters. 1690 else: 1691 raise RelaxUnknownParamCombError('geometric parameter set', geo_params) 1692 1693 1694 # Orientational parameters. 1695 ########################### 1696 1697 # A single orientational parameter. 1698 if len(orient_params) == 1: 1699 # The single parameter alpha. 1700 if orient_params[0] == 'alpha': 1701 if errors: 1702 tensor.alpha_err = orient_values[orient_params.index('alpha')] 1703 else: 1704 tensor.alpha = orient_values[orient_params.index('alpha')] 1705 1706 # The single parameter beta. 1707 elif orient_params[0] == 'beta': 1708 if errors: 1709 tensor.beta_err = orient_values[orient_params.index('beta')] 1710 else: 1711 tensor.beta = orient_values[orient_params.index('beta')] 1712 1713 # The single parameter gamma. 1714 elif orient_params[0] == 'gamma': 1715 if errors: 1716 tensor.gamma_err = orient_values[orient_params.index('gamma')] 1717 else: 1718 tensor.gamma = orient_values[orient_params.index('gamma')] 1719 1720 # Two orientational parameters. 1721 elif len(orient_params) == 2: 1722 # The orientational parameter set {alpha, beta}. 1723 if orient_params.count('alpha') == 1 and orient_params.count('beta') == 1: 1724 if errors: 1725 tensor.alpha_err = orient_values[orient_params.index('alpha')] 1726 tensor.beta_err = orient_values[orient_params.index('beta')] 1727 else: 1728 tensor.alpha = orient_values[orient_params.index('alpha')] 1729 tensor.beta = orient_values[orient_params.index('beta')] 1730 1731 # The orientational parameter set {alpha, gamma}. 1732 if orient_params.count('alpha') == 1 and orient_params.count('gamma') == 1: 1733 if errors: 1734 tensor.alpha_err = orient_values[orient_params.index('alpha')] 1735 tensor.gamma_err = orient_values[orient_params.index('gamma')] 1736 else: 1737 tensor.alpha = orient_values[orient_params.index('alpha')] 1738 tensor.gamma = orient_values[orient_params.index('gamma')] 1739 1740 # The orientational parameter set {beta, gamma}. 1741 if orient_params.count('beta') == 1 and orient_params.count('gamma') == 1: 1742 if errors: 1743 tensor.beta_err = orient_values[orient_params.index('beta')] 1744 tensor.gamma_err = orient_values[orient_params.index('gamma')] 1745 else: 1746 tensor.beta = orient_values[orient_params.index('beta')] 1747 tensor.gamma = orient_values[orient_params.index('gamma')] 1748 1749 # Unknown parameter combination. 1750 else: 1751 raise RelaxUnknownParamCombError('orientational parameter set', orient_params) 1752 1753 # Three orientational parameters. 1754 elif len(orient_params) == 3: 1755 # The orientational parameter set {alpha, beta, gamma}. 1756 if orient_params.count('alpha') == 1 and orient_params.count('beta') == 1: 1757 if errors: 1758 tensor.alpha_err = orient_values[orient_params.index('alpha')] 1759 tensor.beta_err = orient_values[orient_params.index('beta')] 1760 tensor.gamma_err = orient_values[orient_params.index('gamma')] 1761 else: 1762 tensor.alpha = orient_values[orient_params.index('alpha')] 1763 tensor.beta = orient_values[orient_params.index('beta')] 1764 tensor.gamma = orient_values[orient_params.index('gamma')] 1765 1766 # Unknown parameter combination. 1767 else: 1768 raise RelaxUnknownParamCombError('orientational parameter set', orient_params) 1769 1770 # More than three orientational parameters. 1771 elif len(orient_params) > 3: 1772 raise RelaxUnknownParamCombError('orientational parameter set', orient_params) 1773 1774 1775 # Fold the angles in. 1776 ##################### 1777 1778 if orient_params: 1779 fold_angles()
1780 1781 # User function documentation. 1782 __set_prompt_doc__ = """ 1783 Alignment tensor set details 1784 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 1785 1786 If the alignment tensor has not been setup, use the more powerful function 1787 'alignment_tensor.init' to initialise the tensor parameters. 1788 1789 The alignment tensor parameters can only be set when the data pipe corresponds to model-free 1790 analysis. The units of the parameters are: 1791 1792 Unitless for Sxx, Syy, Szz, Sxxyy, Sxy, Sxz, Syz. 1793 Unitless for Axx, Ayy, Azz, Axxyy, Axy, Axz, Ayz. 1794 Unitless for Pxx, Pyy, Pzz, Pxxyy, Pxy, Pxz, Pyz. 1795 Radians for all angles (alpha, beta, gamma). 1796 1797 If a single geometric parameter is supplied, it must be one of Bxx, Byy, Bxy, Bxz, Byz, where B 1798 is one of S, A, or P. For the parameters Bzz and Bxxyy, it is not possible to determine how to 1799 use the currently set values together with the supplied value to calculate the new internal 1800 parameters. When supplying multiple geometric parameters, the set must belong to one of 1801 1802 {Sxx, Syy, Sxy, Sxz, Syz}, 1803 {Szz, Sxxyy, Sxy, Sxz, Syz}. 1804 {Axx, Ayy, Axy, Axz, Ayz}, 1805 {Azz, Axxyy, Axy, Axz, Ayz}. 1806 {Pxx, Pyy, Pxy, Pxz, Pyz}, 1807 {Pzz, Pxxyy, Pxy, Pxz, Pyz}. 1808 """ 1809 1810
1811 -def set_domain(tensor=None, domain=None):
1812 """Set the domain label for the given tensor. 1813 1814 @param tensor: The alignment tensor label. 1815 @type tensor: str 1816 @param domain: The domain label. 1817 @type domain: str 1818 """ 1819 1820 # Test if alignment tensor data exists. 1821 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'): 1822 raise RelaxNoTensorError('alignment') 1823 1824 # Loop over the tensors. 1825 match = False 1826 for tensor_cont in cdp.align_tensors: 1827 # Find the matching tensor and then store the domain label. 1828 if tensor_cont.name == tensor: 1829 tensor_cont.domain = domain 1830 match = True 1831 1832 # The tensor label doesn't exist. 1833 if not match: 1834 raise RelaxNoTensorError('alignment', tensor)
1835 1836
1837 -def svd(basis_set=0, tensors=None):
1838 """Function for calculating the singular values of all the loaded tensors. 1839 1840 The matrix on which SVD will be performed is:: 1841 1842 | Sxx1 Syy1 Sxy1 Sxz1 Syz1 | 1843 | Sxx2 Syy2 Sxy2 Sxz2 Syz2 | 1844 | Sxx3 Syy3 Sxy3 Sxz3 Syz3 | 1845 | . . . . . | 1846 | . . . . . | 1847 | . . . . . | 1848 | SxxN SyyN SxyN SxzN SyzN | 1849 1850 This is the default unitary basis set (selected when basis_set is 0). Alternatively a geometric 1851 basis set consisting of the stretching and skewing parameters Szz and Sxx-yy respectively 1852 replacing Sxx and Syy can be chosen by setting basis_set to 1. The matrix in this case is:: 1853 1854 | Szz1 Sxxyy1 Sxy1 Sxz1 Syz1 | 1855 | Szz2 Sxxyy2 Sxy2 Sxz2 Syz2 | 1856 | Szz3 Sxxyy3 Sxy3 Sxz3 Syz3 | 1857 | . . . . . | 1858 | . . . . . | 1859 | . . . . . | 1860 | SzzN SxxyyN SxyN SxzN SyzN | 1861 1862 The relationships between the geometric and unitary basis sets are:: 1863 1864 Szz = - Sxx - Syy, 1865 Sxxyy = Sxx - Syy, 1866 1867 The SVD values and condition number are dependendent upon the basis set chosen. 1868 1869 1870 @param basis_set: The basis set to create the 5 by n matrix on which to perform SVD. 1871 @type basis_set: int 1872 @param tensors: An array of tensors to apply SVD to. If None, all tensors will be used. 1873 @type tensors: None or array of str 1874 """ 1875 1876 # Test that alignment tensor data exists. 1877 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'): 1878 raise RelaxNoTensorError('alignment') 1879 1880 # Count the number of tensors used in the SVD. 1881 tensor_num = 0 1882 for tensor in cdp.align_tensors: 1883 if tensors and tensor.name not in tensors: 1884 continue 1885 tensor_num = tensor_num + 1 1886 1887 # Create the matrix to apply SVD on. 1888 matrix = zeros((tensor_num, 5), float64) 1889 1890 # Pack the elements. 1891 i = 0 1892 for tensor in cdp.align_tensors: 1893 # Skip tensors. 1894 if tensors and tensor.name not in tensors: 1895 continue 1896 1897 # Unitary basis set. 1898 if basis_set == 0: 1899 matrix[i, 0] = tensor.Sxx 1900 matrix[i, 1] = tensor.Syy 1901 matrix[i, 2] = tensor.Sxy 1902 matrix[i, 3] = tensor.Sxz 1903 matrix[i, 4] = tensor.Syz 1904 1905 # Geometric basis set. 1906 elif basis_set == 1: 1907 matrix[i, 0] = tensor.Szz 1908 matrix[i, 1] = tensor.Sxxyy 1909 matrix[i, 2] = tensor.Sxy 1910 matrix[i, 3] = tensor.Sxz 1911 matrix[i, 4] = tensor.Syz 1912 1913 # Increment the index. 1914 i = i + 1 1915 1916 # SVD. 1917 u, s, vh = linalg.svd(matrix) 1918 1919 # Store the singular values. 1920 cdp.align_tensors.singular_vals = s 1921 1922 # Calculate and store the condition number. 1923 cdp.align_tensors.cond_num = s[0] / s[-1] 1924 1925 # Print out. 1926 print(("\nData pipe: " + repr(pipes.cdp_name()))) 1927 print("\nSingular values:") 1928 for val in s: 1929 print(("\t" + repr(val))) 1930 print(("\nCondition number: " + repr(cdp.align_tensors.cond_num)))
1931