data.align_tensor

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

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

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

73 74

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

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

106 107

108 -def calc_A_diag(A):

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

148 149

150 -def calc_Aa(A_diag):

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

168 169

170 -def calc_Ar(A_diag):

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

188 189

190 -def calc_Axxyy(Axx, Ayy):

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

207 208

209 -def calc_Azz(Axx, Ayy):

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

226 227

228 -def calc_eigvals(A):

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

251 252

253 -def calc_eta(A_diag):

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

275 276

277 -def calc_euler(rotation):

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

287 288

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

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

323 324

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

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

356 357

358 -def calc_S_diag(tensor):

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

398 399

400 -def calc_Sxx(Axx):

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

414 415

416 -def calc_Sxxyy(Sxx, Syy):

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

433 434

435 -def calc_Sxy(Axy):

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

449 450

451 -def calc_Sxz(Axz):

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

465 466

467 -def calc_Syy(Ayy):

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

481 482

483 -def calc_Syz(Ayz):

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

497 498

499 -def calc_Szz(Sxx, Syy):

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

516 517

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

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

555 556

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

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

588 589

590 -def calc_P_diag(tensor):

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

633 634

635 -def calc_Pxx(Axx):

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

649 650

651 -def calc_Pxxyy(Pxx, Pyy):

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

668 669

670 -def calc_Pxy(Axy):

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

684 685

686 -def calc_Pxz(Axz):

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

700 701

702 -def calc_Pyy(Ayy):

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

716 717

718 -def calc_Pyz(Ayz):

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

732 733

734 -def calc_Pzz(Pxx, Pyy):

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

751 752

753 -def calc_R(Aa, Ar):

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

775 776

777 -def calc_rotation(A):

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

818 819

820 -def calc_unit_x(rotation):

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

834 835

836 -def calc_unit_y(rotation):

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

850 851

852 -def calc_unit_z(rotation):

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

866 867

868 -def dependency_generator():

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

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

928 -class AlignTensorList(ListType):

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

934 - def __repr__(self):

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

943 944

945 - def add_item(self, name):

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

949 950

951 - def from_xml(self, align_tensor_super_node):

952 """Recreate the alignment tensor data structure from the XML alignment tensor node. 953 954 @param align_tensor_super_node: The alignment tensor XML nodes. 955 @type align_tensor_super_node: xml.dom.minicompat.Element instance 956 """ 957 958 # Recreate all the alignment tensor data structures. 959 xml_to_object(align_tensor_super_node, self, blacklist=['align_tensor']) 960 961 # Get the individual tensors. 962 align_tensor_nodes = align_tensor_super_node.getElementsByTagName('align_tensor') 963 964 # Loop over the child nodes. 965 for align_tensor_node in align_tensor_nodes: 966 # Add the alignment tensor data container. 967 self.add_item(align_tensor_node.getAttribute('name')) 968 969 # Recreate all the other data structures. 970 xml_to_object(align_tensor_node, self[-1])

971 972

973 - def names(self):

974 """Return a list of the alignment tensor names.""" 975 976 # Loop over the tensors. 977 names = [] 978 for i in xrange(len(self)): 979 names.append(self[i].name) 980 981 # Return the list. 982 return names

983 984

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

986 """Create an XML element for the alignment tensors. 987 988 @param doc: The XML document object. 989 @type doc: xml.dom.minidom.Document instance 990 @param element: The element to add the alignment tensors XML element to. 991 @type element: XML element object 992 """ 993 994 # Create the alignment tensors element and add it to the higher level element. 995 tensor_list_element = doc.createElement('align_tensors') 996 element.appendChild(tensor_list_element) 997 998 # Set the alignment tensor attributes. 999 tensor_list_element.setAttribute('desc', 'Alignment tensor list') 1000 1001 # Add all simple python objects within the PipeContainer to the pipe element. 1002 fill_object_contents(doc, tensor_list_element, object=self, blacklist=list(self.__class__.__dict__.keys() + list.__dict__.keys())) 1003 1004 # Loop over the tensors. 1005 for i in xrange(len(self)): 1006 # Create an XML element for a single tensor. 1007 tensor_element = doc.createElement('align_tensor') 1008 tensor_list_element.appendChild(tensor_element) 1009 tensor_element.setAttribute('index', repr(i)) 1010 tensor_element.setAttribute('desc', 'Alignment tensor') 1011 1012 # Add all simple python objects within the PipeContainer to the pipe element. 1013 fill_object_contents(doc, tensor_element, object=self[i], blacklist=list(self[i].__class__.__dict__.keys()))

1014 1015

1016 -class AlignTensorData(Element):

1017 """An empty data container for the alignment tensor elements.""" 1018 1019 # List of modifiable attributes. 1020 __mod_attr__ = ['name', 1021 'Axx', 'Axx_sim', 'Axx_err', 1022 'Ayy', 'Ayy_sim', 'Ayy_err', 1023 'Axy', 'Axy_sim', 'Axy_err', 1024 'Axz', 'Axz_sim', 'Axz_err', 1025 'Ayz', 'Ayz_sim', 'Ayz_err', 1026 'domain', 1027 'red', 1028 'fixed'] 1029

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

