data.align_tensor

1 ############################################################################### 2 # # 3 # Copyright (C) 2003-2012 Edward d'Auvergne # 4 # # 5 # This file is part of the program relax (http://www.nmr-relax.com). # 6 # # 7 # This program is free software: you can redistribute it and/or modify # 8 # it under the terms of the GNU General Public License as published by # 9 # the Free Software Foundation, either version 3 of the License, or # 10 # (at your option) any later version. # 11 # # 12 # This program is distributed in the hope that it will be useful, # 13 # but WITHOUT ANY WARRANTY; without even the implied warranty of # 14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the # 15 # GNU General Public License for more details. # 16 # # 17 # You should have received a copy of the GNU General Public License # 18 # along with this program. If not, see <http://www.gnu.org/licenses/>. # 19 # # 20 ############################################################################### 21 22 # Python module imports. 23 from re import search 24 from math import cos, sin 25 from numpy import array, dot, eye, float64, identity, transpose, zeros 26 from numpy.linalg import det, eig, eigvals 27 import sys 28 29 # relax module imports. 30 from compat import py_version 31 from data.data_classes import Element 32 from data.relax_xml import fill_object_contents, xml_to_object 33 from float import nan 34 from maths_fns.rotation_matrix import R_to_euler_zyz 35 from relax_errors import RelaxError 36 37

38 -def calc_A(Axx, Ayy, Azz, Axy, Axz, Ayz):

39 """Function for calculating the alignment tensor (in the structural frame). 40 41 @param Axx: The Axx tensor element. 42 @type Axx: float 43 @param Ayy: The Ayy tensor element. 44 @type Ayy: float 45 @param Azz: The Azz tensor element. 46 @type Azz: float 47 @param Axy: The Axy tensor element. 48 @type Axy: float 49 @param Axz: The Axz tensor element. 50 @type Axz: float 51 @param Ayz: The Ayz tensor element. 52 @type Ayz: float 53 @return: The alignment tensor (within the structural frame). 54 @rtype: 3x3 numpy float64 array 55 """ 56 57 # Initialise the tensor. 58 tensor = zeros((3, 3), float64) 59 60 # Populate the diagonal elements. 61 tensor[0, 0] = Axx 62 tensor[1, 1] = Ayy 63 tensor[2, 2] = Azz 64 65 # Populate the off diagonal elements. 66 tensor[0, 1] = tensor[1, 0] = Axy 67 tensor[0, 2] = tensor[2, 0] = Axz 68 tensor[1, 2] = tensor[2, 1] = Ayz 69 70 # Return the tensor. 71 return tensor

72 73

74 -def calc_A_5D(Axx, Ayy, Azz, Axy, Axz, Ayz):

75 """Function for calculating the alignment tensor in the 5D vector notation. 76 77 @param Axx: The Axx tensor element. 78 @type Axx: float 79 @param Ayy: The Ayy tensor element. 80 @type Ayy: float 81 @param Azz: The Azz tensor element. 82 @type Azz: float 83 @param Axy: The Axy tensor element. 84 @type Axy: float 85 @param Axz: The Axz tensor element. 86 @type Axz: float 87 @param Ayz: The Ayz tensor element. 88 @type Ayz: float 89 @return: The alignment 5D tensor (within the structural frame). 90 @rtype: numpy rank-1 5D tensor 91 """ 92 93 # Initialise the tensor. 94 tensor = zeros(5, float64) 95 96 # Populate the tensor. 97 tensor[0] = Axx 98 tensor[1] = Ayy 99 tensor[2] = Axy 100 tensor[3] = Axz 101 tensor[4] = Ayz 102 103 # Return the tensor. 104 return tensor

105 106

107 -def calc_A_diag(A):

108 """Calculate the diagonalised alignment tensor. 109 110 The diagonalised alignment tensor is defined as:: 111 112 | Axx' 0 0 | 113 tensor = | 0 Ayy' 0 |. 114 | 0 0 Azz'| 115 116 The diagonalised alignment tensor is calculated by eigenvalue decomposition. 117 118 119 @param A: The full alignment tensor. 120 @type A: numpy array ((3, 3), float64) 121 @return: The diagonalised alignment tensor. 122 @rtype: numpy array ((3, 3), float64) 123 """ 124 125 # The eigenvalues. 126 vals = eigvals(A) 127 128 # Find the |x| < |y| < |z| indices. 129 abs_vals = abs(vals).tolist() 130 Axx_index = abs_vals.index(min(abs_vals)) 131 Azz_index = abs_vals.index(max(abs_vals)) 132 last_index = list(range(3)) 133 last_index.pop(max(Axx_index, Azz_index)) 134 last_index.pop(min(Axx_index, Azz_index)) 135 Ayy_index = last_index[0] 136 137 # Empty tensor. 138 tensor_diag = zeros((3, 3), float64) 139 140 # Fill the elements. 141 tensor_diag[0, 0] = vals[Axx_index] 142 tensor_diag[1, 1] = vals[Ayy_index] 143 tensor_diag[2, 2] = vals[Azz_index] 144 145 # Return the tensor. 146 return tensor_diag

147 148

149 -def calc_Aa(A_diag):

150 """Calculate the anisotropic parameter Aa. 151 152 This is given by:: 153 154 Aa = 3/2Azz = Szz, 155 156 where Azz and Szz are the eigenvalues. 157 158 159 @param A_diag: The full alignment tensor, diagonalised. 160 @type A_diag: numpy array ((3, 3), float64) 161 @return: The Aa parameter 162 @rtype: float 163 """ 164 165 # Return Aa. 166 return 1.5 * A_diag[2, 2]

