Package target_functions :: Module frame_order :: Class Frame_order
[hide private]
[frames] | no frames]

Class Frame_order

source code

Class containing the target function of the optimisation of Frame Order matrix components.

Instance Methods [hide private]
 
__init__(self, model=None, init_params=None, full_tensors=None, full_in_ref_frame=None, rdcs=None, rdc_errors=None, rdc_weights=None, rdc_vect=None, dip_const=None, pcs=None, pcs_errors=None, pcs_weights=None, atomic_pos=None, temp=None, frq=None, paramag_centre=array([ 0., 0., 0.]), scaling_matrix=None, sobol_max_points=200, sobol_oversample=100, com=None, ave_pos_pivot=array([ 0., 0., 0.]), pivot=None, pivot_opt=False, quad_int=False)
Set up the target functions for the Frame Order theories.
source code
 
_init_tensors(self)
Set up isotropic cone optimisation against the alignment tensor data.
source code
float
func_double_rotor_qr_int(self, params)
Quasi-random Sobol' integration target function for the double rotor model.
source code
float
func_double_rotor_quad_int(self, params)
SciPy quadratic integration target function for the double rotor model.
source code
float
func_free_rotor_qr_int(self, params)
Quasi-random Sobol' integration target function for the free rotor model.
source code
float
func_free_rotor_quad_int(self, params)
SciPy quadratic integration target function for the free rotor model.
source code
float
func_iso_cone_qr_int(self, params)
Quasi-random Sobol' integration target function for the isotropic cone model.
source code
float
func_iso_cone_quad_int(self, params)
SciPy quadratic integration target function for the isotropic cone model.
source code
float
func_iso_cone_free_rotor_qr_int(self, params)
Quasi-random Sobol' integration target function for the free rotor isotropic cone model.
source code
float
func_iso_cone_free_rotor_quad_int(self, params)
SciPy quadratic integration target function for the free rotor isotropic cone model.
source code
float
func_iso_cone_torsionless_qr_int(self, params)
Quasi-random Sobol' integration target function for the torsionless isotropic cone model.
source code
float
func_iso_cone_torsionless_quad_int(self, params)
SciPy quadratic integration target function for the torsionless isotropic cone model.
source code
float
func_pseudo_ellipse_qr_int(self, params)
Quasi-random Sobol' integration target function for the pseudo-ellipse model.
source code
float
func_pseudo_ellipse_quad_int(self, params)
SciPy quadratic integration target function for the pseudo-ellipse model.
source code
float
func_pseudo_ellipse_free_rotor_qr_int(self, params)
Quasi-random Sobol' integration target function for the free-rotor pseudo-ellipse model.
source code
float
func_pseudo_ellipse_free_rotor_quad_int(self, params)
SciPy quadratic integration target function for the free-rotor pseudo-ellipse model.
source code
float
func_pseudo_ellipse_torsionless_qr_int(self, params)
Quasi-random Sobol' integration target function for the torsionless pseudo-ellipse model.
source code
float
func_pseudo_ellipse_torsionless_quad_int(self, params)
SciPy quadratic integration target function for the torsionless pseudo-ellipse model.
source code
float
func_rigid(self, params)
Target function for rigid model optimisation.
source code
float
func_rotor_qr_int(self, params)
Quasi-random Sobol' integration target function for the rotor model.
source code
float
func_rotor_quad_int(self, params)
SciPy quadratic integration target function for rotor model.
source code
 
calc_vectors(self, pivot=None, pivot2=None, R_ave=None, RT_ave=None)
Calculate the pivot to atom and lanthanide to pivot vectors for the target functions.
source code
 
create_sobol_data(self, dims=None)
Create the Sobol' quasi-random data for numerical integration.
source code
 
reduce_and_rot(self, ave_pos_alpha=None, ave_pos_beta=None, ave_pos_gamma=None, daeg=None)
Reduce and rotate the alignments tensors using the frame order matrix and Euler angles.
source code
Method Details [hide private]

__init__(self, model=None, init_params=None, full_tensors=None, full_in_ref_frame=None, rdcs=None, rdc_errors=None, rdc_weights=None, rdc_vect=None, dip_const=None, pcs=None, pcs_errors=None, pcs_weights=None, atomic_pos=None, temp=None, frq=None, paramag_centre=array([ 0., 0., 0.]), scaling_matrix=None, sobol_max_points=200, sobol_oversample=100, com=None, ave_pos_pivot=array([ 0., 0., 0.]), pivot=None, pivot_opt=False, quad_int=False)
(Constructor)

source code 

Set up the target functions for the Frame Order theories.

