pipe_control.align

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 acos, 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 if cdp.align_tensors[i].name == None: 1059 table[0].append(repr(i)) 1060 else: 1061 table[0].append(cdp.align_tensors[i].name) 1062 1063 # First loop over the rows. 1064 for i in range(tensor_num): 1065 # Add the tensor name. 1066 if cdp.align_tensors[i].name == None: 1067 table.append([repr(i)]) 1068 else: 1069 table.append([cdp.align_tensors[i].name]) 1070 1071 # Second loop over the columns. 1072 for j in range(tensor_num): 1073 # The vector angles. 1074 if basis_set in ['unitary 9D', 'unitary 5D', 'geometric 5D']: 1075 # Dot product. 1076 delta = dot(matrix[i], matrix[j]) 1077 1078 # Check. 1079 if delta > 1: 1080 delta = 1 1081 1082 # The angle. 1083 theta = arccos(delta) 1084 1085 # The irreducible complex conjugate angles. 1086 if basis_set in ['irreducible 5D']: 1087 theta = vector_angle_complex_conjugate(v1=matrix[i], v2=matrix[j], v1_conj=matrix_conj[i], v2_conj=matrix_conj[j]) 1088 1089 # The full matrix angle. 1090 elif basis_set in ['matrix']: 1091 # The Euclidean inner product. 1092 nom = inner(matrix[i].flatten(), matrix[j].flatten()) 1093 1094 # The Frobenius norms. 1095 denom = norm(matrix[i]) * norm(matrix[j]) 1096 1097 # The angle. 1098 ratio = nom / denom 1099 if ratio <= 1.0: 1100 theta = arccos(nom / denom) 1101 else: 1102 theta = 0.0 1103 1104 # Store the angle (in rad). 1105 cdp.align_tensors.angles[i, j] = theta 1106 1107 # Add to the table as degrees. 1108 angle = cdp.align_tensors.angles[i, j] 1109 if angle_units == 'deg': 1110 angle = angle * 180.0 / pi 1111 format = "%" + repr(7+precision) + "." + repr(precision) + "f" 1112 table[i+1].append(format % angle) 1113 1114 # Write out the table. 1115 write_data(out=sys.stdout, data=table)

1116 1117

1118 -def num_tensors(skip_fixed=True):

1119 """Count the number of tensors. 1120 1121 @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. 1122 @type skip_fixed: bool 1123 @return: The number of tensors (excluding fixed tensors by default). 1124 @rtype: int 1125 """ 1126 1127 # Init. 1128 count = 0 1129 1130 # Loop over the tensors. 1131 for tensor_cont in cdp.align_tensors: 1132 # Skip fixed tensors. 1133 if skip_fixed and tensor_cont.fixed: 1134 continue 1135 1136 # Increment. 1137 count += 1 1138 1139 # Return the count. 1140 return count

1141 1142

1143 -def opt_uses_align_data(align_id=None):

1144 """Determine if the PCS or RDC data for the given alignment ID is needed for optimisation. 1145 1146 @keyword align_id: The optional alignment ID string. 1147 @type align_id: str 1148 @return: True if alignment data is to be used for optimisation, False otherwise. 1149 @rtype: bool 1150 """ 1151 1152 # No alignment IDs. 1153 if not hasattr(cdp, 'align_ids'): 1154 return False 1155 1156 # Convert the align IDs to an array, or take all IDs. 1157 if align_id: 1158 align_ids = [align_id] 1159 else: 1160 align_ids = cdp.align_ids 1161 1162 # Check the PCS and RDC. 1163 for align_id in align_ids: 1164 if pipe_control.pcs.opt_uses_pcs(align_id) or pipe_control.rdc.opt_uses_rdc(align_id): 1165 return True 1166 1167 # No alignment data is used for optimisation. 1168 return False

1169 1170

1171 -def opt_uses_tensor(tensor):

1172 """Determine if the given tensor is to be optimised. 1173 1174 @param tensor: The alignment tensor to check. 1175 @type tensor: AlignmentTensor object. 1176 @return: True if the tensor is to be optimised, False otherwise. 1177 @rtype: bool 1178 """ 1179 1180 # Combine all RDC and PCS IDs. 1181 ids = [] 1182 if hasattr(cdp, 'rdc_ids'): 1183 ids += cdp.rdc_ids 1184 if hasattr(cdp, 'pcs_ids'): 1185 ids += cdp.pcs_ids 1186 1187 # No RDC or PCS data for the alignment, so skip the tensor as it will not be optimised. 1188 if tensor.align_id not in ids: 1189 return False 1190 1191 # Fixed tensor. 1192 if tensor.fixed: 1193 return False 1194 1195 # The tensor is to be optimised. 1196 return True