167 168

169 -def calc_Ar(A_diag):

170 """Calculate the rhombic parameter Ar. 171 172 This is given by:: 173 174 Ar = Axx - Ayy, 175 176 where Axx and Ayy are the eigenvalues. 177 178 179 @param A_diag: The full alignment tensor, diagonalised. 180 @type A_diag: numpy array ((3, 3), float64) 181 @return: The Ar parameter 182 @rtype: float 183 """ 184 185 # Return Ar. 186 return A_diag[0, 0] - A_diag[1, 1]

187 188

189 -def calc_Axxyy(Axx, Ayy):

190 """Function for calculating the Axx-yy value. 191 192 The equation for calculating the parameter is:: 193 194 Axx-yy = Axx - Ayy. 195 196 @param Axx: The Axx component of the alignment tensor. 197 @type Axx: float 198 @param Ayy: The Ayy component of the alignment tensor. 199 @type Ayy: float 200 @return: The Axx-yy component of the alignment tensor. 201 @rtype: float 202 """ 203 204 # Calculate and return the Axx-yy value. 205 return Axx - Ayy

206 207

208 -def calc_Azz(Axx, Ayy):

209 """Function for calculating the Azz value. 210 211 The equation for calculating the parameter is:: 212 213 Azz = - Axx - Ayy. 214 215 @param Axx: The Axx component of the alignment tensor. 216 @type Axx: float 217 @param Ayy: The Ayy component of the alignment tensor. 218 @type Ayy: float 219 @return: The Azz component of the alignment tensor. 220 @rtype: float 221 """ 222 223 # Calculate and return the Azz value. 224 return - Axx - Ayy

225 226

227 -def calc_eigvals(A):

228 """Calculate the eigenvalues and eigenvectors of the alignment tensor (A). 229 230 @param A: The full alignment tensor. 231 @type A: numpy array ((3, 3), float64) 232 @return: The eigensystem. 233 @rtype: tuple of numpy array (float64) 234 """ 235 236 # The eigenvalues. 237 vals = eigvals(A) 238 239 # Find the |x| < |y| < |z| indices. 240 abs_vals = abs(vals).tolist() 241 x_index = abs_vals.index(min(abs_vals)) 242 z_index = abs_vals.index(max(abs_vals)) 243 last_index = list(range(3)) 244 last_index.pop(max(x_index, z_index)) 245 last_index.pop(min(x_index, z_index)) 246 y_index = last_index[0] 247 248 # Return the sorted eigenvalues. 249 return [vals[x_index], vals[y_index], vals[z_index]]

250 251

252 -def calc_eta(A_diag):

253 """Calculate the asymmetry parameter eta. 254 255 This is given by:: 256 257 eta = (Axx - Ayy) / Azz 258 259 where Aii are the eigenvalues. 260 261 262 @param A_diag: The full alignment tensor, diagonalised. 263 @type A_diag: numpy array ((3, 3), float64) 264 @return: The eta parameter 265 @rtype: float 266 """ 267 268 # Zero Azz value, so return NaN. 269 if A_diag[2, 2] == 0: 270 return nan 271 272 # Return eta. 273 return (A_diag[0, 0] - A_diag[1, 1]) / A_diag[2, 2]

274 275

276 -def calc_euler(rotation):

277 """Calculate the zyz notation Euler angles. 278 279 @param rotation: The rotation matrix. 280 @type rotation: numpy 3D, rank-2 array 281 @return: The Euler angles alpha, beta, and gamma in zyz notation. 282 @rtype: tuple of float 283 """ 284 285 return R_to_euler_zyz(rotation)

286 287

288 -def calc_S(Sxx, Syy, Szz, Sxy, Sxz, Syz):

289 """Function for calculating the alignment tensor (in the structural frame). 290 291 @param Sxx: The Sxx tensor element. 292 @type Sxx: float 293 @param Syy: The Syy tensor element. 294 @type Syy: float 295 @param Szz: The Szz tensor element. 296 @type Szz: float 297 @param Sxy: The Sxy tensor element. 298 @type Sxy: float 299 @param Sxz: The Sxz tensor element. 300 @type Sxz: float 301 @param Syz: The Syz tensor element. 302 @type Syz: float 303 @return: The alignment tensor (within the structural frame). 304 @rtype: 3x3 numpy float64 array 305 """ 306 307 # Initialise the tensor. 308 tensor = zeros((3, 3), float64) 309 310 # Populate the diagonal elements. 311 tensor[0, 0] = Sxx 312 tensor[1, 1] = Syy 313 tensor[2, 2] = Szz 314 315 # Populate the off diagonal elements. 316 tensor[0, 1] = tensor[1, 0] = Sxy 317 tensor[0, 2] = tensor[2, 0] = Sxz 318 tensor[1, 2] = tensor[2, 1] = Syz 319 320 # Return the tensor. 321 return tensor

322 323

324 -def calc_S_5D(Sxx, Syy, Szz, Sxy, Sxz, Syz):

325 """Function for calculating the alignment tensor in the 5D vector notation. 326 327 @param Sxx: The Sxx tensor element. 328 @type Sxx: float 329 @param Syy: The Syy tensor element. 330 @type Syy: float 331 @param Szz: The Szz tensor element. 332 @type Szz: float 333 @param Sxy: The Sxy tensor element. 334 @type Sxy: float 335 @param Sxz: The Sxz tensor element. 336 @type Sxz: float 337 @param Syz: The Syz tensor element. 338 @type Syz: float 339 @return: The alignment 5D tensor (within the structural frame). 340 @rtype: numpy rank-1 5D tensor 341 """ 342 343 # Initialise the tensor. 344 tensor = zeros(5, float64) 345 346 # Populate the tensor. 347 tensor[0] = Sxx 348 tensor[1] = Syy 349 tensor[2] = Sxy 350 tensor[3] = Sxz 351 tensor[4] = Syz 352 353 # Return the tensor. 354 return tensor

