data.align_tensor

1 ############################################################################### 2 # # 3 # Copyright (C) 2003-2013 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 """The alignment tensor objects of the relax data store.""" 24 25 # Python module imports. 26 from re import search 27 from math import cos, sin 28 from numpy import array, dot, eye, float64, identity, transpose, zeros 29 from numpy.linalg import det, eig, eigvals 30 import sys 31 32 # relax module imports. 33 from data.data_classes import Element 34 from data.relax_xml import fill_object_contents, xml_to_object 35 from float import nan 36 from maths_fns.rotation_matrix import R_to_euler_zyz 37 from relax_errors import RelaxError 38 39

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

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

74 75

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

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

107 108

109 -def calc_A_diag(A):

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

149 150

151 -def calc_Aa(A_diag):

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

169 170

171 -def calc_Ar(A_diag):

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

189 190

191 -def calc_Axxyy(Axx, Ayy):

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

208 209

210 -def calc_Azz(Axx, Ayy):

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

227 228

229 -def calc_eigvals(A):

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

252 253

254 -def calc_eta(A_diag):

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

276 277

278 -def calc_euler(rotation):

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

288 289

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

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

324 325

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

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

357 358

359 -def calc_S_diag(tensor):

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

399 400

401 -def calc_Sxx(Axx):

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

415 416

417 -def calc_Sxxyy(Sxx, Syy):

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

434 435

436 -def calc_Sxy(Axy):

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

450 451

452 -def calc_Sxz(Axz):

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

466 467

468 -def calc_Syy(Ayy):

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

482 483

484 -def calc_Syz(Ayz):

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

498 499

500 -def calc_Szz(Sxx, Syy):

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

517 518

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

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

556 557

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

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

589 590

591 -def calc_P_diag(tensor):

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

634 635

636 -def calc_Pxx(Axx):

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

650 651

652 -def calc_Pxxyy(Pxx, Pyy):

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

669 670

671 -def calc_Pxy(Axy):

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

685 686

687 -def calc_Pxz(Axz):

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

701 702

703 -def calc_Pyy(Ayy):

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

717 718

719 -def calc_Pyz(Ayz):

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

733 734

735 -def calc_Pzz(Pxx, Pyy):

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

752 753

754 -def calc_R(Aa, Ar):

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

776 777

778 -def calc_rotation(A):

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

819 820

821 -def calc_unit_x(rotation):

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

835 836

837 -def calc_unit_y(rotation):

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

851 852

853 -def calc_unit_z(rotation):

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

867 868

869 -def dependency_generator():

870 """Generator for the automatic updating the alignment tensor data structures. 871 872 @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. 873 """ 874 875 # Primary objects (only dependant on the modifiable objects). 876 yield ('A', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 877 yield ('A_5D', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 878 yield ('Axxyy', ['Axx', 'Ayy'], ['Axx', 'Ayy']) 879 yield ('Azz', ['Axx', 'Ayy'], ['Axx', 'Ayy']) 880 881 yield ('P', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 882 yield ('P_5D', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Axx', 'Ayy', 'Azz', 'Axy', 'Axz', 'Ayz']) 883 yield ('Pxx', ['Axx'], ['Axx']) 884 yield ('Pxxyy', ['Axx', 'Ayy'], ['Axx', 'Ayy']) 885 yield ('Pxy', ['Axy'], ['Axy']) 886 yield ('Pxz', ['Axz'], ['Axz']) 887 yield ('Pyy', ['Ayy'], ['Ayy']) 888 yield ('Pyz', ['Ayz'], ['Ayz']) 889 890 yield ('Sxx', ['Axx'], ['Axx']) 891 yield ('Sxy', ['Axy'], ['Axy']) 892 yield ('Sxz', ['Axz'], ['Axz']) 893 yield ('Syy', ['Ayy'], ['Ayy']) 894 yield ('Syz', ['Ayz'], ['Ayz']) 895 896 # Secondary objects (dependant on the primary objects). 897 yield ('A_diag', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A']) 898 yield ('eigvals', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A']) 899 yield ('rotation', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A']) 900 901 yield ('P_diag', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['P']) 902 yield ('Pzz', ['Axx', 'Ayy'], ['Pxx', 'Pyy']) 903 904 yield ('S', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Sxx', 'Syy', 'Szz', 'Sxy', 'Sxz', 'Syz']) 905 yield ('S_5D', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Sxx', 'Syy', 'Szz', 'Sxy', 'Sxz', 'Syz']) 906 yield ('S_diag', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['S']) 907 yield ('Sxxyy', ['Axx', 'Ayy'], ['Sxx', 'Syy']) 908 yield ('Szz', ['Axx', 'Ayy'], ['Sxx', 'Syy']) 909 910 # Tertiary objects (dependant on the secondary objects). 911 yield ('Aa', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A_diag']) 912 yield ('Ar', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A_diag']) 913 yield ('eta', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['A_diag']) 914 915 yield ('unit_x', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 916 yield ('unit_y', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 917 yield ('unit_z', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 918 919 yield ('euler', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['rotation']) 920 921 # Quaternary objects (dependant on the tertiary objects). 922 yield ('R', ['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], ['Aa', 'Ar'])

923 924 925 926 # Alignment tensor specific data. 927 ################################# 928

929 -class AlignTensorList(list):

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

935 - def __repr__(self):

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

944 945

946 - def add_item(self, name):

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

950 951

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

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

1018 1019

1020 - def names(self):

1021 """Return a list of the alignment tensor names.""" 1022 1023 # Loop over the tensors. 1024 names = [] 1025 for i in range(len(self)): 1026 names.append(self[i].name) 1027 1028 # Return the list. 1029 return names

1030 1031

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

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

1067 1068

1069 -class AlignTensorData(Element):

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

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

1086 """Set up the tensor data. 1087 1088 @param name: The tensor ID string. 1089 @type name: str 1090 @keyword fixed: The optimisation flag. 1091 @type fixed: bool 1092 """ 1093 1094 # Store the values. 1095 self.__dict__['name'] = name 1096 self.__dict__['fixed'] = fixed 1097 1098 # The number of simulations. 1099 self.__dict__['_sim_num'] = None

1100 1101

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

1103 """Make this object read-only.""" 1104 1105 raise RelaxError("The alignment tensor is a read-only object. The alignment tensor set() method must be used instead.")

1106 1107

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

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

1243 1244

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

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

1300 1301

1302 - def set_fixed(self, flag):

1303 """Set if the alignment tensor should be fixed during optimisation or not. 1304 1305 @param flag: The fixed flag. 1306 @type flag: bool 1307 """ 1308 1309 self.__dict__['fixed'] = flag

1310 1311

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

1313 """Set the number of Monte Carlo simulations for the construction of the simulation structures. 1314 1315 @keyword sim_number: The number of Monte Carlo simulations. 1316 @type sim_number: int 1317 """ 1318 1319 # Store the value. 1320 self.__dict__['_sim_num'] = sim_number

1321 1322 1323

1324 -class AlignTensorSimList(list):

1325 """Empty data container for Monte Carlo simulation alignment tensor data.""" 1326

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

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

1337 1338

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

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

1343 1344

1345 - def _append(self, value):

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

1354 1355

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

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

1367 1368

1369 - def append(self, value):

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

1373

Source Code for Module data.align_tensor