1197 1198

1199 -def reduction(full_tensor=None, red_tensor=None):

1200 """Specify which tensor is a reduction of which other tensor. 1201 1202 @param full_tensor: The full alignment tensor. 1203 @type full_tensor: str 1204 @param red_tensor: The reduced alignment tensor. 1205 @type red_tensor: str 1206 """ 1207 1208 # Tensor information. 1209 match_full = False 1210 match_red = False 1211 i = 0 1212 for tensor_cont in cdp.align_tensors: 1213 # Test the tensor names. 1214 if tensor_cont.name == full_tensor: 1215 match_full = True 1216 index_full = i 1217 if tensor_cont.name == red_tensor: 1218 match_red = True 1219 index_red = i 1220 1221 # Increment. 1222 i = i + 1 1223 1224 # No match. 1225 if not match_full: 1226 raise RelaxNoTensorError('alignment', full_tensor) 1227 if not match_red: 1228 raise RelaxNoTensorError('alignment', red_tensor) 1229 1230 # Store. 1231 if not hasattr(cdp.align_tensors, 'reduction'): 1232 cdp.align_tensors.reduction = [] 1233 cdp.align_tensors.reduction.append([index_full, index_red])

1234 1235

1236 -def return_tensor(index, skip_fixed=True):

1237 """Return the tensor container for the given index, skipping fixed tensors if required. 1238 1239 @param index: The index of the tensor (if skip_fixed is True, then fixed tensors are not included in the index count). 1240 @type index: int 1241 @keyword skip_fixed: A flag which if True will exclude fixed tensors from the indexation. 1242 @type skip_fixed: bool 1243 @return: The tensor corresponding to the index. 1244 @rtype: data.align_tensor.AlignTensorData instance 1245 """ 1246 1247 # Init. 1248 count = 0 1249 1250 # Loop over the tensors. 1251 for tensor_cont in cdp.align_tensors: 1252 # Skip fixed tensors. 1253 if skip_fixed and tensor_cont.fixed: 1254 continue 1255 1256 # Index match, so return the container. 1257 if index == count: 1258 return tensor_cont 1259 1260 # Increment. 1261 count += 1 1262 1263 # Return False if the container was not found. 1264 return False

1265 1266

1267 -def set(tensor=None, value=None, param=None, errors=False):