355 356

357 -def calc_S_diag(tensor):

358 """Calculate the diagonalised alignment tensor. 359 360 The diagonalised alignment tensor is defined as:: 361 362 | Sxx' 0 0 | 363 tensor = | 0 Syy' 0 |. 364 | 0 0 Szz'| 365 366 The diagonalised alignment tensor is calculated by eigenvalue decomposition. 367 368 369 @param tensor: The full alignment tensor in its eigenframe. 370 @type tensor: numpy array ((3, 3), float64) 371 @return: The diagonalised alignment tensor. 372 @rtype: numpy array ((3, 3), float64) 373 """ 374 375 # The eigenvalues. 376 vals = eigvals(tensor) 377 378 # Find the |x| < |y| < |z| indices. 379 abs_vals = abs(vals).tolist() 380 Sxx_index = abs_vals.index(min(abs_vals)) 381 Szz_index = abs_vals.index(max(abs_vals)) 382 last_index = list(range(3)) 383 last_index.pop(max(Sxx_index, Szz_index)) 384 last_index.pop(min(Sxx_index, Szz_index)) 385 Syy_index = last_index[0] 386 387 # Empty tensor. 388 tensor_diag = zeros((3, 3), float64) 389 390 # Fill the elements. 391 tensor_diag[0, 0] = vals[Sxx_index] 392 tensor_diag[1, 1] = vals[Syy_index] 393 tensor_diag[2, 2] = vals[Szz_index] 394 395 # Return the tensor. 396 return tensor_diag

397 398

399 -def calc_Sxx(Axx):

400 """Function for calculating the Axx value. 401 402 The equation for calculating the parameter is:: 403 404 Sxx = 3/2 Axx. 405 406 @param Axx: The Axx component of the alignment tensor. 407 @type Axx: float 408 @rtype: float 409 """ 410 411 # Calculate and return the Axx value. 412 return 3.0/2.0 * Axx

413 414

415 -def calc_Sxxyy(Sxx, Syy):

416 """Function for calculating the Sxx-yy value. 417 418 The equation for calculating the parameter is:: 419 420 Sxx-yy = Sxx - Syy. 421 422 @param Sxx: The Sxx component of the Saupe order matrix. 423 @type Sxx: float 424 @param Syy: The Syy component of the Saupe order matrix. 425 @type Syy: float 426 @return: The Sxx-yy component of the Saupe order matrix. 427 @rtype: float 428 """ 429 430 # Calculate and return the Sxx-yy value. 431 return Sxx - Syy

432 433

434 -def calc_Sxy(Axy):

435 """Function for calculating the Axy value. 436 437 The equation for calculating the parameter is:: 438 439 Sxy = 3/2 Axy. 440 441 @param Axy: The Axy component of the alignment tensor. 442 @type Axy: float 443 @rtype: float 444 """ 445 446 # Calculate and return the Axy value. 447 return 3.0/2.0 * Axy

448 449

450 -def calc_Sxz(Axz):

451 """Function for calculating the Axz value. 452 453 The equation for calculating the parameter is:: 454 455 Sxz = 3/2 Axz. 456 457 @param Axz: The Axz component of the alignment tensor. 458 @type Axz: float 459 @rtype: float 460 """ 461 462 # Calculate and return the Axz value. 463 return 3.0/2.0 * Axz

464 465

466 -def calc_Syy(Ayy):

467 """Function for calculating the Ayy value. 468 469 The equation for calculating the parameter is:: 470 471 Syy = 3/2 Ayy. 472 473 @param Ayy: The Ayy component of the alignment tensor. 474 @type Ayy: float 475 @rtype: float 476 """ 477 478 # Calculate and return the Ayy value. 479 return 3.0/2.0 * Ayy

480 481

482 -def calc_Syz(Ayz):

483 """Function for calculating the Ayz value. 484 485 The equation for calculating the parameter is:: 486 487 Syz = 3/2 Ayz. 488 489 @param Ayz: The Ayz component of the alignment tensor. 490 @type Ayz: float 491 @rtype: float 492 """ 493 494 # Calculate and return the Ayz value. 495 return 3.0/2.0 * Ayz

496 497

498 -def calc_Szz(Sxx, Syy):

499 """Function for calculating the Szz value. 500 501 The equation for calculating the parameter is:: 502 503 Szz = - Sxx - Syy. 504 505 @param Sxx: The Sxx component of the Saupe order matrix. 506 @type Sxx: float 507 @param Syy: The Syy component of the Saupe order matrix. 508 @type Syy: float 509 @return: The Szz component of the Saupe order matrix. 510 @rtype: float 511 """ 512 513 # Calculate and return the Szz value. 514 return - Sxx - Syy

515 516

517 -def calc_P(Axx, Ayy, Azz, Axy, Axz, Ayz):