Parameters:
  • model (str) - The name of the Frame Order model.
  • init_params (numpy float64 array) - The initial parameter values.
  • full_tensors (numpy nx5D, rank-1 float64 array) - An array of the {Axx, Ayy, Axy, Axz, Ayz} values for all full alignment tensors. The format is [Axx1, Ayy1, Axy1, Axz1, Ayz1, Axx2, Ayy2, Axy2, Axz2, Ayz2, ..., Axxn, Ayyn, Axyn, Axzn, Ayzn].
  • full_in_ref_frame (numpy rank-1 array) - An array of flags specifying if the tensor in the reference frame is the full or reduced tensor.
  • rdcs (numpy rank-2 array) - The RDC lists. The first index must correspond to the different alignment media i and the second index to the spin systems j.
  • rdc_errors (numpy rank-2 array) - The RDC error lists. The dimensions of this argument are the same as for 'rdcs'.
  • rdc_weights (numpy rank-2 array) - The RDC weight lists. The dimensions of this argument are the same as for 'rdcs'.
  • rdc_vect (numpy rank-2 array) - The unit XH vector lists corresponding to the RDC values. The first index must correspond to the spin systems and the second index to the x, y, z elements.
  • dip_const (numpy rank-1 array) - The dipolar constants for each RDC. The indices correspond to the spin systems j.
  • pcs (numpy rank-2 array) - The PCS lists. The first index must correspond to the different alignment media i and the second index to the spin systems j.
  • pcs_errors (numpy rank-2 array) - The PCS error lists. The dimensions of this argument are the same as for 'pcs'.
  • pcs_weights (numpy rank-2 array) - The PCS weight lists. The dimensions of this argument are the same as for 'pcs'.
  • atomic_pos (numpy rank-3 array) - The atomic positions of all spins for the PCS and PRE data. The first index is the spin systems j and the second is the structure or state c.
  • temp (numpy rank-1 array) - The temperature of each PCS data set.
  • frq (numpy rank-1 array) - The frequency of each PCS data set.
  • paramag_centre (numpy rank-1, 3D array or rank-2, Nx3 array) - The paramagnetic centre position (or positions).
  • scaling_matrix (numpy rank-2 array) - The square and diagonal scaling matrix.
  • sobol_max_points (int) - The maximum number of Sobol' points to use for the numerical PCS integration technique.
  • sobol_oversample (int) - The oversampling factor Ov used for the total number of points N * Ov * 10**M, where N is the maximum number of Sobol' points and M is the number of dimensions or torsion-tilt angles for the system.
  • com (numpy 3D rank-1 array) - The centre of mass of the system. This is used for defining the rotor model systems.
  • ave_pos_pivot (numpy 3D rank-1 array) - The pivot point to rotate all atoms about to the average domain position. In most cases this will be the centre of mass of the moving domain. This pivot is shifted by the translation vector.
  • pivot (numpy rank-1, 3D array or None) - The pivot point for the ball-and-socket joint motion. This is needed if PCS or PRE values are used.
  • pivot_opt (bool) - A flag which if True will allow the pivot point of the motion to be optimised.
  • quad_int (bool) - A flag which if True will perform high precision numerical integration via the scipy.integrate quad(), dblquad() and tplquad() integration methods rather than the rough quasi-random numerical integration.

func_double_rotor_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the double rotor model.

This function optimises the model parameters using the RDC and PCS base data. Quasi-random, Sobol' sequence based, numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values. These can include {pivot_x, pivot_y, pivot_z, pivot_disp, ave_pos_x, ave_pos_y, ave_pos_z, ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_sigma_max, cone_sigma_max_2}.
Returns: float
The chi-squared or SSE value.

func_double_rotor_quad_int(self, params)

source code 

SciPy quadratic integration target function for the double rotor model.

This function optimises the model parameters using the RDC and PCS base data. Quasi-random, Sobol' sequence based, numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values. These can include {pivot_x, pivot_y, pivot_z, pivot_disp, ave_pos_x, ave_pos_y, ave_pos_z, ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_sigma_max, cone_sigma_max_2}.
Returns: float
The chi-squared or SSE value.

func_free_rotor_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the free rotor model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Simple numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma, theta, phi}.
Returns: float
The chi-squared or SSE value.

func_free_rotor_quad_int(self, params)

source code 

SciPy quadratic integration target function for the free rotor model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma, theta, phi}.
Returns: float
The chi-squared or SSE value.

func_iso_cone_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the isotropic cone model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Simple numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, theta, phi, cone_theta, sigma_max} where the first 3 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, cone_theta is the cone opening half angle, and sigma_max is the torsion angle.
Returns: float
The chi-squared or SSE value.

func_iso_cone_quad_int(self, params)

source code 

SciPy quadratic integration target function for the isotropic cone model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {beta, gamma, theta, phi, s1} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and s1 is the isotropic cone order parameter.
Returns: float
The chi-squared or SSE value.

func_iso_cone_free_rotor_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the free rotor isotropic cone model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Simple numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {beta, gamma, theta, phi, s1} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and s1 is the isotropic cone order parameter.
Returns: float
The chi-squared or SSE value.