1268 """Set the tensor. 1269 1270 @keyword tensor: The alignment tensor object. 1271 @type tensor: AlignTensorData instance 1272 @keyword value: The list of values to set the parameters to. 1273 @type value: list of float 1274 @keyword param: The list of parameter names. 1275 @type param: list of str 1276 @keyword errors: A flag which determines if the alignment tensor data or its errors are being 1277 input. 1278 @type errors: bool 1279 """ 1280 1281 # Initialise. 1282 geo_params = [] 1283 geo_values = [] 1284 orient_params = [] 1285 orient_values = [] 1286 1287 # Loop over the parameters. 1288 for i in range(len(param)): 1289 # Unknown parameter. 1290 if param[i] == None: 1291 raise RelaxUnknownParamError("alignment tensor", 'None') 1292 1293 # Default value. 1294 if value[i] == None: 1295 value[i] = 0.0 1296 1297 # Geometric parameter. 1298 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']: 1299 geo_params.append(param[i]) 1300 geo_values.append(value[i]) 1301 1302 # Orientational parameter. 1303 if param[i] in ['alpha', 'beta', 'gamma']: 1304 orient_params.append(param[i]) 1305 orient_values.append(value[i]) 1306 1307 1308 # Geometric parameters. 1309 ####################### 1310 1311 # A single geometric parameter. 1312 if len(geo_params) == 1: 1313 # Saupe order matrix. 1314 ##################### 1315 1316 # The single parameter Sxx. 1317 if geo_params[0] == 'Sxx': 1318 if errors: 1319 tensor.set(param='Sxx', value=geo_values[0], category='err') 1320 else: 1321 tensor.set(param='Sxx', value=geo_values[0]) 1322 1323 # The single parameter Syy. 1324 elif geo_params[0] == 'Syy': 1325 if errors: 1326 tensor.set(param='Syy', value=geo_values[0], category='err') 1327 else: 1328 tensor.set(param='Syy', value=geo_values[0]) 1329 1330 # The single parameter Sxy. 1331 elif geo_params[0] == 'Sxy': 1332 if errors: 1333 tensor.set(param='Sxy', value=geo_values[0], category='err') 1334 else: 1335 tensor.set(param='Sxy', value=geo_values[0]) 1336 1337 # The single parameter Sxz. 1338 elif geo_params[0] == 'Sxz': 1339 if errors: 1340 tensor.set(param='Sxz', value=geo_values[0], category='err') 1341 else: 1342 tensor.set(param='Sxz', value=geo_values[0]) 1343 1344 # The single parameter Syz. 1345 elif geo_params[0] == 'Syz': 1346 if errors: 1347 tensor.set(param='Syz', value=geo_values[0], category='err') 1348 else: 1349 tensor.set(param='Syz', value=geo_values[0]) 1350 1351 1352 # Alignment tensor. 1353 ################### 1354 1355 # The single parameter Axx. 1356 elif geo_params[0] == 'Axx': 1357 if errors: 1358 tensor.set(param='Sxx', value=3.0/2.0 * geo_values[0], category='err') 1359 else: 1360 tensor.set(param='Sxx', value=3.0/2.0 * geo_values[0]) 1361 1362 # The single parameter Ayy. 1363 elif geo_params[0] == 'Ayy': 1364 if errors: 1365 tensor.set(param='Syy', value=3.0/2.0 * geo_values[0], category='err') 1366 else: 1367 tensor.set(param='Syy', value=3.0/2.0 * geo_values[0]) 1368 1369 # The single parameter Axy. 1370 elif geo_params[0] == 'Axy': 1371 if errors: 1372 tensor.set(param='Sxy', value=3.0/2.0 * geo_values[0], category='err') 1373 else: 1374 tensor.set(param='Sxy', value=3.0/2.0 * geo_values[0]) 1375 1376 # The single parameter Axz. 1377 elif geo_params[0] == 'Axz': 1378 if errors: 1379 tensor.set(param='Sxz', value=3.0/2.0 * geo_values[0], category='err') 1380 else: 1381 tensor.set(param='Sxz', value=3.0/2.0 * geo_values[0]) 1382 1383 # The single parameter Ayz. 1384 elif geo_params[0] == 'Ayz': 1385 if errors: 1386 tensor.set(param='Syz', value=3.0/2.0 * geo_values[0], category='err') 1387 else: 1388 tensor.set(param='Syz', value=3.0/2.0 * geo_values[0]) 1389 1390 1391 # Probability tensor. 1392 ##################### 1393 1394 # The single parameter Pxx. 1395 elif geo_params[0] == 'Pxx': 1396 if errors: 1397 tensor.set(param='Sxx', value=3.0/2.0 * (geo_values[0] - 1.0/3.0), category='err') 1398 else: 1399 tensor.set(param='Sxx', value=3.0/2.0 * (geo_values[0] - 1.0/3.0)) 1400 1401 # The single parameter Pyy. 1402 elif geo_params[0] == 'Pyy': 1403 if errors: 1404 tensor.set(param='Syy', value=3.0/2.0 * (geo_values[0] - 1.0/3.0), category='err') 1405 else: 1406 tensor.set(param='Syy', value=3.0/2.0 * (geo_values[0] - 1.0/3.0)) 1407 1408 # The single parameter Pxy. 1409 elif geo_params[0] == 'Pxy': 1410 if errors: 1411 tensor.set(param='Sxy', value=3.0/2.0 * geo_values[0], category='err') 1412 else: 1413 tensor.set(param='Sxy', value=3.0/2.0 * geo_values[0]) 1414 1415 # The single parameter Pxz. 1416 elif geo_params[0] == 'Pxz': 1417 if errors: 1418 tensor.set(param='Sxz', value=3.0/2.0 * geo_values[0], category='err') 1419 else: 1420 tensor.set(param='Sxz', value=3.0/2.0 * geo_values[0]) 1421 1422 # The single parameter Pyz. 1423 elif geo_params[0] == 'Pyz': 1424 if errors: 1425 tensor.set(param='Syz', value=3.0/2.0 * geo_values[0], category='err') 1426 else: 1427 tensor.set(param='Syz', value=3.0/2.0 * geo_values[0]) 1428 1429 # Cannot set the single parameter. 1430 else: 1431 raise RelaxError("The geometric alignment parameter " + repr(geo_params[0]) + " cannot be set.") 1432 1433 # 5 geometric parameters. 1434 elif len(geo_params) == 5: 1435 # The geometric parameter set {Sxx, Syy, Sxy, Sxz, Syz}. 1436 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: 1437 # The parameters. 1438 Sxx = geo_values[geo_params.index('Sxx')] 1439 Syy = geo_values[geo_params.index('Syy')] 1440 Sxy = geo_values[geo_params.index('Sxy')] 1441 Sxz = geo_values[geo_params.index('Sxz')] 1442 Syz = geo_values[geo_params.index('Syz')] 1443 1444 # Set the internal parameter values. 1445 if errors: 1446 tensor.set(param='Axx', value=2.0/3.0 * Sxx, category='err') 1447 tensor.set(param='Ayy', value=2.0/3.0 * Syy, category='err') 1448 tensor.set(param='Axy', value=2.0/3.0 * Sxy, category='err') 1449 tensor.set(param='Axz', value=2.0/3.0 * Sxz, category='err') 1450 tensor.set(param='Ayz', value=2.0/3.0 * Syz, category='err') 1451 else: 1452 tensor.set(param='Axx', value=2.0/3.0 * Sxx) 1453 tensor.set(param='Ayy', value=2.0/3.0 * Syy) 1454 tensor.set(param='Axy', value=2.0/3.0 * Sxy) 1455 tensor.set(param='Axz', value=2.0/3.0 * Sxz) 1456 tensor.set(param='Ayz', value=2.0/3.0 * Syz) 1457 1458 # The geometric parameter set {Szz, Sxxyy, Sxy, Sxz, Syz}. 1459 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: 1460 # The parameters. 1461 Szz = geo_values[geo_params.index('Szz')] 1462 Sxxyy = geo_values[geo_params.index('Sxxyy')] 1463 Sxy = geo_values[geo_params.index('Sxy')] 1464 Sxz = geo_values[geo_params.index('Sxz')] 1465 Syz = geo_values[geo_params.index('Syz')] 1466 1467 # Set the internal parameter values. 1468 if errors: 1469 tensor.set(param='Axx', value=2.0/3.0 * -0.5*(Szz-Sxxyy), category='err') 1470 tensor.set(param='Ayy', value=2.0/3.0 * -0.5*(Szz+Sxxyy), category='err') 1471 tensor.set(param='Axy', value=2.0/3.0 * Sxy, category='err') 1472 tensor.set(param='Axz', value=2.0/3.0 * Sxz, category='err') 1473 tensor.set(param='Ayz', value=2.0/3.0 * Syz, category='err') 1474 else: 1475 tensor.set(param='Axx', value=2.0/3.0 * -0.5*(Szz-Sxxyy)) 1476 tensor.set(param='Ayy', value=2.0/3.0 * -0.5*(Szz+Sxxyy)) 1477 tensor.set(param='Axy', value=2.0/3.0 * Sxy) 1478 tensor.set(param='Axz', value=2.0/3.0 * Sxz) 1479 tensor.set(param='Ayz', value=2.0/3.0 * Syz) 1480 1481 # The geometric parameter set {Axx, Ayy, Axy, Axz, Ayz}. 1482 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: 1483 # The parameters. 1484 Axx = geo_values[geo_params.index('Axx')] 1485 Ayy = geo_values[geo_params.index('Ayy')] 1486 Axy = geo_values[geo_params.index('Axy')] 1487 Axz = geo_values[geo_params.index('Axz')] 1488 Ayz = geo_values[geo_params.index('Ayz')] 1489 1490 # Set the internal parameter values. 1491 if errors: 1492 tensor.set(param='Axx', value=Axx, category='err') 1493 tensor.set(param='Ayy', value=Ayy, category='err') 1494 tensor.set(param='Axy', value=Axy, category='err') 1495 tensor.set(param='Axz', value=Axz, category='err') 1496 tensor.set(param='Ayz', value=Ayz, category='err') 1497 else: 1498 tensor.set(param='Axx', value=Axx) 1499 tensor.set(param='Ayy', value=Ayy) 1500 tensor.set(param='Axy', value=Axy) 1501 tensor.set(param='Axz', value=Axz) 1502 tensor.set(param='Ayz', value=Ayz) 1503 1504 # The geometric parameter set {Azz, Axxyy, Axy, Axz, Ayz}. 1505 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: 1506 # The parameters. 1507 Azz = geo_values[geo_params.index('Azz')] 1508 Axxyy = geo_values[geo_params.index('Axxyy')] 1509 Axy = geo_values[geo_params.index('Axy')] 1510 Axz = geo_values[geo_params.index('Axz')] 1511 Ayz = geo_values[geo_params.index('Ayz')] 1512 1513 # Set the internal parameter values. 1514 if errors: 1515 tensor.set(param='Axx', value=-0.5*(Azz-Axxyy), category='err') 1516 tensor.set(param='Ayy', value=-0.5*(Azz+Axxyy), category='err') 1517 tensor.set(param='Axy', value=Axy, category='err') 1518 tensor.set(param='Axz', value=Axz, category='err') 1519 tensor.set(param='Ayz', value=Ayz, category='err') 1520 else: 1521 tensor.set(param='Axx', value=-0.5*(Azz-Axxyy)) 1522 tensor.set(param='Ayy', value=-0.5*(Azz+Axxyy)) 1523 tensor.set(param='Axy', value=Axy) 1524 tensor.set(param='Axz', value=Axz) 1525 tensor.set(param='Ayz', value=Ayz) 1526 1527 # The geometric parameter set {Pxx, Pyy, Pxy, Pxz, Pyz}. 1528 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: 1529 # The parameters. 1530 Pxx = geo_values[geo_params.index('Pxx')] 1531 Pyy = geo_values[geo_params.index('Pyy')] 1532 Pxy = geo_values[geo_params.index('Pxy')] 1533 Pxz = geo_values[geo_params.index('Pxz')] 1534 Pyz = geo_values[geo_params.index('Pyz')] 1535 1536 # Set the internal parameter values. 1537 if errors: 1538 tensor.set(param='Axx', value=Pxx - 1.0/3.0, category='err') 1539 tensor.set(param='Ayy', value=Pyy - 1.0/3.0, category='err') 1540 tensor.set(param='Axy', value=Pxy, category='err') 1541 tensor.set(param='Axz', value=Pxz, category='err') 1542 tensor.set(param='Ayz', value=Pyz, category='err') 1543 else: 1544 tensor.set(param='Axx', value=Pxx - 1.0/3.0) 1545 tensor.set(param='Ayy', value=Pyy - 1.0/3.0) 1546 tensor.set(param='Axy', value=Pxy) 1547 tensor.set(param='Axz', value=Pxz) 1548 tensor.set(param='Ayz', value=Pyz) 1549 1550 # The geometric parameter set {Pzz, Pxxyy, Pxy, Pxz, Pyz}. 1551 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: 1552 # The parameters. 1553 Pzz = geo_values[geo_params.index('Pzz')] 1554 Pxxyy = geo_values[geo_params.index('Pxxyy')] 1555 Pxy = geo_values[geo_params.index('Pxy')] 1556 Pxz = geo_values[geo_params.index('Pxz')] 1557 Pyz = geo_values[geo_params.index('Pyz')] 1558 1559 # Set the internal parameter values. 1560 if errors: 1561 tensor.set(param='Axx', value=-0.5*(Pzz-Pxxyy) - 1.0/3.0, category='err') 1562 tensor.set(param='Ayy', value=-0.5*(Pzz+Pxxyy) - 1.0/3.0, category='err') 1563 tensor.set(param='Axy', value=Pxy, category='err') 1564 tensor.set(param='Axz', value=Pxz, category='err') 1565 tensor.set(param='Ayz', value=Pyz, category='err') 1566 else: 1567 tensor.set(param='Axx', value=-0.5*(Pzz-Pxxyy) - 1.0/3.0) 1568 tensor.set(param='Ayy', value=-0.5*(Pzz+Pxxyy) - 1.0/3.0) 1569 tensor.set(param='Axy', value=Pxy) 1570 tensor.set(param='Axz', value=Pxz) 1571 tensor.set(param='Ayz', value=Pyz) 1572 1573 # Unknown parameter combination. 1574 else: 1575 raise RelaxUnknownParamCombError('geometric parameter set', geo_params) 1576 1577 1578 # Unknown geometric parameters. 1579 else: 1580 raise RelaxUnknownParamCombError('geometric parameter set', geo_params) 1581 1582 1583 # Orientational parameters. 1584 ########################### 1585 1586 # A single orientational parameter. 1587 if len(orient_params) == 1: 1588 # The single parameter alpha. 1589 if orient_params[0] == 'alpha': 1590 if errors: 1591 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')], category='err') 1592 else: 1593 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')]) 1594 1595 # The single parameter beta. 1596 elif orient_params[0] == 'beta': 1597 if errors: 1598 tensor.set(param='beta', value=orient_values[orient_params.index('beta')], category='err') 1599 else: 1600 tensor.set(param='beta', value=orient_values[orient_params.index('beta')]) 1601 1602 # The single parameter gamma. 1603 elif orient_params[0] == 'gamma': 1604 if errors: 1605 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')], category='err') 1606 else: 1607 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')]) 1608 1609 # Two orientational parameters. 1610 elif len(orient_params) == 2: 1611 # The orientational parameter set {alpha, beta}. 1612 if orient_params.count('alpha') == 1 and orient_params.count('beta') == 1: 1613 if errors: 1614 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')], category='err') 1615 tensor.set(param='beta', value=orient_values[orient_params.index('beta')], category='err') 1616 else: 1617 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')]) 1618 tensor.set(param='beta', value=orient_values[orient_params.index('beta')]) 1619 1620 # The orientational parameter set {alpha, gamma}. 1621 if orient_params.count('alpha') == 1 and orient_params.count('gamma') == 1: 1622 if errors: 1623 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')], category='err') 1624 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')], category='err') 1625 else: 1626 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')]) 1627 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')]) 1628 1629 # The orientational parameter set {beta, gamma}. 1630 if orient_params.count('beta') == 1 and orient_params.count('gamma') == 1: 1631 if errors: 1632 tensor.set(param='beta', value=orient_values[orient_params.index('beta')], category='err') 1633 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')], category='err') 1634 else: 1635 tensor.set(param='beta', value=orient_values[orient_params.index('beta')]) 1636 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')]) 1637 1638 # Unknown parameter combination. 1639 else: 1640 raise RelaxUnknownParamCombError('orientational parameter set', orient_params) 1641 1642 # Three orientational parameters. 1643 elif len(orient_params) == 3: 1644 # The orientational parameter set {alpha, beta, gamma}. 1645 if orient_params.count('alpha') == 1 and orient_params.count('beta') == 1: 1646 if errors: 1647 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')], category='err') 1648 tensor.set(param='beta', value=orient_values[orient_params.index('beta')], category='err') 1649 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')], category='err') 1650 else: 1651 tensor.set(param='alpha', value=orient_values[orient_params.index('alpha')]) 1652 tensor.set(param='beta', value=orient_values[orient_params.index('beta')]) 1653 tensor.set(param='gamma', value=orient_values[orient_params.index('gamma')]) 1654 1655 # Unknown parameter combination. 1656 else: 1657 raise RelaxUnknownParamCombError('orientational parameter set', orient_params) 1658 1659 # More than three orientational parameters. 1660 elif len(orient_params) > 3: 1661 raise RelaxUnknownParamCombError('orientational parameter set', orient_params) 1662 1663 1664 # Fold the angles in. 1665 ##################### 1666 1667 if orient_params: 1668 fold_angles()