518 """Function for calculating the alignment tensor (in the structural frame). 519 520 @param Axx: The Axx tensor element. 521 @type Axx: float 522 @param Ayy: The Ayy tensor element. 523 @type Ayy: float 524 @param Azz: The Azz tensor element. 525 @type Azz: float 526 @param Axy: The Axy tensor element. 527 @type Axy: float 528 @param Axz: The Axz tensor element. 529 @type Axz: float 530 @param Ayz: The Ayz tensor element. 531 @type Ayz: float 532 @return: The alignment tensor (within the structural frame). 533 @rtype: 3x3 numpy float64 array 534 """ 535 536 # Initialise the tensor. 537 tensor = zeros((3, 3), float64) 538 539 # Populate the diagonal elements. 540 tensor[0, 0] = Axx 541 tensor[1, 1] = Ayy 542 tensor[2, 2] = Azz 543 544 # Populate the off diagonal elements. 545 tensor[0, 1] = tensor[1, 0] = Axy 546 tensor[0, 2] = tensor[2, 0] = Axz 547 tensor[1, 2] = tensor[2, 1] = Ayz 548 549 # Add 1/3 the identity matrix. 550 tensor = tensor + eye(3)/3.0 551 552 # Return the tensor. 553 return tensor

554 555

556 -def calc_P_5D(Axx, Ayy, Azz, Axy, Axz, Ayz):

557 """Function for calculating the alignment tensor in the 5D vector notation. 558 559 @param Axx: The Axx tensor element. 560 @type Axx: float 561 @param Ayy: The Ayy tensor element. 562 @type Ayy: float 563 @param Azz: The Azz tensor element. 564 @type Azz: float 565 @param Axy: The Axy tensor element. 566 @type Axy: float 567 @param Axz: The Axz tensor element. 568 @type Axz: float 569 @param Ayz: The Ayz tensor element. 570 @type Ayz: float 571 @return: The alignment 5D tensor (within the structural frame). 572 @rtype: numpy rank-1 5D tensor 573 """ 574 575 # Initialise the tensor. 576 tensor = zeros(5, float64) 577 578 # Populate the tensor. 579 tensor[0] = Axx + 1.0/3.0 580 tensor[1] = Ayy + 1.0/3.0 581 tensor[2] = Axy 582 tensor[3] = Axz 583 tensor[4] = Ayz 584 585 # Return the tensor. 586 return tensor

587 588

589 -def calc_P_diag(tensor):

590 """Calculate the diagonalised alignment tensor. 591 592 The diagonalised alignment tensor is defined as:: 593 594 | Pxx' 0 0 | 595 tensor = | 0 Pyy' 0 |. 596 | 0 0 Pzz'| 597 598 The diagonalised alignment tensor is calculated by eigenvalue decomposition. 599 600 601 @param tensor: The full alignment tensor in its eigenframe. 602 @type tensor: numpy array ((3, 3), float64) 603 @return: The diagonalised alignment tensor. 604 @rtype: numpy array ((3, 3), float64) 605 """ 606 607 # The eigenvalues. 608 vals = eigvals(tensor) 609 610 # Find the |x| < |y| < |z| indices. 611 abs_vals = abs(vals).tolist() 612 Pxx_index = abs_vals.index(min(abs_vals)) 613 Pzz_index = abs_vals.index(max(abs_vals)) 614 last_index = list(range(3)) 615 last_index.pop(max(Pxx_index, Pzz_index)) 616 last_index.pop(min(Pxx_index, Pzz_index)) 617 Pyy_index = last_index[0] 618 619 # Empty tensor. 620 tensor_diag = zeros((3, 3), float64) 621 622 # Fill the elements. 623 tensor_diag[0, 0] = vals[Pxx_index] 624 tensor_diag[1, 1] = vals[Pyy_index] 625 tensor_diag[2, 2] = vals[Pzz_index] 626 627 # Add 1/3 the identity matrix. 628 tensor = tensor + eye(3)/3.0 629 630 # Return the tensor. 631 return tensor_diag

632 633

634 -def calc_Pxx(Axx):

635 """Function for calculating the Pxx value. 636 637 The equation for calculating the parameter is:: 638 639 Pxx = Axx + 1/3. 640 641 @param Axx: The Axx component of the alignment tensor. 642 @type Axx: float 643 @rtype: float 644 """ 645 646 # Calculate and return the Pxx value. 647 return Axx + 1.0/3.0

648 649

650 -def calc_Pxxyy(Pxx, Pyy):

651 """Function for calculating the Pxx-yy value. 652 653 The equation for calculating the parameter is:: 654 655 Pxx-yy = Pxx - Pyy. 656 657 @param Pxx: The Pxx component of the alignment tensor. 658 @type Pxx: float 659 @param Pyy: The Pyy component of the alignment tensor. 660 @type Pyy: float 661 @return: The Pxx-yy component of the alignment tensor. 662 @rtype: float 663 """ 664 665 # Calculate and return the Pxx-yy value. 666 return Pxx - Pyy

667 668

669 -def calc_Pxy(Axy):

670 """Function for calculating the Pxy value. 671 672 The equation for calculating the parameter is:: 673 674 Pxy = Axy. 675 676 @param Axy: The Axy component of the alignment tensor. 677 @type Axy: float 678 @rtype: float 679 """ 680 681 # Calculate and return the Pxy value. 682 return Axy

683 684

685 -def calc_Pxz(Axz):

686 """Function for calculating the Pxz value. 687 688 The equation for calculating the parameter is:: 689 690 Pxz = Axz. 691 692 @param Axz: The Axz component of the alignment tensor. 693 @type Axz: float 694 @rtype: float 695 """ 696 697 # Calculate and return the Pxz value. 698 return Axz

699 700

701 -def calc_Pyy(Ayy):

