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

Source Code for Module pipe_control.align_tensor

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