1031 """Set up the tensor data. 1032 1033 @param name: The tensor ID string. 1034 @type name: str 1035 @keyword fixed: The optimisation flag. 1036 @type fixed: bool 1037 """ 1038 1039 # Store the values. 1040 self.name = name 1041 self.fixed = fixed

1042 1043

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

1045 """Function for calculating the parameters, unit vectors, and tensors on the fly. 1046 1047 @param name: The name of the object to set. 1048 @type name: str 1049 @param value: The value to set the object corresponding to the name argument to. 1050 @type value: Any Python object type 1051 """ 1052 1053 # Get the base parameter name and determine the object category ('val', 'err', or 'sim'). 1054 if search('_err$', name): 1055 category = 'err' 1056 param_name = name[:-4] 1057 elif search('_sim$', name): 1058 category = 'sim' 1059 param_name = name[:-4] 1060 else: 1061 category = 'val' 1062 param_name = name 1063 1064 # Test if the attribute that is trying to be set is modifiable. 1065 if not param_name in self.__mod_attr__: 1066 raise RelaxError("The object " + repr(name) + " is not modifiable.") 1067 1068 # Set the attribute normally. 1069 self.__dict__[name] = value 1070 1071 # Update the data structures. 1072 for target, update_if_set, depends in dependency_generator(): 1073 self.__update_object(param_name, target, update_if_set, depends, category)

1074 1075

1076 - def __update_sim_append(self, param_name, index):

1077 """Update the Monte Carlo simulation data lists when a simulation value is appended. 1078 1079 @param param_name: The MC sim parameter name which is being appended to. 1080 @type param_name: str 1081 @param index: The index of the Monte Carlo simulation which was set. 1082 @type index: int 1083 """ 1084 1085 # Loop over the targets. 1086 for target, update_if_set, depends in dependency_generator(): 1087 # Only update if the parameter name is within the 'update_if_set' list. 1088 if not param_name in update_if_set: 1089 continue 1090 1091 # Get the function for calculating the value. 1092 fn = globals()['calc_'+target] 1093 1094 # Get all the dependencies if possible. 1095 missing_dep = 0 1096 deps = () 1097 for dep_name in depends: 1098 # Modify the dependency name. 1099 if dep_name != 'type': 1100 dep_name = dep_name+'_sim' 1101 1102 # Test if the MC sim object exists. 1103 if not hasattr(self, dep_name): 1104 missing_dep = 1 1105 break 1106 1107 # Get the MC dependency. 1108 dep_obj = getattr(self, dep_name) 1109 1110 # The alignment tensor type. 1111 if dep_name == 'type': 1112 deps = deps+(dep_obj,) 1113 continue 1114 1115 # Test if the MC sim dependency is long enough. 1116 if len(dep_obj) <= index: 1117 missing_dep = 1 1118 break 1119 1120 # Place the value corresponding to the index into the 'deps' array. 1121 deps = deps+(dep_obj[index],) 1122 1123 # Only update the MC simulation object if its dependencies exist. 1124 if not missing_dep: 1125 # Get the target object. 1126 target_obj = getattr(self, target+'_sim') 1127 1128 # Calculate and set the value. 1129 target_obj.append_untouchable_item(fn(*deps))

1130 1131

1132 - def __update_sim_set(self, param_name, index):

1133 """Update the Monte Carlo simulation data lists when a simulation value is set. 1134 1135 @param param_name: The MC sim parameter name which is being set. 1136 @type param_name: str 1137 @param index: The index of the Monte Carlo simulation which was set. 1138 @type index: int 1139 """ 1140 1141 # Loop over the targets. 1142 for target, update_if_set, depends in dependency_generator(): 1143 # Only update if the parameter name is within the 'update_if_set' list. 1144 if not param_name in update_if_set: 1145 continue 1146 1147 # Get the function for calculating the value. 1148 fn = globals()['calc_'+target] 1149 1150 # Get all the dependencies if possible. 1151 missing_dep = 0 1152 deps = () 1153 for dep_name in depends: 1154 # Modify the dependency name. 1155 if dep_name != 'type': 1156 dep_name = dep_name+'_sim' 1157 1158 # Test if the MC sim object exists. 1159 if not hasattr(self, dep_name): 1160 missing_dep = 1 1161 break 1162 1163 # Get the MC dependency. 1164 dep_obj = getattr(self, dep_name) 1165 1166 # The alignment tensor type. 1167 if dep_name == 'type': 1168 deps = deps+(dep_obj,) 1169 continue 1170 1171 # Test if the MC sim dependency is long enough. 1172 if len(dep_obj) <= index: 1173 missing_dep = 1 1174 break 1175 1176 # Place the value corresponding to the index into the 'deps' array. 1177 deps = deps+(dep_obj[index],) 1178 1179 # Only update the MC simulation object if its dependencies exist. 1180 if not missing_dep: 1181 # Get the target object. 1182 target_obj = getattr(self, target+'_sim') 1183 1184 # Calculate and set the value. 1185 target_obj.set_untouchable_item(index, fn(*deps))