702 """Function for calculating the Pyy value. 703 704 The equation for calculating the parameter is:: 705 706 Pyy = Ayy + 1/3. 707 708 @param Ayy: The Ayy component of the alignment tensor. 709 @type Ayy: float 710 @rtype: float 711 """ 712 713 # Calculate and return the Pyy value. 714 return Ayy + 1.0/3.0

715 716

717 -def calc_Pyz(Ayz):

718 """Function for calculating the Pyz value. 719 720 The equation for calculating the parameter is:: 721 722 Pyz = Ayz. 723 724 @param Ayz: The Ayz component of the alignment tensor. 725 @type Ayz: float 726 @rtype: float 727 """ 728 729 # Calculate and return the Pyz value. 730 return Ayz

731 732

733 -def calc_Pzz(Pxx, Pyy):

734 """Function for calculating the Pzz value. 735 736 The equation for calculating the parameter is:: 737 738 Pzz = 1 - Pxx - Pyy. 739 740 @param Pxx: The Pxx component of the alignment tensor. 741 @type Pxx: float 742 @param Pyy: The Pyy component of the alignment tensor. 743 @type Pyy: float 744 @return: The Pzz component of the alignment tensor. 745 @rtype: float 746 """ 747 748 # Calculate and return the Pzz value. 749 return 1.0 - Pxx - Pyy

750 751

752 -def calc_R(Aa, Ar):

753 """Calculate the rhombicity parameter R. 754 755 This is given by:: 756 757 R = Ar / Aa. 758 759 760 @param Aa: The Aa parameter. 761 @type Aa: float 762 @param Ar: The Ar parameter. 763 @type Ar: float 764 @return: The R parameter. 765 @rtype: float 766 """ 767 768 # Zero Aa value, so return NaN. 769 if Aa == 0: 770 return nan 771 772 # Return R. 773 return Ar / Aa

774 775

776 -def calc_rotation(A):

777 """Calculate the rotation matrix from the molecular frame to the tensor frame. 778 779 This is defined by:: 780 781 | Azz | >= | Ayy | >= | Axx |. 782 783 784 @param A: The full alignment tensor. 785 @type A: numpy array ((3, 3), float64) 786 @return: The array of x, y, and z indices. 787 @rtype: list 788 """ 789 790 # The eigenvalues. 791 vals, rot = eig(A) 792 793 # Find the |x| < |y| < |z| indices. 794 abs_vals = abs(vals).tolist() 795 x_index = abs_vals.index(min(abs_vals)) 796 z_index = abs_vals.index(max(abs_vals)) 797 last_index = list(range(3)) 798 last_index.pop(max(x_index, z_index)) 799 last_index.pop(min(x_index, z_index)) 800 y_index = last_index[0] 801 802 # Empty rotation matrix for index permutations. 803 rot_perm = zeros((3, 3), float64) 804 805 # Permute the rotation matrix. 806 perm = [x_index, y_index, z_index] 807 for i in range(3): 808 for j in range(3): 809 rot_perm[i, j] = rot[i, perm[j]] 810 811 # Switch from the left handed to right handed universe if required. 812 if abs(det(rot_perm) - 1.0) > 1e-7: 813 rot_perm[:, 0] = -rot_perm[:, 0] 814 815 # Return the permuted rotation matrix. 816 return rot_perm

817 818

819 -def calc_unit_x(rotation):

820 """Calculate the x unit vector. 821 822 This is given by the eigenvalue decomposition. 823 824 825 @param rotation: The rotation matrix. 826 @type rotation: numpy 3D, rank-2 array 827 @return: The x unit vector. 828 @rtype: numpy array (float64) 829 """ 830 831 # Return the x unit vector. 832 return rotation[:, 0]

833 834

835 -def calc_unit_y(rotation):

836 """Calculate the y unit vector. 837 838 This is given by the eigenvalue decomposition. 839 840 841 @param rotation: The rotation matrix. 842 @type rotation: numpy 3D, rank-2 array 843 @return: The y unit vector. 844 @rtype: numpy array (float64) 845 """ 846 847 # Return the y unit vector. 848 return rotation[:, 1]

849 850

851 -def calc_unit_z(rotation):

852 """Calculate the z unit vector. 853 854 This is given by the eigenvalue decomposition. 855 856 857 @param rotation: The rotation matrix. 858 @type rotation: numpy 3D, rank-2 array 859 @return: The z unit vector. 860 @rtype: numpy array (float64) 861 """ 862 863 # Return the z unit vector. 864 return rotation[:, 2]

865 866

867 -def dependency_generator():

