Source Code for Module maths_fns.frame_order_matrix_ops

   1  ############################################################################### 
   2  #                                                                             # 
   3  # Copyright (C) 2009-2011 Edward d'Auvergne                                   # 
   4  #                                                                             # 
   5  # This file is part of the program relax (http://www.nmr-relax.com).          # 
   6  #                                                                             # 
   7  # This program is free software: you can redistribute it and/or modify        # 
   8  # it under the terms of the GNU General Public License as published by        # 
   9  # the Free Software Foundation, either version 3 of the License, or           # 
  10  # (at your option) any later version.                                         # 
  11  #                                                                             # 
  12  # This program is distributed in the hope that it will be useful,             # 
  13  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
  14  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
  15  # GNU General Public License for more details.                                # 
  16  #                                                                             # 
  17  # You should have received a copy of the GNU General Public License           # 
  18  # along with this program.  If not, see <http://www.gnu.org/licenses/>.       # 
  19  #                                                                             # 
  20  ############################################################################### 
  21   
  22  # Module docstring. 
  23  """Module for the handling of Frame Order.""" 
  24   
  25  # Dependency check module. 
  26  import dep_check 
  27   
  28  # Python module imports. 
  29  from math import cos, pi, sin, sqrt 
  30  from numpy import cross, dot, sinc, transpose 
  31  from numpy.linalg import norm 
  32  if dep_check.scipy_module: 
  33      from scipy.integrate import quad 
  34   
  35  # relax module imports. 
  36  from float import isNaN 
  37  from maths_fns import order_parameters 
  38  from maths_fns.coord_transform import spherical_to_cartesian 
  39  from maths_fns.kronecker_product import kron_prod, transpose_23 
  40  from maths_fns.pseudo_ellipse import pec 
  41  from maths_fns.rotation_matrix import euler_to_R_zyz, two_vect_to_R 
  42   
  43   
  44 -def compile_1st_matrix_pseudo_ellipse(matrix, theta_x, theta_y, sigma_max): 
  45      """Generate the 1st degree Frame Order matrix for the pseudo-ellipse. 
  46   
  47      @param matrix:      The Frame Order matrix, 1st degree to be populated. 
  48      @type matrix:       numpy 3D, rank-2 array 
  49      @param theta_x:     The cone opening angle along x. 
  50      @type theta_x:      float 
  51      @param theta_y:     The cone opening angle along y. 
  52      @type theta_y:      float 
  53      @param sigma_max:   The maximum torsion angle. 
  54      @type sigma_max:    float 
  55      """ 
  56   
  57      # The surface area normalisation factor. 
  58      fact = 1.0 / (2.0 * sigma_max * pec(theta_x, theta_y)) 
  59   
  60      # Numerical integration of phi of each element. 
  61      matrix[0, 0] = fact * quad(part_int_daeg1_pseudo_ellipse_xx, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0] 
  62      matrix[1, 1] = fact * quad(part_int_daeg1_pseudo_ellipse_yy, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0] 
  63      matrix[2, 2] = fact * quad(part_int_daeg1_pseudo_ellipse_zz, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0] 
  64   
  65   
  66 -def compile_2nd_matrix_free_rotor(matrix, R, z_axis, axis, theta_axis, phi_axis): 
  67      """Generate the rotated 2nd degree Frame Order matrix for the free rotor model. 
  68   
  69      The rotor axis is assumed to be parallel to the z-axis in the eigenframe. 
  70   
  71   
  72      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
  73      @type matrix:       numpy 9D, rank-2 array 
  74      @param R:           The rotation matrix to be populated. 
  75      @type R:            numpy 3D, rank-2 array 
  76      @param z_axis:      The molecular frame z-axis from which the rotor axis is rotated from. 
  77      @type z_axis:       numpy 3D, rank-1 array 
  78      @param axis:        The storage structure for the axis. 
  79      @type axis:         numpy 3D, rank-1 array 
  80      @param theta_axis:  The axis polar angle. 
  81      @type theta_axis:   float 
  82      @param phi_axis:    The axis azimuthal angle. 
  83      @type phi_axis:     float 
  84      """ 
  85   
  86      # Zeros. 
  87      for i in range(9): 
  88          for j in range(9): 
  89              matrix[i, j] = 0.0 
  90   
  91      # Diagonal. 
  92      matrix[0, 0] = 0.5 
  93      matrix[1, 1] = 0.5 
  94      matrix[3, 3] = 0.5 
  95      matrix[4, 4] = 0.5 
  96      matrix[8, 8] = 1 
  97   
  98      # Off diagonal set 1. 
  99      matrix[0, 4] = matrix[4, 0] = 0.5 
 100   
 101      # Off diagonal set 2. 
 102      matrix[1, 3] = matrix[3, 1] = -0.5 
 103   
 104      # Generate the cone axis from the spherical angles. 
 105      spherical_to_cartesian([1.0, theta_axis, phi_axis], axis) 
 106   
 107      # Average position rotation. 
 108      two_vect_to_R(z_axis, axis, R) 
 109   
 110      # Rotate and return the frame order matrix. 
 111      return rotate_daeg(matrix, R) 
 112   
 113   
 114 -def compile_2nd_matrix_iso_cone(matrix, R, eigen_alpha, eigen_beta, eigen_gamma, cone_theta, sigma_max): 
 115      """Generate the rotated 2nd degree Frame Order matrix for the isotropic cone. 
 116   
 117      The cone axis is assumed to be parallel to the z-axis in the eigenframe. 
 118   
 119      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 120      @type matrix:       numpy 9D, rank-2 array 
 121      @param R:           The rotation matrix to be populated. 
 122      @type R:            numpy 3D, rank-2 array 
 123      @param eigen_alpha: The eigenframe rotation alpha Euler angle. 
 124      @type eigen_alpha:  float 
 125      @param eigen_beta:  The eigenframe rotation beta Euler angle. 
 126      @type eigen_beta:   float 
 127      @param eigen_gamma: The eigenframe rotation gamma Euler angle. 
 128      @type eigen_gamma:  float 
 129      @param cone_theta:  The cone opening angle. 
 130      @type cone_theta:   float 
 131      @param sigma_max:   The maximum torsion angle. 
 132      @type sigma_max:    float 
 133      """ 
 134   
 135      # Populate the Frame Order matrix in the eigenframe. 
 136      populate_2nd_eigenframe_iso_cone(matrix, cone_theta, sigma_max) 
 137   
 138      # Average position rotation. 
 139      euler_to_R_zyz(eigen_alpha, eigen_beta, eigen_gamma, R) 
 140   
 141      # Rotate and return the frame order matrix. 
 142      return rotate_daeg(matrix, R) 
 143   
 144   
 145 -def compile_2nd_matrix_iso_cone_free_rotor(matrix, R, z_axis, cone_axis, theta_axis, phi_axis, s1): 
 146      """Generate the rotated 2nd degree Frame Order matrix for the free rotor isotropic cone. 
 147   
 148      The cone axis is assumed to be parallel to the z-axis in the eigenframe.  In this model, the three order parameters are defined as:: 
 149   
 150          S1 = S2, 
 151          S3 = 0 
 152   
 153   
 154      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 155      @type matrix:       numpy 9D, rank-2 array 
 156      @param R:           The rotation matrix to be populated. 
 157      @type R:            numpy 3D, rank-2 array 
 158      @param z_axis:      The molecular frame z-axis from which the cone axis is rotated from. 
 159      @type z_axis:       numpy 3D, rank-1 array 
 160      @param cone_axis:   The storage structure for the cone axis. 
 161      @type cone_axis:    numpy 3D, rank-1 array 
 162      @param theta_axis:  The cone axis polar angle. 
 163      @type theta_axis:   float 
 164      @param phi_axis:    The cone axis azimuthal angle. 
 165      @type phi_axis:     float 
 166      @param s1:          The cone order parameter. 
 167      @type s1:           float 
 168      """ 
 169   
 170      # Generate the cone axis from the spherical angles. 
 171      spherical_to_cartesian([1.0, theta_axis, phi_axis], cone_axis) 
 172   
 173      # Populate the Frame Order matrix in the eigenframe. 
 174      populate_2nd_eigenframe_iso_cone_free_rotor(matrix, s1) 
 175   
 176      # Average position rotation. 
 177      two_vect_to_R(z_axis, cone_axis, R) 
 178   
 179      # Rotate and return the frame order matrix. 
 180      return rotate_daeg(matrix, R) 
 181   
 182   
 183 -def compile_2nd_matrix_iso_cone_torsionless(matrix, R, z_axis, cone_axis, theta_axis, phi_axis, cone_theta): 
 184      """Generate the rotated 2nd degree Frame Order matrix for the torsionless isotropic cone. 
 185   
 186      The cone axis is assumed to be parallel to the z-axis in the eigenframe. 
 187   
 188   
 189      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 190      @type matrix:       numpy 9D, rank-2 array 
 191      @param R:           The rotation matrix to be populated. 
 192      @type R:            numpy 3D, rank-2 array 
 193      @param z_axis:      The molecular frame z-axis from which the cone axis is rotated from. 
 194      @type z_axis:       numpy 3D, rank-1 array 
 195      @param cone_axis:   The storage structure for the cone axis. 
 196      @type cone_axis:    numpy 3D, rank-1 array 
 197      @param theta_axis:  The cone axis polar angle. 
 198      @type theta_axis:   float 
 199      @param phi_axis:    The cone axis azimuthal angle. 
 200      @type phi_axis:     float 
 201      @param cone_theta:  The cone opening angle. 
 202      @type cone_theta:   float 
 203      """ 
 204   
 205      # Zeros. 
 206      for i in range(9): 
 207          for j in range(9): 
 208              matrix[i, j] = 0.0 
 209   
 210      # Repetitive trig calculations. 
 211      cos_tmax = cos(cone_theta) 
 212      cos_tmax2 = cos_tmax**2 
 213   
 214      # Diagonal. 
 215      matrix[0, 0] = (3.0*cos_tmax2 + 6.0*cos_tmax + 15.0) / 24.0 
 216      matrix[1, 1] = (cos_tmax2 + 10.0*cos_tmax + 13.0) / 24.0 
 217      matrix[2, 2] = (4.0*cos_tmax2 + 10.0*cos_tmax + 10.0) / 24.0 
 218      matrix[3, 3] = matrix[1, 1] 
 219      matrix[4, 4] = matrix[0, 0] 
 220      matrix[5, 5] = matrix[2, 2] 
 221      matrix[6, 6] = matrix[2, 2] 
 222      matrix[7, 7] = matrix[2, 2] 
 223      matrix[8, 8] = (cos_tmax2 + cos_tmax + 1.0) / 3.0 
 224   
 225      # Off diagonal set 1. 
 226      matrix[0, 4] = matrix[4, 0] = (cos_tmax2 - 2.0*cos_tmax + 1.0) / 24.0 
 227      matrix[0, 8] = matrix[8, 0] = -(cos_tmax2 + cos_tmax - 2.0) / 6.0 
 228      matrix[4, 8] = matrix[8, 4] = matrix[0, 8] 
 229   
 230      # Off diagonal set 2. 
 231      matrix[1, 3] = matrix[3, 1] = matrix[0, 4] 
 232      matrix[2, 6] = matrix[6, 2] = -matrix[0, 8] 
 233      matrix[5, 7] = matrix[7, 5] = -matrix[0, 8] 
 234   
 235      # Generate the cone axis from the spherical angles. 
 236      spherical_to_cartesian([1.0, theta_axis, phi_axis], cone_axis) 
 237   
 238      # Average position rotation. 
 239      two_vect_to_R(z_axis, cone_axis, R) 
 240   
 241      # Rotate and return the frame order matrix. 
 242      return rotate_daeg(matrix, R) 
 243   
 244   
 245 -def compile_2nd_matrix_pseudo_ellipse(matrix, R, eigen_alpha, eigen_beta, eigen_gamma, theta_x, theta_y, sigma_max): 
 246      """Generate the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 247   
 248      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 249      @type matrix:       numpy 9D, rank-2 array 
 250      @param R:           The rotation matrix to be populated. 
 251      @type R:            numpy 3D, rank-2 array 
 252      @param eigen_alpha: The eigenframe rotation alpha Euler angle. 
 253      @type eigen_alpha:  float 
 254      @param eigen_beta:  The eigenframe rotation beta Euler angle. 
 255      @type eigen_beta:   float 
 256      @param eigen_gamma: The eigenframe rotation gamma Euler angle. 
 257      @type eigen_gamma:  float 
 258      @param theta_x:     The cone opening angle along x. 
 259      @type theta_x:      float 
 260      @param theta_y:     The cone opening angle along y. 
 261      @type theta_y:      float 
 262      @param sigma_max:   The maximum torsion angle. 
 263      @type sigma_max:    float 
 264      """ 
 265   
 266      # The surface area normalisation factor. 
 267      fact = 12.0 * pec(theta_x, theta_y) 
 268   
 269      # Invert. 
 270      if fact == 0.0: 
 271          fact = 1e100 
 272      else: 
 273          fact = 1.0 / fact 
 274   
 275      # Sigma_max part. 
 276      if sigma_max == 0.0: 
 277          fact2 = 1e100 
 278      else: 
 279          fact2 = fact / (2.0 * sigma_max) 
 280   
 281      # Diagonal. 
 282      matrix[0, 0] = fact * (4.0*pi*(sinc(2.0*sigma_max/pi) + 2.0) + quad(part_int_daeg2_pseudo_ellipse_00, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 283      matrix[1, 1] = fact * (4.0*pi*sinc(2.0*sigma_max/pi) + quad(part_int_daeg2_pseudo_ellipse_11, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 284      matrix[2, 2] = fact * 2.0*sinc(sigma_max/pi) * (5.0*pi - quad(part_int_daeg2_pseudo_ellipse_22, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 285      matrix[3, 3] = matrix[1, 1] 
 286      matrix[4, 4] = fact * (4.0*pi*(sinc(2.0*sigma_max/pi) + 2.0) + quad(part_int_daeg2_pseudo_ellipse_44, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 287      matrix[5, 5] = fact * 2.0*sinc(sigma_max/pi) * (5.0*pi - quad(part_int_daeg2_pseudo_ellipse_55, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 288      matrix[6, 6] = matrix[2, 2] 
 289      matrix[7, 7] = matrix[5, 5] 
 290      matrix[8, 8] = 4.0 * fact * (2.0*pi - quad(part_int_daeg2_pseudo_ellipse_88, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 291   
 292      # Off diagonal set 1. 
 293      matrix[0, 4] = fact * (4.0*pi*(2.0 - sinc(2.0*sigma_max/pi)) + quad(part_int_daeg2_pseudo_ellipse_04, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 294      matrix[4, 0] = fact * (4.0*pi*(2.0 - sinc(2.0*sigma_max/pi)) + quad(part_int_daeg2_pseudo_ellipse_40, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 295      matrix[0, 8] = 4.0 * fact * (2.0*pi - quad(part_int_daeg2_pseudo_ellipse_08, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 296      matrix[8, 0] = fact * (8.0*pi + quad(part_int_daeg2_pseudo_ellipse_80, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 297      matrix[4, 8] = 4.0 * fact * (2.0*pi - quad(part_int_daeg2_pseudo_ellipse_48, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 298      matrix[8, 4] = fact * (8.0*pi - quad(part_int_daeg2_pseudo_ellipse_84, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 299   
 300      # Off diagonal set 2. 
 301      matrix[1, 3] = matrix[3, 1] = fact * (4.0*pi*sinc(2.0*sigma_max/pi) + quad(part_int_daeg2_pseudo_ellipse_13, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 302      matrix[2, 6] = matrix[6, 2] = -fact * 4.0 * sinc(sigma_max/pi) * (2.0*pi + quad(part_int_daeg2_pseudo_ellipse_26, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 303      matrix[5, 7] = matrix[7, 5] = -fact * 4.0 * sinc(sigma_max/pi) * (2.0*pi + quad(part_int_daeg2_pseudo_ellipse_57, -pi, pi, args=(theta_x, theta_y, sigma_max), full_output=1)[0]) 
 304   
 305      # Average position rotation. 
 306      euler_to_R_zyz(eigen_alpha, eigen_beta, eigen_gamma, R) 
 307   
 308      # Rotate and return the frame order matrix. 
 309      return rotate_daeg(matrix, R) 
 310   
 311   
 312 -def compile_2nd_matrix_pseudo_ellipse_free_rotor(matrix, R, eigen_alpha, eigen_beta, eigen_gamma, theta_x, theta_y): 
 313      """Generate the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
 314   
 315      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 316      @type matrix:       numpy 9D, rank-2 array 
 317      @param R:           The rotation matrix to be populated. 
 318      @type R:            numpy 3D, rank-2 array 
 319      @param eigen_alpha: The eigenframe rotation alpha Euler angle. 
 320      @type eigen_alpha:  float 
 321      @param eigen_beta:  The eigenframe rotation beta Euler angle. 
 322      @type eigen_beta:   float 
 323      @param eigen_gamma: The eigenframe rotation gamma Euler angle. 
 324      @type eigen_gamma:  float 
 325      @param theta_x:     The cone opening angle along x. 
 326      @type theta_x:      float 
 327      @param theta_y:     The cone opening angle along y. 
 328      @type theta_y:      float 
 329      """ 
 330   
 331      # The surface area normalisation factor. 
 332      fact = 1.0 / (6.0 * pec(theta_x, theta_y)) 
 333   
 334      # Diagonal. 
 335      matrix[0, 0] = fact * (4.0*pi - quad(part_int_daeg2_pseudo_ellipse_free_rotor_00, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 336      matrix[1, 1] = matrix[3, 3] = fact * 3.0/2.0 * quad(part_int_daeg2_pseudo_ellipse_free_rotor_11, -pi, pi, args=(theta_x, theta_y), full_output=1)[0] 
 337      matrix[4, 4] = fact * (4.0*pi - quad(part_int_daeg2_pseudo_ellipse_free_rotor_44, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 338      matrix[8, 8] = fact * (4.0*pi - 2.0*quad(part_int_daeg2_pseudo_ellipse_free_rotor_88, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 339   
 340      # Off diagonal set 1. 
 341      matrix[0, 4] = matrix[0, 0] 
 342      matrix[4, 0] = matrix[4, 4] 
 343      matrix[0, 8] = fact * (4.0*pi + quad(part_int_daeg2_pseudo_ellipse_free_rotor_08, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 344      matrix[8, 0] = fact * (4.0*pi + quad(part_int_daeg2_pseudo_ellipse_free_rotor_80, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 345      matrix[4, 8] = fact * (4.0*pi + quad(part_int_daeg2_pseudo_ellipse_free_rotor_48, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 346      matrix[8, 4] = matrix[8, 0] 
 347   
 348      # Off diagonal set 2. 
 349      matrix[1, 3] = matrix[3, 1] = -matrix[1, 1] 
 350   
 351      # Average position rotation. 
 352      euler_to_R_zyz(eigen_alpha, eigen_beta, eigen_gamma, R) 
 353   
 354      # Rotate and return the frame order matrix. 
 355      return rotate_daeg(matrix, R) 
 356   
 357   
 358 -def compile_2nd_matrix_pseudo_ellipse_torsionless(matrix, R, eigen_alpha, eigen_beta, eigen_gamma, theta_x, theta_y): 
 359      """Generate the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
 360   
 361      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 362      @type matrix:       numpy 9D, rank-2 array 
 363      @param R:           The rotation matrix to be populated. 
 364      @type R:            numpy 3D, rank-2 array 
 365      @param eigen_alpha: The eigenframe rotation alpha Euler angle. 
 366      @type eigen_alpha:  float 
 367      @param eigen_beta:  The eigenframe rotation beta Euler angle. 
 368      @type eigen_beta:   float 
 369      @param eigen_gamma: The eigenframe rotation gamma Euler angle. 
 370      @type eigen_gamma:  float 
 371      @param theta_x:     The cone opening angle along x. 
 372      @type theta_x:      float 
 373      @param theta_y:     The cone opening angle along y. 
 374      @type theta_y:      float 
 375      """ 
 376   
 377      # The surface area normalisation factor. 
 378      fact = 1.0 / (6.0 * pec(theta_x, theta_y)) 
 379   
 380      # Diagonal. 
 381      matrix[0, 0] = fact * (6.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_00, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 382      matrix[1, 1] = fact * (2.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_11, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 383      matrix[2, 2] = fact * (5.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_22, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 384      matrix[3, 3] = matrix[1, 1] 
 385      matrix[4, 4] = fact * (6.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_44, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 386      matrix[5, 5] = fact * (5.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_55, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 387      matrix[6, 6] = matrix[2, 2] 
 388      matrix[7, 7] = matrix[5, 5] 
 389      matrix[8, 8] = fact * quad(part_int_daeg2_pseudo_ellipse_torsionless_88, -pi, pi, args=(theta_x, theta_y), full_output=1)[0] 
 390   
 391      # Off diagonal set 1. 
 392      matrix[0, 4] = matrix[4, 0] = fact * (2.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_04, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 393      matrix[0, 8] = matrix[8, 0] = fact * (4.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_08, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 394      matrix[4, 8] = matrix[8, 4] = fact * (4.0*pi + quad(part_int_daeg2_pseudo_ellipse_torsionless_48, -pi, pi, args=(theta_x, theta_y), full_output=1)[0]) 
 395   
 396      # Off diagonal set 2. 
 397      matrix[1, 3] = matrix[3, 1] = matrix[0, 4] 
 398      matrix[2, 6] = matrix[6, 2] = -matrix[0, 8] 
 399      matrix[5, 7] = matrix[7, 5] = -matrix[4, 8] 
 400   
 401      # Average position rotation. 
 402      euler_to_R_zyz(eigen_alpha, eigen_beta, eigen_gamma, R) 
 403   
 404      # Rotate and return the frame order matrix. 
 405      return rotate_daeg(matrix, R) 
 406   
 407   
 408 -def compile_2nd_matrix_rotor(matrix, R, z_axis, axis, theta_axis, phi_axis, smax): 
 409      """Generate the rotated 2nd degree Frame Order matrix for the rotor model. 
 410   
 411      The cone axis is assumed to be parallel to the z-axis in the eigenframe. 
 412   
 413   
 414      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
 415      @type matrix:       numpy 9D, rank-2 array 
 416      @param R:           The rotation matrix to be populated. 
 417      @type R:            numpy 3D, rank-2 array 
 418      @param z_axis:      The molecular frame z-axis from which the rotor axis is rotated from. 
 419      @type z_axis:       numpy 3D, rank-1 array 
 420      @param axis:        The storage structure for the axis. 
 421      @type axis:         numpy 3D, rank-1 array 
 422      @param theta_axis:  The axis polar angle. 
 423      @type theta_axis:   float 
 424      @param phi_axis:    The axis azimuthal angle. 
 425      @type phi_axis:     float 
 426      @param smax:        The maximum torsion angle. 
 427      @type smax:         float 
 428      """ 
 429   
 430      # Zeros. 
 431      for i in range(9): 
 432          for j in range(9): 
 433              matrix[i, j] = 0.0 
 434   
 435      # Repetitive trig calculations. 
 436      sinc_smax = sinc(smax/pi) 
 437      sinc_2smax = sinc(2.0*smax/pi) 
 438   
 439      # Diagonal. 
 440      matrix[0, 0] = (sinc_2smax + 1.0) / 2.0 
 441      matrix[1, 1] = matrix[0, 0] 
 442      matrix[2, 2] = sinc_smax 
 443      matrix[3, 3] = matrix[0, 0] 
 444      matrix[4, 4] = matrix[0, 0] 
 445      matrix[5, 5] = matrix[2, 2] 
 446      matrix[6, 6] = matrix[2, 2] 
 447      matrix[7, 7] = matrix[2, 2] 
 448      matrix[8, 8] = 1.0 
 449   
 450      # Off diagonal set 1. 
 451      matrix[0, 4] = matrix[4, 0] = -(sinc_2smax - 1.0) / 2.0 
 452   
 453      # Off diagonal set 2. 
 454      matrix[1, 3] = matrix[3, 1] = -matrix[0, 4] 
 455   
 456      # Generate the cone axis from the spherical angles. 
 457      spherical_to_cartesian([1.0, theta_axis, phi_axis], axis) 
 458   
 459      # Average position rotation. 
 460      two_vect_to_R(z_axis, axis, R) 
 461   
 462      # Rotate and return the frame order matrix. 
 463      return rotate_daeg(matrix, R) 
 464   
 465   
 466 -def daeg_to_rotational_superoperator(daeg, Rsuper): 
 467      """Convert the frame order matrix (daeg) to the rotational superoperator. 
 468   
 469      @param daeg:    The second degree frame order matrix, daeg.  This must be in the Kronecker product layout. 
 470      @type daeg:     numpy 9D, rank-2 array or numpy 3D, rank-4 array 
 471      @param Rsuper:  The rotational superoperator structure to be populated. 
 472      @type Rsuper:   numpy 5D, rank-2 array 
 473      """ 
 474   
 475      # First perform the T23 transpose. 
 476      transpose_23(daeg) 
 477   
 478      # Convert to rank-4. 
 479      orig_shape = daeg.shape 
 480      daeg.shape = (3, 3, 3, 3) 
 481   
 482      # First column of the superoperator. 
 483      Rsuper[0, 0] = daeg[0, 0, 0, 0] - daeg[2, 0, 2, 0] 
 484      Rsuper[1, 0] = daeg[0, 1, 0, 1] - daeg[2, 1, 2, 1] 
 485      Rsuper[2, 0] = daeg[0, 0, 0, 1] - daeg[2, 0, 2, 1] 
 486      Rsuper[3, 0] = daeg[0, 0, 0, 2] - daeg[2, 0, 2, 2] 
 487      Rsuper[4, 0] = daeg[0, 1, 0, 2] - daeg[2, 1, 2, 2] 
 488   
 489      # Second column of the superoperator. 
 490      Rsuper[0, 1] = daeg[1, 0, 1, 0] - daeg[2, 0, 2, 0] 
 491      Rsuper[1, 1] = daeg[1, 1, 1, 1] - daeg[2, 1, 2, 1] 
 492      Rsuper[2, 1] = daeg[1, 0, 1, 1] - daeg[2, 0, 2, 1] 
 493      Rsuper[3, 1] = daeg[1, 0, 1, 2] - daeg[2, 0, 2, 2] 
 494      Rsuper[4, 1] = daeg[1, 1, 1, 2] - daeg[2, 1, 2, 2] 
 495   
 496      # Third column of the superoperator. 
 497      Rsuper[0, 2] = daeg[0, 0, 1, 0] + daeg[1, 0, 0, 0] 
 498      Rsuper[1, 2] = daeg[0, 1, 1, 1] + daeg[1, 1, 0, 1] 
 499      Rsuper[2, 2] = daeg[0, 0, 1, 1] + daeg[1, 0, 0, 1] 
 500      Rsuper[3, 2] = daeg[0, 0, 1, 2] + daeg[1, 0, 0, 2] 
 501      Rsuper[4, 2] = daeg[0, 1, 1, 2] + daeg[1, 1, 0, 2] 
 502   
 503      # Fourth column of the superoperator. 
 504      Rsuper[0, 3] = daeg[0, 0, 2, 0] + daeg[2, 0, 0, 0] 
 505      Rsuper[1, 3] = daeg[0, 1, 2, 1] + daeg[2, 1, 0, 1] 
 506      Rsuper[2, 3] = daeg[0, 0, 2, 1] + daeg[2, 0, 0, 1] 
 507      Rsuper[3, 3] = daeg[0, 0, 2, 2] + daeg[2, 0, 0, 2] 
 508      Rsuper[4, 3] = daeg[0, 1, 2, 2] + daeg[2, 1, 0, 2] 
 509   
 510      # Fifth column of the superoperator. 
 511      Rsuper[0, 4] = daeg[1, 0, 2, 0] + daeg[2, 0, 1, 0] 
 512      Rsuper[1, 4] = daeg[1, 1, 2, 1] + daeg[2, 1, 1, 1] 
 513      Rsuper[2, 4] = daeg[1, 0, 2, 1] + daeg[2, 0, 1, 1] 
 514      Rsuper[3, 4] = daeg[1, 0, 2, 2] + daeg[2, 0, 1, 2] 
 515      Rsuper[4, 4] = daeg[1, 1, 2, 2] + daeg[2, 1, 1, 2] 
 516   
 517      # Revert the shape. 
 518      daeg.shape = orig_shape 
 519   
 520      # Undo the T23 transpose. 
 521      transpose_23(daeg) 
 522   
 523   
 524 -def part_int_daeg1_pseudo_ellipse_xx(phi, x, y, smax): 
 525      """The theta-sigma partial integral of the 1st degree Frame Order matrix element xx for the pseudo-ellipse. 
 526   
 527      @param phi:     The azimuthal tilt-torsion angle. 
 528      @type phi:      float 
 529      @param x:       The cone opening angle along x. 
 530      @type x:        float 
 531      @param y:       The cone opening angle along y. 
 532      @type y:        float 
 533      @param smax:    The maximum torsion angle. 
 534      @type smax:     float 
 535      @return:        The theta-sigma partial integral. 
 536      @rtype:         float 
 537      """ 
 538   
 539      # Theta max. 
 540      tmax = tmax_pseudo_ellipse(phi, x, y) 
 541   
 542      # The theta-sigma integral. 
 543      return sin(smax) * (2 * (1 - cos(tmax)) * sin(phi)**2 + cos(phi)**2 * sin(tmax)**2) 
 544   
 545   
 546 -def part_int_daeg1_pseudo_ellipse_yy(phi, x, y, smax): 
 547      """The theta-sigma partial integral of the 1st degree Frame Order matrix element yy for the pseudo-ellipse. 
 548   
 549      @param phi:     The azimuthal tilt-torsion angle. 
 550      @type phi:      float 
 551      @param x:       The cone opening angle along x. 
 552      @type x:        float 
 553      @param y:       The cone opening angle along y. 
 554      @type y:        float 
 555      @param smax:    The maximum torsion angle. 
 556      @type smax:     float 
 557      @return:        The theta-sigma partial integral. 
 558      @rtype:         float 
 559      """ 
 560   
 561      # Theta max. 
 562      tmax = tmax_pseudo_ellipse(phi, x, y) 
 563   
 564      # The theta-sigma integral. 
 565      return sin(smax) * (2.0 * cos(phi)**2 * (1.0 - cos(tmax)) + sin(phi)**2 * sin(tmax)**2) 
 566   
 567   
 568 -def part_int_daeg1_pseudo_ellipse_zz(phi, x, y, smax): 
 569      """The theta-sigma partial integral of the 1st degree Frame Order matrix element zz for the pseudo-ellipse. 
 570   
 571      @param phi:     The azimuthal tilt-torsion angle. 
 572      @type phi:      float 
 573      @param x:       The cone opening angle along x. 
 574      @type x:        float 
 575      @param y:       The cone opening angle along y. 
 576      @type y:        float 
 577      @param smax:    The maximum torsion angle. 
 578      @type smax:     float 
 579      @return:        The theta-sigma partial integral. 
 580      @rtype:         float 
 581      """ 
 582   
 583      # Theta max. 
 584      tmax = tmax_pseudo_ellipse(phi, x, y) 
 585   
 586      # The theta-sigma integral. 
 587      return smax * sin(tmax)**2 
 588   
 589   
 590 -def part_int_daeg2_pseudo_ellipse_00(phi, x, y, smax): 
 591      """The theta-sigma partial integral of the 2nd degree Frame Order matrix element 11 for the pseudo-ellipse. 
 592   
 593      @param phi:     The azimuthal tilt-torsion angle. 
 594      @type phi:      float 
 595      @param x:       The cone opening angle along x. 
 596      @type x:        float 
 597      @param y:       The cone opening angle along y. 
 598      @type y:        float 
 599      @param smax:    The maximum torsion angle. 
 600      @type smax:     float 
 601      @return:        The theta-sigma partial integral. 
 602      @rtype:         float 
 603      """ 
 604   
 605      # Theta max. 
 606      tmax = tmax_pseudo_ellipse(phi, x, y) 
 607   
 608      # Repetitive trig. 
 609      cos_tmax = cos(tmax) 
 610      sin_tmax2 = sin(tmax)**2 
 611      cos_phi2 = cos(phi)**2 
 612   
 613      # Components. 
 614      a = sinc(2.0*smax/pi) * (2.0 * sin_tmax2 * cos_phi2 * ((2.0*cos_phi2 - 1.0)*cos(tmax)  -  6.0*(cos_phi2 - 1.0)) - 2.0*cos_tmax * (2.0*cos_phi2*(4.0*cos_phi2 - 5.0) + 3.0)) 
 615      b = 2.0*cos_phi2*cos_tmax*(sin_tmax2 + 2.0) - 6.0*cos_tmax 
 616   
 617      # Return the theta-sigma integral. 
 618      return a + b 
 619   
 620   
 621 -def part_int_daeg2_pseudo_ellipse_04(phi, x, y, smax): 
 622      """The theta-sigma partial integral of the 2nd degree Frame Order matrix element 22 for the pseudo-ellipse. 
 623   
 624      @param phi:     The azimuthal tilt-torsion angle. 
 625      @type phi:      float 
 626      @param x:       The cone opening angle along x. 
 627      @type x:        float 
 628      @param y:       The cone opening angle along y. 
 629      @type y:        float 
 630      @param smax:    The maximum torsion angle. 
 631      @type smax:     float 
 632      @return:        The theta-sigma partial integral. 
 633      @rtype:         float 
 634      """ 
 635   
 636      # Theta max. 
 637      tmax = tmax_pseudo_ellipse(phi, x, y) 
 638   
 639      # Repetitive trig. 
 640      cos_tmax = cos(tmax) 
 641      sin_tmax2 = sin(tmax)**2 
 642      cos_phi2 = cos(phi)**2 
 643      sin_phi2 = sin(phi)**2 
 644   
 645      # Components. 
 646      a = sinc(2.0*smax/pi) * (2.0 * sin_tmax2 * cos_phi2 * ((2.0*sin_phi2 - 1.0)*cos(tmax)  -  6.0*sin_phi2) + 2.0*cos_tmax * (2.0*cos_phi2*(4.0*cos_phi2 - 5.0) + 3.0)) 
 647      b = 2.0*cos_phi2*cos_tmax*(sin_tmax2 + 2.0) - 6.0*cos_tmax 
 648   
 649      # Return the theta-sigma integral. 
 650      return a + b 
 651   
 652   
 653 -def part_int_daeg2_pseudo_ellipse_08(phi, x, y, smax): 
 654      """The theta-sigma partial integral of the 2nd degree Frame Order matrix element 33 for the pseudo-ellipse. 
 655   
 656      @param phi:     The azimuthal tilt-torsion angle. 
 657      @type phi:      float 
 658      @param x:       The cone opening angle along x. 
 659      @type x:        float 
 660      @param y:       The cone opening angle along y. 
 661      @type y:        float 
 662      @param smax:    The maximum torsion angle. 
 663      @type smax:     float 
 664      @return:        The theta-sigma partial integral. 
 665      @rtype:         float 
 666      """ 
 667   
 668      # Theta max. 
 669      tmax = tmax_pseudo_ellipse(phi, x, y) 
 670   
 671      # The theta-sigma integral. 
 672      return cos(tmax) * cos(phi)**2 * (sin(tmax)**2 + 2.0) 
 673   
 674   
 675 -def part_int_daeg2_pseudo_ellipse_11(phi, x, y, smax): 
 676      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 677   
 678      @param phi:     The azimuthal tilt-torsion angle. 
 679      @type phi:      float 
 680      @param x:       The cone opening angle along x. 
 681      @type x:        float 
 682      @param y:       The cone opening angle along y. 
 683      @type y:        float 
 684      @param smax:    The maximum torsion angle. 
 685      @type smax:     float 
 686      @return:        The theta-sigma partial integral. 
 687      @rtype:         float 
 688      """ 
 689   
 690      # Theta max. 
 691      tmax = tmax_pseudo_ellipse(phi, x, y) 
 692   
 693      # Repetitive trig. 
 694      cos_tmax = cos(tmax) 
 695      cos_phi2 = cos(phi)**2 
 696      sin_phi2 = sin(phi)**2 
 697   
 698      # The integral. 
 699      a = sinc(2.0*smax/pi) * ((4.0*cos_phi2*((1.0 - cos_phi2)*cos_tmax + 3.0*(cos_phi2-1)) + 3.0)*sin(tmax)**2 - 16.0*cos_phi2*sin_phi2*cos_tmax) + 3.0*sin(tmax)**2 
 700   
 701      # The theta-sigma integral. 
 702      return a 
 703   
 704   
 705 -def part_int_daeg2_pseudo_ellipse_13(phi, x, y, smax): 
 706      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 707   
 708      @param phi:     The azimuthal tilt-torsion angle. 
 709      @type phi:      float 
 710      @param x:       The cone opening angle along x. 
 711      @type x:        float 
 712      @param y:       The cone opening angle along y. 
 713      @type y:        float 
 714      @param smax:    The maximum torsion angle. 
 715      @type smax:     float 
 716      @return:        The theta-sigma partial integral. 
 717      @rtype:         float 
 718      """ 
 719   
 720      # Theta max. 
 721      tmax = tmax_pseudo_ellipse(phi, x, y) 
 722   
 723      # Repetitive trig. 
 724      cos_tmax = cos(tmax) 
 725      sin_tmax2 = sin(tmax)**2 
 726      sinc_2smax = sinc(2.0*smax/pi) 
 727      cos_sin_phi2 = cos(phi)**2*sin(phi)**2 
 728   
 729      # The theta-sigma integral. 
 730      return sinc_2smax * (sin_tmax2 * (4*cos_sin_phi2*cos_tmax - 12*cos_sin_phi2 + 3) - 16*cos_sin_phi2*cos_tmax) - 3.0*sin_tmax2 
 731   
 732   
 733 -def part_int_daeg2_pseudo_ellipse_22(phi, x, y, smax): 
 734      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 735   
 736      @param phi:     The azimuthal tilt-torsion angle. 
 737      @type phi:      float 
 738      @param x:       The cone opening angle along x. 
 739      @type x:        float 
 740      @param y:       The cone opening angle along y. 
 741      @type y:        float 
 742      @param smax:    The maximum torsion angle. 
 743      @type smax:     float 
 744      @return:        The theta-sigma partial integral. 
 745      @rtype:         float 
 746      """ 
 747   
 748      # Theta max. 
 749      tmax = tmax_pseudo_ellipse(phi, x, y) 
 750   
 751      # Repetitive trig. 
 752      cos_tmax = cos(tmax) 
 753      cos_phi2 = cos(phi)**2 
 754   
 755      # Components. 
 756      a = 2.0*cos_phi2*cos_tmax**3 + 3.0*(1.0 - cos_phi2)*cos_tmax**2 
 757   
 758      # Return the theta-sigma integral. 
 759      return a 
 760   
 761   
 762 -def part_int_daeg2_pseudo_ellipse_26(phi, x, y, smax): 
 763      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 764   
 765      @param phi:     The azimuthal tilt-torsion angle. 
 766      @type phi:      float 
 767      @param x:       The cone opening angle along x. 
 768      @type x:        float 
 769      @param y:       The cone opening angle along y. 
 770      @type y:        float 
 771      @param smax:    The maximum torsion angle. 
 772      @type smax:     float 
 773      @return:        The theta-sigma partial integral. 
 774      @rtype:         float 
 775      """ 
 776   
 777      # Theta max. 
 778      tmax = tmax_pseudo_ellipse(phi, x, y) 
 779   
 780      # The theta-sigma integral. 
 781      return cos(phi)**2 * (cos(tmax)**3 - 3.0*cos(tmax)) 
 782   
 783   
 784 -def part_int_daeg2_pseudo_ellipse_40(phi, x, y, smax): 
 785      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 786   
 787      @param phi:     The azimuthal tilt-torsion angle. 
 788      @type phi:      float 
 789      @param x:       The cone opening angle along x. 
 790      @type x:        float 
 791      @param y:       The cone opening angle along y. 
 792      @type y:        float 
 793      @param smax:    The maximum torsion angle. 
 794      @type smax:     float 
 795      @return:        The theta-sigma partial integral. 
 796      @rtype:         float 
 797      """ 
 798   
 799      # Theta max. 
 800      tmax = tmax_pseudo_ellipse(phi, x, y) 
 801   
 802      # Repetitive trig. 
 803      cos_tmax = cos(tmax) 
 804      sin_tmax2 = sin(tmax)**2 
 805      cos_phi2 = cos(phi)**2 
 806      sin_phi2 = sin(phi)**2 
 807   
 808      # Components. 
 809      a = sinc(2.0*smax/pi) * (2.0 * sin_tmax2 * sin_phi2 * ((2.0*cos_phi2 - 1.0)*cos(tmax) - 6.0*cos_phi2) + 2.0*cos_tmax * (2.0*sin_phi2*(4.0*sin_phi2 - 5.0) + 3.0)) 
 810      b = 2.0*sin_phi2*cos_tmax*(sin_tmax2 + 2.0) - 6.0*cos_tmax 
 811   
 812      # Return the theta-sigma integral. 
 813      return a + b 
 814   
 815   
 816 -def part_int_daeg2_pseudo_ellipse_44(phi, x, y, smax): 
 817      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 818   
 819      @param phi:     The azimuthal tilt-torsion angle. 
 820      @type phi:      float 
 821      @param x:       The cone opening angle along x. 
 822      @type x:        float 
 823      @param y:       The cone opening angle along y. 
 824      @type y:        float 
 825      @param smax:    The maximum torsion angle. 
 826      @type smax:     float 
 827      @return:        The theta-sigma partial integral. 
 828      @rtype:         float 
 829      """ 
 830   
 831      # Theta max. 
 832      tmax = tmax_pseudo_ellipse(phi, x, y) 
 833   
 834      # Repetitive trig. 
 835      cos_tmax = cos(tmax) 
 836      sin_tmax2 = sin(tmax)**2 
 837      cos_phi2 = cos(phi)**2 
 838      sin_phi2 = sin(phi)**2 
 839   
 840      # Components. 
 841      a = sinc(2.0*smax/pi) * (2.0 * sin_tmax2 * sin_phi2 * ((2.0*sin_phi2 - 1.0)*cos(tmax)  +  6.0*cos_phi2) - 2.0*cos_tmax * (2.0*sin_phi2*(4.0*sin_phi2 - 5.0) + 3.0)) 
 842      b = 2.0*sin_phi2*cos_tmax*(sin_tmax2 + 2.0) - 6.0*cos_tmax 
 843   
 844      # Return the theta-sigma integral. 
 845      return a + b 
 846   
 847   
 848 -def part_int_daeg2_pseudo_ellipse_48(phi, x, y, smax): 
 849      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 850   
 851      @param phi:     The azimuthal tilt-torsion angle. 
 852      @type phi:      float 
 853      @param x:       The cone opening angle along x. 
 854      @type x:        float 
 855      @param y:       The cone opening angle along y. 
 856      @type y:        float 
 857      @param smax:    The maximum torsion angle. 
 858      @type smax:     float 
 859      @return:        The theta-sigma partial integral. 
 860      @rtype:         float 
 861      """ 
 862   
 863      # Theta max. 
 864      tmax = tmax_pseudo_ellipse(phi, x, y) 
 865   
 866      # The theta-sigma integral. 
 867      return cos(tmax) * sin(phi)**2 * (sin(tmax)**2 + 2.0) 
 868   
 869   
 870 -def part_int_daeg2_pseudo_ellipse_55(phi, x, y, smax): 
 871      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 872   
 873      @param phi:     The azimuthal tilt-torsion angle. 
 874      @type phi:      float 
 875      @param x:       The cone opening angle along x. 
 876      @type x:        float 
 877      @param y:       The cone opening angle along y. 
 878      @type y:        float 
 879      @param smax:    The maximum torsion angle. 
 880      @type smax:     float 
 881      @return:        The theta-sigma partial integral. 
 882      @rtype:         float 
 883      """ 
 884   
 885      # Theta max. 
 886      tmax = tmax_pseudo_ellipse(phi, x, y) 
 887   
 888      # Repetitive trig. 
 889      cos_tmax = cos(tmax) 
 890      sin_phi2 = sin(phi)**2 
 891   
 892      # Return the theta-sigma integral. 
 893      return 2.0*sin_phi2*cos_tmax**3 + 3.0*(1.0 - sin_phi2)*cos_tmax**2 
 894   
 895   
 896 -def part_int_daeg2_pseudo_ellipse_57(phi, x, y, smax): 
 897      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 898   
 899      @param phi:     The azimuthal tilt-torsion angle. 
 900      @type phi:      float 
 901      @param x:       The cone opening angle along x. 
 902      @type x:        float 
 903      @param y:       The cone opening angle along y. 
 904      @type y:        float 
 905      @param smax:    The maximum torsion angle. 
 906      @type smax:     float 
 907      @return:        The theta-sigma partial integral. 
 908      @rtype:         float 
 909      """ 
 910   
 911      # Theta max. 
 912      tmax = tmax_pseudo_ellipse(phi, x, y) 
 913   
 914      # The theta-sigma integral. 
 915      return sin(phi)**2 * (cos(tmax)**3 - 3.0*cos(tmax)) 
 916   
 917   
 918 -def part_int_daeg2_pseudo_ellipse_80(phi, x, y, smax): 
 919      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 920   
 921      @param phi:     The azimuthal tilt-torsion angle. 
 922      @type phi:      float 
 923      @param x:       The cone opening angle along x. 
 924      @type x:        float 
 925      @param y:       The cone opening angle along y. 
 926      @type y:        float 
 927      @param smax:    The maximum torsion angle. 
 928      @type smax:     float 
 929      @return:        The theta-sigma partial integral. 
 930      @rtype:         float 
 931      """ 
 932   
 933      # Theta max. 
 934      tmax = tmax_pseudo_ellipse(phi, x, y) 
 935   
 936      # Repetitive trig. 
 937      cos_tmax = cos(tmax) 
 938      sin_tmax2 = sin(tmax)**2 
 939      cos_phi2 = cos(phi)**2 
 940   
 941      # The theta-sigma integral. 
 942      return sinc(2.0*smax/pi) * (2.0*(1.0 - 2.0*cos_phi2)*cos_tmax*(sin_tmax2 + 2.0)) + 2.0*cos_tmax**3 - 6.0*cos_tmax 
 943   
 944   
 945 -def part_int_daeg2_pseudo_ellipse_84(phi, x, y, smax): 
 946      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 947   
 948      @param phi:     The azimuthal tilt-torsion angle. 
 949      @type phi:      float 
 950      @param x:       The cone opening angle along x. 
 951      @type x:        float 
 952      @param y:       The cone opening angle along y. 
 953      @type y:        float 
 954      @param smax:    The maximum torsion angle. 
 955      @type smax:     float 
 956      @return:        The theta-sigma partial integral. 
 957      @rtype:         float 
 958      """ 
 959   
 960      # Theta max. 
 961      tmax = tmax_pseudo_ellipse(phi, x, y) 
 962   
 963      # Repetitive trig. 
 964      cos_tmax = cos(tmax) 
 965      sin_tmax2 = sin(tmax)**2 
 966      cos_phi2 = cos(phi)**2 
 967   
 968      # The theta-sigma integral. 
 969      return sinc(2.0*smax/pi) * (2.0*(1.0 - 2.0*cos_phi2)*cos_tmax*(sin_tmax2 + 2.0)) - 2.0*cos_tmax**3 + 6.0*cos_tmax 
 970   
 971   
 972 -def part_int_daeg2_pseudo_ellipse_88(phi, x, y, smax): 
 973      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the pseudo-ellipse. 
 974   
 975      @param phi:     The azimuthal tilt-torsion angle. 
 976      @type phi:      float 
 977      @param x:       The cone opening angle along x. 
 978      @type x:        float 
 979      @param y:       The cone opening angle along y. 
 980      @type y:        float 
 981      @param smax:    The maximum torsion angle. 
 982      @type smax:     float 
 983      @return:        The theta-sigma partial integral. 
 984      @rtype:         float 
 985      """ 
 986   
 987      # Theta max. 
 988      tmax = tmax_pseudo_ellipse(phi, x, y) 
 989   
 990      # The theta-sigma integral. 
 991      return cos(tmax)**3 
 992   
 993   
 994 -def part_int_daeg2_pseudo_ellipse_free_rotor_00(phi, x, y): 
 995      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
 996   
 997      @param phi:     The azimuthal tilt-torsion angle. 
 998      @type phi:      float 
 999      @param x:       The cone opening angle along x. 
1000      @type x:        float 
1001      @param y:       The cone opening angle along y. 
1002      @type y:        float 
1003      @return:        The theta-sigma partial integral. 
1004      @rtype:         float 
1005      """ 
1006   
1007      # Theta max. 
1008      tmax = tmax_pseudo_ellipse(phi, x, y) 
1009   
1010      # The theta-sigma integral. 
1011      return cos(phi)**2*cos(tmax)**3 + 3.0*sin(phi)**2*cos(tmax) 
1012   
1013   
1014 -def part_int_daeg2_pseudo_ellipse_free_rotor_08(phi, x, y): 
1015      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
1016   
1017      @param phi:     The azimuthal tilt-torsion angle. 
1018      @type phi:      float 
1019      @param x:       The cone opening angle along x. 
1020      @type x:        float 
1021      @param y:       The cone opening angle along y. 
1022      @type y:        float 
1023      @return:        The theta-sigma partial integral. 
1024      @rtype:         float 
1025      """ 
1026   
1027      # Theta max. 
1028      tmax = tmax_pseudo_ellipse(phi, x, y) 
1029   
1030      # The theta-sigma integral. 
1031      return cos(phi)**2*(2.0*cos(tmax)**3 - 6.0*cos(tmax)) 
1032   
1033   
1034 -def part_int_daeg2_pseudo_ellipse_free_rotor_11(phi, x, y): 
1035      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
1036   
1037      @param phi:     The azimuthal tilt-torsion angle. 
1038      @type phi:      float 
1039      @param x:       The cone opening angle along x. 
1040      @type x:        float 
1041      @param y:       The cone opening angle along y. 
1042      @type y:        float 
1043      @return:        The theta-sigma partial integral. 
1044      @rtype:         float 
1045      """ 
1046   
1047      # Theta max. 
1048      tmax = tmax_pseudo_ellipse(phi, x, y) 
1049   
1050      # The theta-sigma integral. 
1051      return sin(tmax)**2 
1052   
1053   
1054 -def part_int_daeg2_pseudo_ellipse_free_rotor_44(phi, x, y): 
1055      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
1056   
1057      @param phi:     The azimuthal tilt-torsion angle. 
1058      @type phi:      float 
1059      @param x:       The cone opening angle along x. 
1060      @type x:        float 
1061      @param y:       The cone opening angle along y. 
1062      @type y:        float 
1063      @return:        The theta-sigma partial integral. 
1064      @rtype:         float 
1065      """ 
1066   
1067      # Theta max. 
1068      tmax = tmax_pseudo_ellipse(phi, x, y) 
1069   
1070      # The theta-sigma integral. 
1071      return sin(phi)**2*cos(tmax)**3 + 3*cos(phi)**2*cos(tmax) 
1072   
1073   
1074 -def part_int_daeg2_pseudo_ellipse_free_rotor_48(phi, x, y): 
1075      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
1076   
1077      @param phi:     The azimuthal tilt-torsion angle. 
1078      @type phi:      float 
1079      @param x:       The cone opening angle along x. 
1080      @type x:        float 
1081      @param y:       The cone opening angle along y. 
1082      @type y:        float 
1083      @return:        The theta-sigma partial integral. 
1084      @rtype:         float 
1085      """ 
1086   
1087      # Theta max. 
1088      tmax = tmax_pseudo_ellipse(phi, x, y) 
1089   
1090      # The theta-sigma integral. 
1091      return sin(phi)**2*(2.0*cos(tmax)**3 - 6.0*cos(tmax)) 
1092   
1093   
1094 -def part_int_daeg2_pseudo_ellipse_free_rotor_80(phi, x, y): 
1095      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
1096   
1097      @param phi:     The azimuthal tilt-torsion angle. 
1098      @type phi:      float 
1099      @param x:       The cone opening angle along x. 
1100      @type x:        float 
1101      @param y:       The cone opening angle along y. 
1102      @type y:        float 
1103      @return:        The theta-sigma partial integral. 
1104      @rtype:         float 
1105      """ 
1106   
1107      # Theta max. 
1108      tmax = tmax_pseudo_ellipse(phi, x, y) 
1109   
1110      # The theta-sigma integral. 
1111      return cos(tmax)**3 - 3.0*cos(tmax) 
1112   
1113   
1114 -def part_int_daeg2_pseudo_ellipse_free_rotor_88(phi, x, y): 
1115      """The theta-sigma partial integral of the 2nd degree Frame Order matrix for the free rotor pseudo-ellipse. 
1116   
1117      @param phi:     The azimuthal tilt-torsion angle. 
1118      @type phi:      float 
1119      @param x:       The cone opening angle along x. 
1120      @type x:        float 
1121      @param y:       The cone opening angle along y. 
1122      @type y:        float 
1123      @return:        The theta-sigma partial integral. 
1124      @rtype:         float 
1125      """ 
1126   
1127      # Theta max. 
1128      tmax = tmax_pseudo_ellipse(phi, x, y) 
1129   
1130      # The theta-sigma integral. 
1131      return cos(tmax)**3 
1132   
1133   
1134 -def part_int_daeg2_pseudo_ellipse_torsionless_00(phi, x, y): 
1135      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1136   
1137      @param phi:     The azimuthal tilt-torsion angle. 
1138      @type phi:      float 
1139      @param x:       The cone opening angle along x. 
1140      @type x:        float 
1141      @param y:       The cone opening angle along y. 
1142      @type y:        float 
1143      @return:        The theta partial integral. 
1144      @rtype:         float 
1145      """ 
1146   
1147      # Theta max. 
1148      tmax = tmax_pseudo_ellipse(phi, x, y) 
1149   
1150      # The theta integral. 
1151      return (2*cos(phi)**4*cos(tmax) + 6*cos(phi)**2*sin(phi)**2)*sin(tmax)**2 - (6*sin(phi)**4 + 2*cos(phi)**4)*cos(tmax) 
1152   
1153   
1154 -def part_int_daeg2_pseudo_ellipse_torsionless_04(phi, x, y): 
1155      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1156   
1157      @param phi:     The azimuthal tilt-torsion angle. 
1158      @type phi:      float 
1159      @param x:       The cone opening angle along x. 
1160      @type x:        float 
1161      @param y:       The cone opening angle along y. 
1162      @type y:        float 
1163      @return:        The theta partial integral. 
1164      @rtype:         float 
1165      """ 
1166   
1167      # Theta max. 
1168      tmax = tmax_pseudo_ellipse(phi, x, y) 
1169   
1170      # The theta integral. 
1171      return (2*cos(phi)**2*sin(phi)**2*cos(tmax) - 6*cos(phi)**2*sin(phi)**2)*sin(tmax)**2 - 8*cos(phi)**2*sin(phi)**2*cos(tmax) 
1172   
1173   
1174 -def part_int_daeg2_pseudo_ellipse_torsionless_08(phi, x, y): 
1175      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1176   
1177      @param phi:     The azimuthal tilt-torsion angle. 
1178      @type phi:      float 
1179      @param x:       The cone opening angle along x. 
1180      @type x:        float 
1181      @param y:       The cone opening angle along y. 
1182      @type y:        float 
1183      @return:        The theta partial integral. 
1184      @rtype:         float 
1185      """ 
1186   
1187      # Theta max. 
1188      tmax = tmax_pseudo_ellipse(phi, x, y) 
1189   
1190      # The theta integral. 
1191      return 2*cos(phi)**2*cos(tmax)**3 - 6*cos(phi)**2*cos(tmax) 
1192   
1193   
1194 -def part_int_daeg2_pseudo_ellipse_torsionless_11(phi, x, y): 
1195      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1196   
1197      @param phi:     The azimuthal tilt-torsion angle. 
1198      @type phi:      float 
1199      @param x:       The cone opening angle along x. 
1200      @type x:        float 
1201      @param y:       The cone opening angle along y. 
1202      @type y:        float 
1203      @return:        The theta partial integral. 
1204      @rtype:         float 
1205      """ 
1206   
1207      # Theta max. 
1208      tmax = tmax_pseudo_ellipse(phi, x, y) 
1209   
1210      # The theta integral. 
1211      return (2*cos(phi)**2*sin(phi)**2*cos(tmax) + 3*sin(phi)**4 + 3*cos(phi)**4)*sin(tmax)**2 - 8*cos(phi)**2*sin(phi)**2*cos(tmax) 
1212   
1213   
1214 -def part_int_daeg2_pseudo_ellipse_torsionless_22(phi, x, y): 
1215      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1216   
1217      @param phi:     The azimuthal tilt-torsion angle. 
1218      @type phi:      float 
1219      @param x:       The cone opening angle along x. 
1220      @type x:        float 
1221      @param y:       The cone opening angle along y. 
1222      @type y:        float 
1223      @return:        The theta partial integral. 
1224      @rtype:         float 
1225      """ 
1226   
1227      # Theta max. 
1228      tmax = tmax_pseudo_ellipse(phi, x, y) 
1229   
1230      # The theta integral. 
1231      return (2*sin(phi)**2 - 2)*cos(tmax)**3 - 3*sin(phi)**2*cos(tmax)**2 
1232   
1233   
1234 -def part_int_daeg2_pseudo_ellipse_torsionless_44(phi, x, y): 
1235      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1236   
1237      @param phi:     The azimuthal tilt-torsion angle. 
1238      @type phi:      float 
1239      @param x:       The cone opening angle along x. 
1240      @type x:        float 
1241      @param y:       The cone opening angle along y. 
1242      @type y:        float 
1243      @return:        The theta partial integral. 
1244      @rtype:         float 
1245      """ 
1246   
1247      # Theta max. 
1248      tmax = tmax_pseudo_ellipse(phi, x, y) 
1249   
1250      # The theta integral. 
1251      return (2*sin(phi)**4*cos(tmax) + 6*cos(phi)**2*sin(phi)**2)*sin(tmax)**2 - (2*sin(phi)**4 + 6*cos(phi)**4)*cos(tmax) 
1252   
1253   
1254 -def part_int_daeg2_pseudo_ellipse_torsionless_48(phi, x, y): 
1255      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1256   
1257      @param phi:     The azimuthal tilt-torsion angle. 
1258      @type phi:      float 
1259      @param x:       The cone opening angle along x. 
1260      @type x:        float 
1261      @param y:       The cone opening angle along y. 
1262      @type y:        float 
1263      @return:        The theta partial integral. 
1264      @rtype:         float 
1265      """ 
1266   
1267      # Theta max. 
1268      tmax = tmax_pseudo_ellipse(phi, x, y) 
1269   
1270      # The theta integral. 
1271      return 2*sin(phi)**2*cos(tmax)**3 - 6*sin(phi)**2*cos(tmax) 
1272   
1273   
1274 -def part_int_daeg2_pseudo_ellipse_torsionless_55(phi, x, y): 
1275      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1276   
1277      @param phi:     The azimuthal tilt-torsion angle. 
1278      @type phi:      float 
1279      @param x:       The cone opening angle along x. 
1280      @type x:        float 
1281      @param y:       The cone opening angle along y. 
1282      @type y:        float 
1283      @return:        The theta partial integral. 
1284      @rtype:         float 
1285      """ 
1286   
1287      # Theta max. 
1288      tmax = tmax_pseudo_ellipse(phi, x, y) 
1289   
1290      # The theta integral. 
1291      return (2*cos(phi)**2 - 2)*cos(tmax)**3 - 3*cos(phi)**2*cos(tmax)**2 
1292   
1293   
1294 -def part_int_daeg2_pseudo_ellipse_torsionless_88(phi, x, y): 
1295      """The theta partial integral of the 2nd degree Frame Order matrix for the torsionless pseudo-ellipse. 
1296   
1297      @param phi:     The azimuthal tilt-torsion angle. 
1298      @type phi:      float 
1299      @param x:       The cone opening angle along x. 
1300      @type x:        float 
1301      @param y:       The cone opening angle along y. 
1302      @type y:        float 
1303      @return:        The theta partial integral. 
1304      @rtype:         float 
1305      """ 
1306   
1307      # Theta max. 
1308      tmax = tmax_pseudo_ellipse(phi, x, y) 
1309   
1310      # The theta integral. 
1311      return 2 - 2*cos(tmax)**3 
1312   
1313   
1314 -def populate_1st_eigenframe_iso_cone(matrix, angle): 
1315      """Populate the 1st degree Frame Order matrix in the eigenframe for an isotropic cone. 
1316   
1317      The cone axis is assumed to be parallel to the z-axis in the eigenframe. 
1318   
1319      @param matrix:  The Frame Order matrix, 1st degree. 
1320      @type matrix:   numpy 3D, rank-2 array 
1321      @param angle:   The cone angle. 
1322      @type angle:    float 
1323      """ 
1324   
1325      # Zeros. 
1326      for i in range(3): 
1327          for j in range(3): 
1328              matrix[i, j] = 0.0 
1329   
1330      # The c33 element. 
1331      matrix[2, 2] = (cos(angle) + 1.0) / 2.0 
1332   
1333   
1334 -def populate_2nd_eigenframe_iso_cone(matrix, tmax, smax): 
1335      """Populate the 2nd degree Frame Order matrix in the eigenframe for the isotropic cone. 
1336   
1337      The cone axis is assumed to be parallel to the z-axis in the eigenframe. 
1338   
1339   
1340      @param matrix:  The Frame Order matrix, 2nd degree. 
1341      @type matrix:   numpy 9D, rank-2 array 
1342      @param tmax:    The cone opening angle. 
1343      @type tmax:     float 
1344      @param smax:    The maximum torsion angle. 
1345      @type smax:     float 
1346      """ 
1347   
1348      # Zeros. 
1349      for i in range(9): 
1350          for j in range(9): 
1351              matrix[i, j] = 0.0 
1352   
1353      # Repetitive trig calculations. 
1354      sinc_smax = sinc(smax/pi) 
1355      sinc_2smax = sinc(2.0*smax/pi) 
1356      cos_tmax = cos(tmax) 
1357      cos_tmax2 = cos_tmax**2 
1358   
1359      # Larger factors. 
1360      fact_sinc_2smax = sinc_2smax*(cos_tmax2 + 4.0*cos_tmax + 7.0) / 24.0 
1361      fact_cos_tmax2 = (cos_tmax2 + cos_tmax + 4.0) / 12.0 
1362      fact_cos_tmax = (cos_tmax + 1.0) / 4.0 
1363   
1364      # Diagonal. 
1365      matrix[0, 0] = fact_sinc_2smax  +  fact_cos_tmax2 
1366      matrix[1, 1] = fact_sinc_2smax  +  fact_cos_tmax 
1367      matrix[2, 2] = sinc_smax * (2.0*cos_tmax2 + 5.0*cos_tmax + 5.0) / 12.0 
1368      matrix[3, 3] = matrix[1, 1] 
1369      matrix[4, 4] = matrix[0, 0] 
1370      matrix[5, 5] = matrix[2, 2] 
1371      matrix[6, 6] = matrix[2, 2] 
1372      matrix[7, 7] = matrix[2, 2] 
1373      matrix[8, 8] = (cos_tmax2 + cos_tmax + 1.0) / 3.0 
1374   
1375      # Off diagonal set 1. 
1376      matrix[0, 4] = matrix[4, 0] = -fact_sinc_2smax  +  fact_cos_tmax2 
1377      matrix[0, 8] = matrix[8, 0] = -(cos_tmax2 + cos_tmax - 2.0) / 6.0 
1378      matrix[4, 8] = matrix[8, 4] = matrix[0, 8] 
1379   
1380      # Off diagonal set 2. 
1381      matrix[1, 3] = matrix[3, 1] = fact_sinc_2smax  -  fact_cos_tmax 
1382      matrix[2, 6] = matrix[6, 2] = sinc_smax * (cos_tmax2 + cos_tmax - 2.0) / 6.0 
1383      matrix[5, 7] = matrix[7, 5] = matrix[2, 6] 
1384   
1385   
1386 -def populate_2nd_eigenframe_iso_cone_free_rotor(matrix, s1): 
1387      """Populate the 2nd degree Frame Order matrix in the eigenframe for the free rotor isotropic cone. 
1388   
1389      The cone axis is assumed to be parallel to the z-axis in the eigenframe.  In this model, the three order parameters are defined as:: 
1390   
1391          S1 = S2, 
1392          S3 = 0 
1393   
1394      This is in the Kronecker product form. 
1395   
1396   
1397      @param matrix:  The Frame Order matrix, 2nd degree. 
1398      @type matrix:   numpy 9D, rank-2 array 
1399      @param s1:      The cone order parameter. 
1400      @type s1:       float 
1401      """ 
1402   
1403      # Zeros. 
1404      for i in range(9): 
1405          for j in range(9): 
1406              matrix[i, j] = 0.0 
1407   
1408      # The c11^2, c22^2, c12^2, and c21^2 elements. 
1409      matrix[0, 0] = matrix[4, 4] = (s1 + 2.0) / 6.0 
1410      matrix[0, 4] = matrix[4, 0] = matrix[0, 0] 
1411   
1412      # The c33^2 element. 
1413      matrix[8, 8] = (2.0*s1 + 1.0) / 3.0 
1414   
1415      # The c13^2, c31^2, c23^2, c32^2 elements. 
1416      matrix[0, 8] = matrix[8, 0] = (1.0 - s1) / 3.0 
1417      matrix[4, 8] = matrix[8, 4] = matrix[0, 8] 
1418   
1419      # Calculate the cone angle. 
1420      theta = order_parameters.iso_cone_S_to_theta(s1) 
1421   
1422      # The c11.c22 and c12.c21 elements. 
1423      matrix[1, 1] = matrix[3, 3] = (cos(theta) + 1.0) / 4.0 
1424      matrix[1, 3] = matrix[3, 1] = -matrix[1, 1] 
1425   
1426   
1427 -def reduce_alignment_tensor(D, A, red_tensor): 
1428      """Calculate the reduction in the alignment tensor caused by the Frame Order matrix. 
1429   
1430      This is both the forward rotation notation and Kronecker product arrangement. 
1431   
1432      @param D:           The Frame Order matrix, 2nd degree to be populated. 
1433      @type D:            numpy 9D, rank-2 array 
1434      @param A:           The full alignment tensor in {Axx, Ayy, Axy, Axz, Ayz} notation. 
1435      @type A:            numpy 5D, rank-1 array 
1436      @param red_tensor:  The structure in {Axx, Ayy, Axy, Axz, Ayz} notation to place the reduced 
1437                          alignment tensor. 
1438      @type red_tensor:   numpy 5D, rank-1 array 
1439      """ 
1440   
1441      # The reduced tensor element A0. 
1442      red_tensor[0] =                 (D[0, 0] - D[0, 8])*A[0] 
1443      red_tensor[0] = red_tensor[0] + (D[0, 4] - D[0, 8])*A[1] 
1444      red_tensor[0] = red_tensor[0] + (D[0, 1] + D[0, 3])*A[2] 
1445      red_tensor[0] = red_tensor[0] + (D[0, 2] + D[0, 6])*A[3] 
1446      red_tensor[0] = red_tensor[0] + (D[0, 5] + D[0, 7])*A[4] 
1447                                                        
1448      # The reduced tensor element A1.                  
1449      red_tensor[1] =                 (D[4, 0] - D[4, 8])*A[0] 
1450      red_tensor[1] = red_tensor[1] + (D[4, 4] - D[4, 8])*A[1] 
1451      red_tensor[1] = red_tensor[1] + (D[4, 1] + D[4, 3])*A[2] 
1452      red_tensor[1] = red_tensor[1] + (D[4, 2] + D[4, 6])*A[3] 
1453      red_tensor[1] = red_tensor[1] + (D[4, 5] + D[4, 7])*A[4] 
1454                                                        
1455      # The reduced tensor element A2.                  
1456      red_tensor[2] =                 (D[1, 0] - D[1, 8])*A[0] 
1457      red_tensor[2] = red_tensor[2] + (D[1, 4] - D[1, 8])*A[1] 
1458      red_tensor[2] = red_tensor[2] + (D[1, 1] + D[1, 3])*A[2] 
1459      red_tensor[2] = red_tensor[2] + (D[1, 2] + D[1, 6])*A[3] 
1460      red_tensor[2] = red_tensor[2] + (D[1, 5] + D[1, 7])*A[4] 
1461                                                        
1462      # The reduced tensor element A3.                  
1463      red_tensor[3] =                 (D[2, 0] - D[2, 8])*A[0] 
1464      red_tensor[3] = red_tensor[3] + (D[2, 4] - D[2, 8])*A[1] 
1465      red_tensor[3] = red_tensor[3] + (D[2, 1] + D[2, 3])*A[2] 
1466      red_tensor[3] = red_tensor[3] + (D[2, 2] + D[2, 6])*A[3] 
1467      red_tensor[3] = red_tensor[3] + (D[2, 5] + D[2, 7])*A[4] 
1468                                                        
1469      # The reduced tensor element A4.                  
1470      red_tensor[4] =                 (D[5, 0] - D[5, 8])*A[0] 
1471      red_tensor[4] = red_tensor[4] + (D[5, 4] - D[5, 8])*A[1] 
1472      red_tensor[4] = red_tensor[4] + (D[5, 1] + D[5, 3])*A[2] 
1473      red_tensor[4] = red_tensor[4] + (D[5, 2] + D[5, 6])*A[3] 
1474      red_tensor[4] = red_tensor[4] + (D[5, 5] + D[5, 7])*A[4] 
1475   
1476   
1477 -def reduce_alignment_tensor_symmetric(D, A, red_tensor): 
1478      """Calculate the reduction in the alignment tensor caused by the Frame Order matrix. 
1479   
1480      This is both the forward rotation notation and Kronecker product arrangement.  This simplification is due to the symmetry in motion of the pseudo-elliptic and isotropic cones.  All element of the frame order matrix where an index appears only once are zero. 
1481   
1482      @param D:           The Frame Order matrix, 2nd degree to be populated. 
1483      @type D:            numpy 9D, rank-2 array 
1484      @param A:           The full alignment tensor in {Axx, Ayy, Axy, Axz, Ayz} notation. 
1485      @type A:            numpy 5D, rank-1 array 
1486      @param red_tensor:  The structure in {Axx, Ayy, Axy, Axz, Ayz} notation to place the reduced 
1487                          alignment tensor. 
1488      @type red_tensor:   numpy 5D, rank-1 array 
1489      """ 
1490   
1491      # The reduced tensor elements. 
1492      red_tensor[0] = (D[0, 0] - D[0, 8])*A[0]  +  (D[0, 4] - D[0, 8])*A[1] 
1493      red_tensor[1] = (D[4, 0] - D[4, 8])*A[0]  +  (D[4, 4] - D[4, 8])*A[1] 
1494      red_tensor[2] = (D[1, 1] + D[1, 3])*A[2] 
1495      red_tensor[3] = (D[2, 2] + D[2, 6])*A[3] 
1496      red_tensor[4] = (D[5, 5] + D[5, 7])*A[4] 
1497   
1498   
1499 -def rotate_daeg(matrix, R): 
1500      """Rotate the given frame order matrix. 
1501   
1502      It is assumed that the frame order matrix is in the Kronecker product form. 
1503   
1504   
1505      @param matrix:      The Frame Order matrix, 2nd degree to be populated. 
1506      @type matrix:       numpy 9D, rank-2 array 
1507      @param R:           The rotation matrix to be populated. 
1508      @type R:            numpy 3D, rank-2 array 
1509      """ 
1510   
1511      # The outer product of R. 
1512      R_kron = kron_prod(R, R) 
1513   
1514      # Rotate. 
1515      matrix_rot = dot(R_kron, dot(matrix, transpose(R_kron))) 
1516   
1517      # Return the matrix. 
1518      return matrix_rot 
1519   
1520   
1521 -def tmax_pseudo_ellipse(phi, theta_x, theta_y): 
1522      """The pseudo-ellipse tilt-torsion polar angle. 
1523   
1524      @param phi:     The azimuthal tilt-torsion angle. 
1525      @type phi:      float 
1526      @param theta_x: The cone opening angle along x. 
1527      @type theta_x:  float 
1528      @param theta_y: The cone opening angle along y. 
1529      @type theta_y:  float 
1530      @return:        The theta max angle for the given phi angle. 
1531      @rtype:         float 
1532      """ 
1533   
1534      # Zero points. 
1535      if theta_x == 0.0: 
1536          return 0.0 
1537      elif theta_y == 0.0: 
1538          return 0.0 
1539   
1540      # Return the maximum angle. 
1541      return theta_x * theta_y / sqrt((cos(phi)*theta_y)**2 + (sin(phi)*theta_x)**2) 
1542