func_iso_cone_free_rotor_quad_int(self, params)

source code 

SciPy quadratic integration target function for the free rotor isotropic cone model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {beta, gamma, theta, phi, s1} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and s1 is the isotropic cone order parameter.
Returns: float
The chi-squared or SSE value.

func_iso_cone_torsionless_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the torsionless isotropic cone model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Simple numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {beta, gamma, theta, phi, cone_theta} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and cone_theta is cone opening angle.
Returns: float
The chi-squared or SSE value.

func_iso_cone_torsionless_quad_int(self, params)

source code 

SciPy quadratic integration target function for the torsionless isotropic cone model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {beta, gamma, theta, phi, cone_theta} where the first 2 are the tensor rotation Euler angles, the next two are the polar and azimuthal angles of the cone axis, and cone_theta is cone opening angle.
Returns: float
The chi-squared or SSE value.

func_pseudo_ellipse_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the pseudo-ellipse model.

This function optimises the model parameters using the RDC and PCS base data. Quasi-random, Sobol' sequence based, numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 3 are the pseudo-elliptic cone geometric parameters.
Returns: float
The chi-squared or SSE value.

func_pseudo_ellipse_quad_int(self, params)

source code 

SciPy quadratic integration target function for the pseudo-ellipse model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 3 are the pseudo-elliptic cone geometric parameters.
Returns: float
The chi-squared or SSE value.

func_pseudo_ellipse_free_rotor_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the free-rotor pseudo-ellipse model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Simple numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 2 are the free_rotor pseudo-elliptic cone geometric parameters.
Returns: float
The chi-squared or SSE value.

func_pseudo_ellipse_free_rotor_quad_int(self, params)

source code 

SciPy quadratic integration target function for the free-rotor pseudo-ellipse model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 2 are the free_rotor pseudo-elliptic cone geometric parameters.
Returns: float
The chi-squared or SSE value.

func_pseudo_ellipse_torsionless_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the torsionless pseudo-ellipse model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Simple numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 2 are the torsionless pseudo-elliptic cone geometric parameters.
Returns: float
The chi-squared or SSE value.

func_pseudo_ellipse_torsionless_quad_int(self, params)

source code 

SciPy quadratic integration target function for the torsionless pseudo-ellipse model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values {alpha, beta, gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y} where the first 3 are the average position rotation Euler angles, the next 3 are the Euler angles defining the eigenframe, and the last 2 are the torsionless pseudo-elliptic cone geometric parameters.
Returns: float
The chi-squared or SSE value.

func_rigid(self, params)

source code 

Target function for rigid model optimisation.

This function optimises the isotropic cone model parameters using the RDC and PCS base data.

Parameters:
  • params (list of float) - The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma}.
Returns: float
The chi-squared or SSE value.

func_rotor_qr_int(self, params)

source code 

Quasi-random Sobol' integration target function for the rotor model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Quasi-random, Sobol' sequence based, numerical integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma, theta, phi, sigma_max}.
Returns: float
The chi-squared or SSE value.

func_rotor_quad_int(self, params)

source code 

SciPy quadratic integration target function for rotor model.

This function optimises the isotropic cone model parameters using the RDC and PCS base data. Scipy quadratic integration is used for the PCS.

Parameters:
  • params (list of float) - The vector of parameter values. These are the tensor rotation angles {alpha, beta, gamma, theta, phi, sigma_max}.
Returns: float
The chi-squared or SSE value.

calc_vectors(self, pivot=None, pivot2=None, R_ave=None, RT_ave=None)

source code 

Calculate the pivot to atom and lanthanide to pivot vectors for the target functions.

Parameters:
  • pivot (numpy rank-1, 3D array) - The pivot point.
  • pivot2 (numpy rank-1, 3D array) - The 2nd pivot point.
  • R_ave (numpy rank-2, 3D array) - The rotation matrix for rotating from the reference frame to the average position.
  • RT_ave (numpy rank-2, 3D array) - The transpose of R_ave.

create_sobol_data(self, dims=None)

source code 

Create the Sobol' quasi-random data for numerical integration.

This uses the external sobol_lib module to create the data. The algorithm is that modified by Antonov and Saleev.

Parameters:
  • dims (list of str) - The list of parameters.

reduce_and_rot(self, ave_pos_alpha=None, ave_pos_beta=None, ave_pos_gamma=None, daeg=None)

source code 

Reduce and rotate the alignments tensors using the frame order matrix and Euler angles.

Parameters:
  • ave_pos_alpha (float) - The alpha Euler angle describing the average domain position, the tensor rotation.
  • ave_pos_beta (float) - The beta Euler angle describing the average domain position, the tensor rotation.
  • ave_pos_gamma (float) - The gamma Euler angle describing the average domain position, the tensor rotation.
  • daeg (rank-2, 9D array) - The 2nd degree frame order matrix.