1186 1187

1188 - def __update_object(self, param_name, target, update_if_set, depends, category):

1189 """Function for updating the target object, its error, and the MC simulations. 1190 1191 If the base name of the object is not within the 'update_if_set' list, this function returns 1192 without doing anything (to avoid wasting time). Dependant upon the category the object 1193 (target), its error (target+'_err'), or all Monte Carlo simulations (target+'_sim') are 1194 updated. 1195 1196 @param param_name: The parameter name which is being set in the __setattr__() function. 1197 @type param_name: str 1198 @param target: The name of the object to update. 1199 @type target: str 1200 @param update_if_set: If the parameter being set by the __setattr__() function is not 1201 within this list of parameters, don't waste time updating the 1202 target. 1203 @param depends: An array of names objects that the target is dependent upon. 1204 @type depends: array of str 1205 @param category: The category of the object to update (one of 'val', 'err', or 1206 'sim'). 1207 @type category: str 1208 @return: None 1209 """ 1210 1211 # Only update if the parameter name is within the 'update_if_set' list. 1212 if not param_name in update_if_set: 1213 return 1214 1215 # Get the function for calculating the value. 1216 fn = globals()['calc_'+target] 1217 1218 1219 # The value. 1220 ############ 1221 1222 if category == 'val': 1223 # Get all the dependencies if possible. 1224 missing_dep = 0 1225 deps = () 1226 for dep_name in depends: 1227 # Test if the object exists. 1228 if not hasattr(self, dep_name): 1229 missing_dep = 1 1230 break 1231 1232 # Get the object and place it into the 'deps' tuple. 1233 deps = deps+(getattr(self, dep_name),) 1234 1235 # Only update the object if its dependencies exist. 1236 if not missing_dep: 1237 # Calculate the value. 1238 value = fn(*deps) 1239 1240 # Set the attribute. 1241 self.__dict__[target] = value 1242 1243 1244 # The error. 1245 ############ 1246 1247 if category == 'err': 1248 # Get all the dependencies if possible. 1249 missing_dep = 0 1250 deps = () 1251 for dep_name in depends: 1252 # Test if the error object exists. 1253 if not hasattr(self, dep_name+'_err'): 1254 missing_dep = 1 1255 break 1256 1257 # Get the object and place it into the 'deps' tuple. 1258 deps = deps+(getattr(self, dep_name+'_err'),) 1259 1260 # Only update the error object if its dependencies exist. 1261 if not missing_dep: 1262 # Calculate the value. 1263 value = fn(*deps) 1264 1265 # Set the attribute. 1266 self.__dict__[target+'_err'] = value 1267 1268 1269 # The Monte Carlo simulations. 1270 ############################## 1271 1272 if category == 'sim': 1273 # Get all the dependencies if possible. 1274 missing_dep = 0 1275 deps = [] 1276 for dep_name in depends: 1277 # Modify the dependency name. 1278 if dep_name != 'type': 1279 dep_name = dep_name+'_sim' 1280 1281 # Test if the MC sim object exists. 1282 if not hasattr(self, dep_name): 1283 missing_dep = 1 1284 break 1285 1286 # Only create the MC simulation object if its dependencies exist. 1287 if not missing_dep: 1288 # Initialise an empty array to store the MC simulation object elements (if it doesn't already exist). 1289 if not target+'_sim' in self.__dict__: 1290 self.__dict__[target+'_sim'] = AlignTensorSimList(target, self)

1291 1292 1293

1294 -class AlignTensorSimList(ListType):

1295 """Empty data container for Monte Carlo simulation alignment tensor data.""" 1296

1297 - def __init__(self, param_name, align_element):

1298 """Initialise the Monte Carlo simulation parameter list. 1299 1300 This function makes the parameter name and parent object accessible to the functions of this 1301 list object. 1302 """ 1303 1304 self.param_name = param_name 1305 self.align_element = align_element

1306 1307

1308 - def __setitem__(self, index, value):

1309 """Set the value.""" 1310 1311 # Set the value. 1312 ListType.__setitem__(self, index, value) 1313 1314 # Then update the other lists. 1315 self.align_element._AlignTensorData__update_sim_set(self.param_name, index)

1316 1317

1318 - def append(self, value):

1319 """Replacement function for the normal self.append() method.""" 1320 1321 # Append the value to the list. 1322 self[len(self):len(self)] = [value] 1323 1324 # Update the other MC lists. 1325 self.align_element._AlignTensorData__update_sim_append(self.param_name, len(self)-1)

1326 1327

1328 - def append_untouchable_item(self, value):

1329 """Append the value for an untouchable MC data structure.""" 1330 1331 # Append the value to the list. 1332 self[len(self):len(self)] = [value]

1333 1334

1335 - def set_untouchable_item(self, index, value):

1336 """Set the value for an untouchable MC data structure.""" 1337 1338 # Set the value. 1339 ListType.__setitem__(self, index, value)

1340

Source Code for Module data.align_tensor