1669 1670

1671 -def set_align_id(tensor=None, align_id=None):

1672 """Set the align ID string for the given tensor. 1673 1674 @keyword tensor: The alignment tensor label. 1675 @type tensor: str 1676 @keyword align_id: The alignment ID string. 1677 @type align_id: str 1678 """ 1679 1680 # Loop over the tensors. 1681 match = False 1682 for tensor_cont in cdp.align_tensors: 1683 # Find the matching tensor and then store the ID. 1684 if tensor_cont.name == tensor: 1685 tensor_cont.align_id = align_id 1686 match = True 1687 1688 # The tensor label doesn't exist. 1689 if not match: 1690 raise RelaxNoTensorError('alignment', tensor)

1691 1692

1693 -def set_domain(tensor=None, domain=None):

1694 """Set the domain label for the given tensor. 1695 1696 @param tensor: The alignment tensor label. 1697 @type tensor: str 1698 @param domain: The domain label. 1699 @type domain: str 1700 """ 1701 1702 # Check that the domain is defined. 1703 if not hasattr(cdp, 'domain') or domain not in cdp.domain: 1704 raise RelaxError("The domain '%s' has not been defined. Please use the domain user function." % domain) 1705 1706 # Loop over the tensors. 1707 match = False 1708 for tensor_cont in cdp.align_tensors: 1709 # Find the matching tensor and then store the domain label. 1710 if tensor_cont.name == tensor: 1711 tensor_cont.set(param='domain', value=domain) 1712 match = True 1713 1714 # The tensor label doesn't exist. 1715 if not match: 1716 raise RelaxNoTensorError('alignment', tensor)