868 """Generator for the automatic updating the alignment tensor data structures. 869 870 @return: This generator successively yields three objects, the target object to update, the list of parameters which if modified cause the target to be updated, and the list of parameters that the target depends upon. 871 """ 872 873 # Primary objects (only dependant on the modifiable objects). 874 yield ('A', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 875 yield ('A_5D', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 876 yield ('Axxyy', ['Axx', 'Ayy'], ['Axx', 'Ayy']) 877 yield ('Azz', ['Axx', 'Ayy'], ['Axx', 'Ayy']) 878 879 yield ('P', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 880 yield ('P_5D', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 881 yield ('Pxx', ['Axx'], ['Axx']) 882 yield ('Pxxyy', ['Axx', 'Ayy'], ['Axx', 'Ayy']) 883 yield ('Pxy', ['Axy'], ['Axy']) 884 yield ('Pxz', ['Axz'], ['Axz']) 885 yield ('Pyy', ['Ayy'], ['Ayy']) 886 yield ('Pyz', ['Ayz'], ['Ayz']) 887 888 yield ('Sxx', ['Axx'], ['Axx']) 889 yield ('Sxy', ['Axy'], ['Axy']) 890 yield ('Sxz', ['Axz'], ['Axz']) 891 yield ('Syy', ['Ayy'], ['Ayy']) 892 yield ('Syz', ['Ayz'], ['Ayz']) 893 894 # Secondary objects (dependant on the primary objects). 895 yield ('A_diag', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A']) 896 yield ('eigvals', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A']) 897 yield ('rotation', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A']) 898 899 yield ('P_diag', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['P']) 900 yield ('Pzz', ['Axx', 'Ayy'], ['Pxx', 'Pyy']) 901 902 yield ('S', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Sxx', 'Syy', 'Szz', 'Sxy', 'Sxz', 'Syz']) 903 yield ('S_5D', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Sxx', 'Syy', 'Szz', 'Sxy', 'Sxz', 'Syz']) 904 yield ('S_diag', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['S']) 905 yield ('Sxxyy', ['Axx', 'Ayy'], ['Sxx', 'Syy']) 906 yield ('Szz', ['Axx', 'Ayy'], ['Sxx', 'Syy']) 907 908 # Tertiary objects (dependant on the secondary objects). 909 yield ('Aa', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A_diag']) 910 yield ('Ar', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A_diag']) 911 yield ('eta', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A_diag']) 912 913 yield ('unit_x', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 914 yield ('unit_y', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 915 yield ('unit_z', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 916 917 yield ('euler', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 918 919 # Quaternary objects (dependant on the tertiary objects). 920 yield ('R', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Aa', 'Ar'])

921 922 923 924 # Alignment tensor specific data. 925 ################################# 926

927 -class AlignTensorList(list):

928 """List type data container for holding all the alignment tensors. 929 930 The elements of the list should be AlignTensorData instances. 931 """ 932

933 - def __repr__(self):

934 """Replacement function for displaying an instance of this class.""" 935 936 text = "Alignment tensors.\n\n" 937 text = text + "%-8s%-20s\n" % ("Index", "Name") 938 for i in range(len(self)): 939 text = text + "%-8i%-20s\n" % (i, self[i].name) 940 text = text + "\nThese can be accessed by typing 'pipe.align_tensor[index]'.\n" 941 return text

942 943

944 - def add_item(self, name):

945 """Function for appending a new AlignTensorData instance to the list.""" 946 947 self.append(AlignTensorData(name))

948 949

950 - def from_xml(self, align_tensor_super_node, file_version=1):

951 """Recreate the alignment tensor data structure from the XML alignment tensor node. 952 953 @param align_tensor_super_node: The alignment tensor XML nodes. 954 @type align_tensor_super_node: xml.dom.minicompat.Element instance 955 @keyword file_version: The relax XML version of the XML file. 956 @type file_version: int 957 """ 958 959 # Recreate all the alignment tensor data structures. 960 xml_to_object(align_tensor_super_node, self, file_version=file_version, blacklist=['align_tensor']) 961 962 # Get the individual tensors. 963 align_tensor_nodes = align_tensor_super_node.getElementsByTagName('align_tensor') 964 965 # Loop over the child nodes. 966 for align_tensor_node in align_tensor_nodes: 967 # Add the alignment tensor data container. 968 self.add_item(align_tensor_node.getAttribute('name')) 969 970 # A temporary object to pack the structures from the XML data into. 971 temp_obj = Element() 972 973 # Recreate all the other data structures (into the temporary object). 974 xml_to_object(align_tensor_node, temp_obj, file_version=file_version) 975 976 # Loop over all modifiable objects in the temporary object and make soft copies of them. 977 for name in self[-1]._mod_attr: 978 # Skip if missing from the object. 979 if not hasattr(temp_obj, name): 980 continue 981 982 # The category. 983 if search('_err$', name): 984 category = 'err' 985 param = name.replace('_err', '') 986 elif search('_sim$', name): 987 category = 'sim' 988 param = name.replace('_sim', '') 989 else: 990 category = 'val' 991 param = name 992 993 # Get the object. 994 value = getattr(temp_obj, name) 995 996 # Normal parameters. 997 if category == 'val': 998 self[-1].set(param=param, value=value) 999 1000 # Errors. 1001 elif category == 'err': 1002 self[-1].set(param=param, value=value, category='err') 1003 1004 # Simulation objects objects. 1005 else: 1006 # Recreate the list elements. 1007 for i in range(len(value)): 1008 self[-1].set(param=param, value=value[i], category='sim', sim_index=i) 1009 1010 # Delete the temporary object. 1011 del temp_obj

1012 1013

1014 - def names(self):

1015 """Return a list of the alignment tensor names.""" 1016 1017 # Loop over the tensors. 1018 names = [] 1019 for i in range(len(self)): 1020 names.append(self[i].name) 1021 1022 # Return the list. 1023 return names

1024 1025

1026 - def to_xml(self, doc, element):

1027 """Create an XML element for the alignment tensors. 1028 1029 @param doc: The XML document object. 1030 @type doc: xml.dom.minidom.Document instance 1031 @param element: The element to add the alignment tensors XML element to. 1032 @type element: XML element object 1033 """ 1034 1035 # Create the alignment tensors element and add it to the higher level element. 1036 tensor_list_element = doc.createElement('align_tensors') 1037 element.appendChild(tensor_list_element) 1038 1039 # Set the alignment tensor attributes. 1040 tensor_list_element.setAttribute('desc', 'Alignment tensor list') 1041 1042 # Add all simple python objects within the PipeContainer to the pipe element. 1043 fill_object_contents(doc, tensor_list_element, object=self, blacklist=list(list(self.__class__.__dict__.keys()) + list(list.__dict__.keys()))) 1044 1045 # Loop over the tensors. 1046 for i in range(len(self)): 1047 # Create an XML element for a single tensor. 1048 tensor_element = doc.createElement('align_tensor') 1049 tensor_list_element.appendChild(tensor_element) 1050 tensor_element.setAttribute('index', repr(i)) 1051 tensor_element.setAttribute('desc', 'Alignment tensor') 1052 1053 # The blacklist. 1054 blacklist = ['type', 'is_empty'] + list(self[i].__class__.__dict__.keys()) 1055 for name in dir(self): 1056 if name not in self[i]._mod_attr: 1057 blacklist.append(name) 1058 1059 # Add all simple python objects within the PipeContainer to the pipe element. 1060 fill_object_contents(doc, tensor_element, object=self[i], blacklist=blacklist)

1061 1062

1063 -class AlignTensorData(Element):

1064 """An empty data container for the alignment tensor elements.""" 1065 1066 # List of modifiable attributes. 1067 _mod_attr = [ 1068 'name', 1069 'Axx', 'Axx_sim', 'Axx_err', 1070 'Ayy', 'Ayy_sim', 'Ayy_err', 1071 'Axy', 'Axy_sim', 'Axy_err', 1072 'Axz', 'Axz_sim', 'Axz_err', 1073 'Ayz', 'Ayz_sim', 'Ayz_err', 1074 'domain', 1075 'red', 1076 'fixed' 1077 ] 1078

1079 - def __init__(self, name, fixed=False):

1080 """Set up the tensor data. 1081 1082 @param name: The tensor ID string. 1083 @type name: str 1084 @keyword fixed: The optimisation flag. 1085 @type fixed: bool 1086 """ 1087 1088 # Store the values. 1089 self.__dict__['name'] = name 1090 self.__dict__['fixed'] = fixed 1091 1092 # The number of simulations. 1093 self.__dict__['_sim_num'] = None

1094 1095

1096 - def __setattr__(self, name, value):

1097 """Make this object read-only.""" 1098 1099 raise RelaxError("The alignment tensor is a read-only object. The alignment tensor set() method must be used instead.")

1100 1101

1102 - def _update_object(self, param_name, target, update_if_set, depends, category):

1103 """Function for updating the target object, its error, and the MC simulations. 1104 1105 If the base name of the object is not within the 'update_if_set' list, this function returns 1106 without doing anything (to avoid wasting time). Dependant upon the category the object 1107 (target), its error (target+'_err'), or all Monte Carlo simulations (target+'_sim') are 1108 updated. 1109 1110 @param param_name: The parameter name which is being set in the __setattr__() function. 1111 @type param_name: str 1112 @param target: The name of the object to update. 1113 @type target: str 1114 @param update_if_set: If the parameter being set by the __setattr__() function is not 1115 within this list of parameters, don't waste time updating the 1116 target. 1117 @param depends: An array of names objects that the target is dependent upon. 1118 @type depends: array of str 1119 @param category: The category of the object to update (one of 'val', 'err', or 1120 'sim'). 1121 @type category: str 1122 @return: None 1123 """ 1124 1125 # Only update if the parameter name is within the 'update_if_set' list. 1126 if not param_name in update_if_set: 1127 return 1128 1129 # Get the function for calculating the value. 1130 fn = globals()['calc_'+target] 1131 1132 1133 # The value. 1134 ############ 1135 1136 if category == 'val': 1137 # Get all the dependencies if possible. 1138 missing_dep = 0 1139 deps = () 1140 for dep_name in depends: 1141 # Test if the object exists. 1142 if not hasattr(self, dep_name): 1143 missing_dep = 1 1144 break 1145 1146 # Get the object and place it into the 'deps' tuple. 1147 deps = deps+(getattr(self, dep_name),) 1148 1149 # Only update the object if its dependencies exist. 1150 if not missing_dep: 1151 # Calculate the value. 1152 value = fn(*deps) 1153 1154 # Set the attribute. 1155 self.__dict__[target] = value 1156 1157 1158 # The error. 1159 ############ 1160 1161 if category == 'err': 1162 # Get all the dependencies if possible. 1163 missing_dep = 0 1164 deps = () 1165 for dep_name in depends: 1166 # Test if the error object exists. 1167 if not hasattr(self, dep_name+'_err'): 1168 missing_dep = 1 1169 break 1170 1171 # Get the object and place it into the 'deps' tuple. 1172 deps = deps+(getattr(self, dep_name+'_err'),) 1173 1174 # Only update the error object if its dependencies exist. 1175 if not missing_dep: 1176 # Calculate the value. 1177 value = fn(*deps) 1178 1179 # Set the attribute. 1180 self.__dict__[target+'_err'] = value 1181 1182 1183 # The Monte Carlo simulations. 1184 ############################## 1185 1186 if category == 'sim': 1187 # Get all the dependencies if possible. 1188 missing_dep = 0 1189 deps = [] 1190 for dep_name in depends: 1191 # Modify the dependency name. 1192 if dep_name != 'type': 1193 dep_name = dep_name+'_sim' 1194 1195 # Test if the MC sim object exists. 1196 if not hasattr(self, dep_name): 1197 missing_dep = 1 1198 break 1199 1200 # Get the object and place it into the 'deps' tuple. 1201 deps.append(getattr(self, dep_name)) 1202 1203 # Only create the MC simulation object if its dependencies exist. 1204 if not missing_dep: 1205 # Initialise an empty array to store the MC simulation object elements (if it doesn't already exist). 1206 if not target+'_sim' in self.__dict__: 1207 self.__dict__[target+'_sim'] = AlignTensorSimList(elements=self._sim_num) 1208 1209 # Repackage the deps structure. 1210 args = [] 1211 skip = False 1212 for i in range(self._sim_num): 1213 args.append(()) 1214 1215 # Loop over the dependent structures. 1216 for j in range(len(deps)): 1217 # None, so skip. 1218 if deps[j] == None or deps[j][i] == None: 1219 skip = True 1220 1221 # String data type. 1222 if isinstance(deps[j], str): 1223 args[-1] = args[-1] + (deps[j],) 1224 1225 # List data type. 1226 else: 1227 args[-1] = args[-1] + (deps[j][i],) 1228 1229 # Loop over the sims and set the values. 1230 if not skip: 1231 for i in range(self._sim_num): 1232 # Calculate the value. 1233 value = fn(*args[i]) 1234 1235 # Set the attribute. 1236 self.__dict__[target+'_sim']._set(value=value, sim_index=i)