1717 1718

1719 -def svd(basis_set='irreducible 5D', tensors=None, precision=1):

1720 """Calculate the singular values of all the loaded tensors. 1721 1722 The basis set can be set to one of: 1723 1724 - 'irreducible 5D', the irreducible 5D basis set {A-2, A-1, A0, A1, A2}. This is a linear map, hence angles are preserved. 1725 - '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. 1726 - 'unitary 5D', the unitary 5D basis set {Sxx, Syy, Sxy, Sxz, Syz}. This is a non-linear map, hence angles are not preserved. 1727 - '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. 1728 1729 If the selected basis set is the default of 'irreducible 5D', the matrix on which SVD will be performed will be:: 1730 1731 | A-2(1) A-1(1) A0(1) A1(1) A2(1) | 1732 | A-2(2) A-1(2) A0(2) A1(2) A2(2) | 1733 | A-2(3) A-1(3) A0(3) A1(3) A2(3) | 1734 | . . . . . | 1735 | . . . . . | 1736 | . . . . . | 1737 | A-2(N) A-1(N) A0(N) A1(N) A2(N) | 1738 1739 If the selected basis set is 'unitary 9D', the matrix on which SVD will be performed will be:: 1740 1741 | Sxx1 Sxy1 Sxz1 Syx1 Syy1 Syz1 Szx1 Szy1 Szz1 | 1742 | Sxx2 Sxy2 Sxz2 Syx2 Syy2 Syz2 Szx2 Szy2 Szz2 | 1743 | Sxx3 Sxy3 Sxz3 Syx3 Syy3 Syz3 Szx3 Szy3 Szz3 | 1744 | . . . . . . . . . | 1745 | . . . . . . . . . | 1746 | . . . . . . . . . | 1747 | SxxN SxyN SxzN SyxN SyyN SyzN SzxN SzyN SzzN | 1748 1749 Otherwise if the selected basis set is 'unitary 5D', the matrix for SVD is:: 1750 1751 | Sxx1 Syy1 Sxy1 Sxz1 Syz1 | 1752 | Sxx2 Syy2 Sxy2 Sxz2 Syz2 | 1753 | Sxx3 Syy3 Sxy3 Sxz3 Syz3 | 1754 | . . . . . | 1755 | . . . . . | 1756 | . . . . . | 1757 | SxxN SyyN SxyN SxzN SyzN | 1758 1759 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:: 1760 1761 | Szz1 Sxxyy1 Sxy1 Sxz1 Syz1 | 1762 | Szz2 Sxxyy2 Sxy2 Sxz2 Syz2 | 1763 | Szz3 Sxxyy3 Sxy3 Sxz3 Syz3 | 1764 | . . . . . | 1765 | . . . . . | 1766 | . . . . . | 1767 | SzzN SxxyyN SxyN SxzN SyzN | 1768 1769 For the irreducible basis set, the Am components are defined as:: 1770 1771 / 4pi \ 1/2 1772 A0 = | --- | Szz , 1773 \ 5 / 1774 1775 / 8pi \ 1/2 1776 A+/-1 = +/- | --- | (Sxz +/- iSyz) , 1777 \ 15 / 1778 1779 / 2pi \ 1/2 1780 A+/-2 = | --- | (Sxx - Syy +/- 2iSxy) . 1781 \ 15 / 1782 1783 The relationships between the geometric and unitary basis sets are:: 1784 1785 Szz = - Sxx - Syy, 1786 Sxxyy = Sxx - Syy, 1787 1788 The SVD values and condition number are dependant upon the basis set chosen. 1789 1790 1791 @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". 1792 @type basis_set: str 1793 @keyword tensors: The list of alignment tensor IDs to calculate inter-matrix angles between. If None, all tensors will be used. 1794 @type tensors: None or list of str 1795 @keyword precision: The precision of the printed out angles. The number corresponds to the number of figures to print after the decimal point. 1796 @type precision: int 1797 """ 1798 1799 # Argument check. 1800 allowed = ['irreducible 5D', 'unitary 9D', 'unitary 5D', 'geometric 5D'] 1801 if basis_set not in allowed: 1802 raise RelaxError("The basis set of '%s' is not one of %s." % (basis_set, allowed)) 1803 1804 # Test that alignment tensor data exists. 1805 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0: 1806 raise RelaxNoTensorError('alignment') 1807 1808 # Count the number of tensors used in the SVD. 1809 tensor_num = 0 1810 for tensor in cdp.align_tensors: 1811 if tensors and tensor.name not in tensors: 1812 continue 1813 tensor_num = tensor_num + 1 1814 1815 # Create the matrix to apply SVD on. 1816 if basis_set in ['unitary 9D']: 1817 matrix = zeros((tensor_num, 9), float64) 1818 elif basis_set in ['irreducible 5D']: 1819 matrix = zeros((tensor_num, 5), complex128) 1820 else: 1821 matrix = zeros((tensor_num, 5), float64) 1822 1823 # Pack the elements. 1824 i = 0 1825 for tensor in cdp.align_tensors: 1826 # Skip tensors. 1827 if tensors and tensor.name not in tensors: 1828 continue 1829 1830 # 5D irreducible basis set. 1831 if basis_set == 'irreducible 5D': 1832 matrix[i, 0] = tensor.Am2 1833 matrix[i, 1] = tensor.Am1 1834 matrix[i, 2] = tensor.A0 1835 matrix[i, 3] = tensor.A1 1836 matrix[i, 4] = tensor.A2 1837 1838 # 5D unitary basis set. 1839 elif basis_set == 'unitary 9D': 1840 matrix[i, 0] = tensor.Sxx 1841 matrix[i, 1] = tensor.Sxy 1842 matrix[i, 2] = tensor.Sxz 1843 matrix[i, 3] = tensor.Sxy 1844 matrix[i, 4] = tensor.Syy 1845 matrix[i, 5] = tensor.Syz 1846 matrix[i, 6] = tensor.Sxz 1847 matrix[i, 7] = tensor.Syz 1848 matrix[i, 8] = tensor.Szz 1849 1850 # 5D unitary basis set. 1851 elif basis_set == 'unitary 5D': 1852 matrix[i, 0] = tensor.Sxx 1853 matrix[i, 1] = tensor.Syy 1854 matrix[i, 2] = tensor.Sxy 1855 matrix[i, 3] = tensor.Sxz 1856 matrix[i, 4] = tensor.Syz 1857 1858 # 5D geometric basis set. 1859 elif basis_set == 'geometric 5D': 1860 matrix[i, 0] = tensor.Szz 1861 matrix[i, 1] = tensor.Sxxyy 1862 matrix[i, 2] = tensor.Sxy 1863 matrix[i, 3] = tensor.Sxz 1864 matrix[i, 4] = tensor.Syz 1865 1866 # Increment the index. 1867 i = i + 1 1868 1869 # SVD. 1870 u, s, vh = linalg.svd(matrix) 1871 1872 # Store the singular values. 1873 cdp.align_tensors.singular_vals = s 1874 1875 # Calculate and store the condition number. 1876 cdp.align_tensors.cond_num = s[0] / s[-1] 1877 1878 # Print out. 1879 if basis_set == 'irreducible 5D': 1880 sys.stdout.write("SVD for the irreducible 5D vectors {A-2, A-1, A0, A1, A2}.\n") 1881 elif basis_set == 'unitary 9D': 1882 sys.stdout.write("SVD for the unitary 9D vectors {Sxx, Sxy, Sxz, Syx, Syy, Syz, Szx, Szy, Szz}.\n") 1883 elif basis_set == 'unitary 5D': 1884 sys.stdout.write("SVD for the unitary 5D vectors {Sxx, Syy, Sxy, Sxz, Syz}.\n") 1885 elif basis_set == 'geometric 5D': 1886 sys.stdout.write("SVD for the geometric 5D vectors {Szz, Sxx-yy, Sxy, Sxz, Syz}.\n") 1887 sys.stdout.write("\nSingular values:\n") 1888 for val in s: 1889 format = " %." + repr(precision) + "e\n" 1890 sys.stdout.write(format % val) 1891 format = "\nCondition number: %." + repr(precision) + "f\n" 1892 sys.stdout.write(format % cdp.align_tensors.cond_num)

1893

Source Code for Module pipe_control.align_tensor