1237 1238

1239 - def set(self, param=None, value=None, category='val', sim_index=None):

1240 """Set a alignment tensor parameter. 1241 1242 @keyword param: The name of the parameter to set. 1243 @type param: str 1244 @keyword value: The parameter value. 1245 @type value: anything 1246 @keyword category: The type of parameter to set. This can be 'val' for the normal parameter, 'err' for the parameter error, or 'sim' for Monte Carlo or other simulated parameters. 1247 @type category: str 1248 @keyword sim_index: The index for a Monte Carlo simulation for simulated parameter. 1249 @type sim_index: int or None 1250 """ 1251 1252 # Check the type. 1253 if category not in ['val', 'err', 'sim']: 1254 raise RelaxError("The category of the parameter '%s' is incorrectly set to %s - it must be one of 'val', 'err' or 'sim'." % (param, category)) 1255 1256 # Test if the attribute that is trying to be set is modifiable. 1257 if not param in self._mod_attr: 1258 raise RelaxError("The object '%s' is not modifiable." % param) 1259 1260 # Set a parameter value. 1261 if category == 'val': 1262 self.__dict__[param] = value 1263 1264 # Set an error. 1265 elif category == 'err': 1266 self.__dict__[param+'_err'] = value 1267 1268 # Set a simulation value. 1269 else: 1270 # Check that the simulation number has been set. 1271 if self._sim_num == None: 1272 raise RelaxError("The alignment tensor simulation number has not yet been specified, therefore a simulation value cannot be set.") 1273 1274 # The simulation parameter name. 1275 sim_param = param+'_sim' 1276 1277 # No object, so create it. 1278 if not hasattr(self, sim_param): 1279 self.__dict__[sim_param] = AlignTensorSimList(elements=self._sim_num) 1280 1281 # The object. 1282 obj = getattr(self, sim_param) 1283 1284 # Set the value. 1285 obj._set(value=value, sim_index=sim_index) 1286 1287 # Skip the updating process for certain objects. 1288 if param in ['type']: 1289 return 1290 1291 # Update the data structures. 1292 for target, update_if_set, depends in dependency_generator(): 1293 self._update_object(param, target, update_if_set, depends, category)

1294 1295

1296 - def set_fixed(self, flag):

1297 """Set if the alignment tensor should be fixed during optimisation or not. 1298 1299 @param flag: The fixed flag. 1300 @type flag: bool 1301 """ 1302 1303 self.__dict__['fixed'] = flag

1304 1305

1306 - def set_sim_num(self, sim_number=None):

1307 """Set the number of Monte Carlo simulations for the construction of the simulation structures. 1308 1309 @keyword sim_number: The number of Monte Carlo simulations. 1310 @type sim_number: int 1311 """ 1312 1313 # Check if not already set. 1314 if self._sim_num != None: 1315 raise RelaxError("The number of simulations has already been set.") 1316 1317 # Store the value. 1318 self.__dict__['_sim_num'] = sim_number

1319 1320 1321

1322 -class AlignTensorSimList(list):

1323 """Empty data container for Monte Carlo simulation alignment tensor data.""" 1324

1325 - def __init__(self, elements=None):

1326 """Initialise the Monte Carlo simulation parameter list. 1327 1328 @keyword elements: The number of elements to initialise the length of the list to. 1329 @type elements: None or int 1330 """ 1331 1332 # Initialise a length. 1333 for i in range(elements): 1334 self._append(None)

1335 1336

1337 - def __setitem__(self, slice_obj, value):

1338 """This is a read-only object!""" 1339 1340 raise RelaxError("The alignment tensor is a read-only object. The alignment tensor set() method must be used instead.")

1341 1342

1343 - def _append(self, value):

1344 """The secret append method. 1345 1346 @param value: The value to append to the list. 1347 @type value: anything 1348 """ 1349 1350 # Execute the base class method. 1351 super(AlignTensorSimList, self).append(value)

1352 1353

1354 - def _set(self, value=None, sim_index=None):

1355 """Replacement secret method for __setitem__(). 1356 1357 @keyword value: The value to set. 1358 @type value: anything 1359 @keyword sim_index: The index of the simulation value to set. 1360 @type sim_index: int 1361 """ 1362 1363 # Execute the base class method. 1364 super(AlignTensorSimList, self).__setitem__(sim_index, value)

1365 1366

1367 - def append(self, value):

1368 """This is a read-only object!""" 1369 1370 raise RelaxError("The alignment tensor is a read-only object. The alignment tensor set() method must be used instead.")

1371

Source Code for Module data.align_tensor