1   
  2   
  3   
  4   
  5   
  6   
  7   
  8   
  9   
 10   
 11   
 12   
 13   
 14   
 15   
 16   
 17   
 18   
 19   
 20   
 21   
 22   
 23   
 24  """Module containing the target functions of the Frame Order theories.""" 
 25   
 26   
 27  from copy import deepcopy 
 28  from numpy import array, dot, float64, ones, transpose, zeros 
 29   
 30   
 31  from generic_fns.frame_order import print_frame_order_2nd_degree 
 32  from maths_fns.alignment_tensor import to_5D, to_tensor 
 33  from maths_fns.chi2 import chi2 
 34  from maths_fns.frame_order_matrix_ops import compile_2nd_matrix_free_rotor, compile_2nd_matrix_iso_cone, compile_2nd_matrix_iso_cone_free_rotor, compile_2nd_matrix_iso_cone_torsionless, compile_2nd_matrix_pseudo_ellipse, compile_2nd_matrix_pseudo_ellipse_free_rotor, compile_2nd_matrix_pseudo_ellipse_torsionless, compile_2nd_matrix_rotor, reduce_alignment_tensor 
 35  from maths_fns.rotation_matrix import euler_to_R_zyz as euler_to_R 
 36  from relax_errors import RelaxError 
 37   
 38   
 40      """Class containing the target function of the optimisation of Frame Order matrix components.""" 
 41   
 42 -    def __init__(self, model=None, init_params=None, full_tensors=None, red_tensors=None, red_errors=None, full_in_ref_frame=None, frame_order_2nd=None): 
  43          """Set up the target functions for the Frame Order theories. 
 44           
 45          @keyword model:             The name of the Frame Order model. 
 46          @type model:                str 
 47          @keyword init_params:       The initial parameter values. 
 48          @type init_params:          numpy float64 array 
 49          @keyword full_tensors:      An array of the {Sxx, Syy, Sxy, Sxz, Syz} values for all full alignment tensors.  The format is [Sxx1, Syy1, Sxy1, Sxz1, Syz1, Sxx2, Syy2, Sxy2, Sxz2, Syz2, ..., Sxxn, Syyn, Sxyn, Sxzn, Syzn]. 
 50          @type full_tensors:         numpy nx5D, rank-1 float64 array 
 51          @keyword red_tensors:       An array of the {Sxx, Syy, Sxy, Sxz, Syz} values for all reduced tensors.  The array format is the same as for full_tensors. 
 52          @type red_tensors:          numpy nx5D, rank-1 float64 array 
 53          @keyword red_errors:        An array of the {Sxx, Syy, Sxy, Sxz, Syz} errors for all reduced tensors.  The array format is the same as for full_tensors. 
 54          @type red_errors:           numpy nx5D, rank-1 float64 array 
 55          @keyword full_in_ref_frame: An array of flags specifying if the tensor in the reference frame is the full or reduced tensor. 
 56          @type full_in_ref_frame:    numpy rank-1 array 
 57          @keyword frame_order_2nd:   The numerical values of the 2nd degree Frame Order matrix.  If supplied, the target functions will optimise directly to these values. 
 58          @type frame_order_2nd:      None or numpy 9D, rank-2 array 
 59          """ 
 60   
 61           
 62          if not model: 
 63              raise RelaxError("The type of Frame Order model must be specified.") 
 64   
 65           
 66          self.params = deepcopy(init_params) 
 67   
 68           
 69          self.model = model 
 70          self.full_tensors = full_tensors 
 71          self.red_tensors = red_tensors 
 72          self.red_errors = red_errors 
 73          self.full_in_ref_frame = full_in_ref_frame 
 74   
 75           
 76          if full_tensors != None and frame_order_2nd != None: 
 77              raise RelaxError("Tensors and Frame Order matrices cannot be supplied together.") 
 78          if full_tensors == None and frame_order_2nd == None: 
 79              raise RelaxError("The arguments have been incorrectly supplied, cannot determine the optimisation mode.") 
 80   
 81           
 82          if full_tensors != None: 
 83              self.__init_tensors() 
 84   
 85           
 86          if model == 'iso cone, free rotor' and frame_order_2nd != None: 
 87              self.__init_iso_cone_elements(frame_order_2nd) 
 88   
 89           
 90          if model == 'pseudo-ellipse': 
 91              self.func = self.func_pseudo_ellipse 
 92          elif model == 'pseudo-ellipse, torsionless': 
 93              self.func = self.func_pseudo_ellipse_torsionless 
 94          elif model == 'pseudo-ellipse, free rotor': 
 95              self.func = self.func_pseudo_ellipse_free_rotor 
 96          elif model == 'iso cone': 
 97              self.func = self.func_iso_cone 
 98          elif model == 'iso cone, torsionless': 
 99              self.func = self.func_iso_cone_torsionless 
100          elif model == 'iso cone, free rotor': 
101              self.func = self.func_iso_cone_free_rotor 
102          elif model == 'line': 
103              self.func = self.func_line 
104          elif model == 'line, torsionless': 
105              self.func = self.func_line_torsionless 
106          elif model == 'line, free rotor': 
107              self.func = self.func_line_free_rotor 
108          elif model == 'rotor': 
109              self.func = self.func_rotor 
110          elif model == 'rigid': 
111              self.func = self.func_rigid 
112          elif model == 'free rotor': 
113              self.func = self.func_free_rotor 
 114   
115   
117          """Set up isotropic cone optimisation against the alignment tensor data.""" 
118   
119           
120          if self.full_tensors == None or not len(self.full_tensors): 
121              raise RelaxError("The full_tensors argument " + repr(self.full_tensors) + " must be supplied.") 
122          if self.red_tensors == None or not len(self.red_tensors): 
123              raise RelaxError("The red_tensors argument " + repr(self.red_tensors) + " must be supplied.") 
124          if self.red_errors == None or not len(self.red_errors): 
125              raise RelaxError("The red_errors argument " + repr(self.red_errors) + " must be supplied.") 
126          if self.full_in_ref_frame == None or not len(self.full_in_ref_frame): 
127              raise RelaxError("The full_in_ref_frame argument " + repr(self.full_in_ref_frame) + " must be supplied.") 
128   
129           
130          self.num_tensors = len(self.full_tensors) / 5 
131          self.red_tensors_bc = zeros(self.num_tensors*5, float64) 
132   
133           
134          self.rot = zeros((3, 3), float64) 
135          self.tensor_3D = zeros((3, 3), float64) 
136   
137           
138          self.cone_axis = zeros(3, float64) 
139          self.z_axis = array([0, 0, 1], float64) 
140   
141           
142          self.frame_order_2nd = zeros((9, 9), float64) 
 143   
144   
146          """Set up isotropic cone optimisation against the 2nd degree Frame Order matrix elements. 
147           
148          @keyword frame_order_2nd:   The numerical values of the 2nd degree Frame Order matrix.  If supplied, the target functions will optimise directly to these values. 
149          @type frame_order_2nd:      numpy 9D, rank-2 array 
150          """ 
151   
152           
153          self.data = frame_order_2nd 
154   
155           
156          self.errors = ones((9, 9), float64) 
157   
158           
159          self.cone_axis = zeros(3, float64) 
160          self.z_axis = array([0, 0, 1], float64) 
161   
162           
163          self.rot = zeros((3, 3), float64) 
164   
165           
166          self.frame_order_2nd = zeros((9, 9), float64) 
167   
168           
169          self.func = self.func_iso_cone_elements 
 170   
171   
173          """Target function for free rotor model optimisation. 
174   
175          This function optimises against alignment tensors.  The Euler angles for the tensor rotation are the first 3 parameters optimised in this model, followed by the polar and azimuthal angles of the cone axis. 
176   
177          @param params:  The vector of parameter values.  These are the tensor rotation angles {alpha, beta, gamma, theta, phi}. 
178          @type params:   list of float 
179          @return:        The chi-squared or SSE value. 
180          @rtype:         float 
181          """ 
182   
183           
184          ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi = params 
185   
186           
187          frame_order_2nd = compile_2nd_matrix_free_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi) 
188   
189           
190          self.reduce_and_rot(0.0, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
191   
192           
193          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 194   
195   
197          """Target function for isotropic cone model optimisation using the Frame Order matrix. 
198   
199          This function optimises by directly matching the elements of the 2nd degree Frame Order 
200          super matrix.  The cone axis spherical angles theta and phi and the cone angle theta are the 
201          3 parameters optimised in this model. 
202   
203          @param params:  The vector of parameter values {theta, phi, theta_cone} where the first two are the polar and azimuthal angles of the cone axis theta_cone is the isotropic cone angle. 
204          @type params:   list of float 
205          @return:        The chi-squared or SSE value. 
206          @rtype:         float 
207          """ 
208   
209           
210          theta, phi, theta_cone = params 
211   
212           
213          self.frame_order_2nd = compile_2nd_matrix_iso_cone_free_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, theta, phi, theta_cone) 
214   
215           
216          self.frame_order_2nd = self.frame_order_2nd.copy() 
217   
218           
219          self.data.shape = (81,) 
220          self.frame_order_2nd.shape = (81,) 
221          self.errors.shape = (81,) 
222   
223           
224          val = chi2(self.data, self.frame_order_2nd, self.errors) 
225   
226           
227          self.data.shape = (9, 9) 
228          self.frame_order_2nd.shape = (9, 9) 
229          self.errors.shape = (9, 9) 
230   
231           
232          return val 
 233   
234   
236          """Target function for isotropic cone model optimisation. 
237   
238          This function optimises against alignment tensors. 
239   
240          @param params:  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. 
241          @type params:   list of float 
242          @return:        The chi-squared or SSE value. 
243          @rtype:         float 
244          """ 
245   
246           
247          ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta, sigma_max = params 
248   
249           
250          frame_order_2nd = compile_2nd_matrix_iso_cone(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta, sigma_max) 
251   
252           
253          self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
254   
255           
256          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 257   
258   
260          """Target function for free rotor isotropic cone model optimisation. 
261   
262          This function optimises against alignment tensors. 
263   
264          @param params:  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. 
265          @type params:   list of float 
266          @return:        The chi-squared or SSE value. 
267          @rtype:         float 
268          """ 
269   
270           
271          ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi, cone_s1 = params 
272   
273           
274          frame_order_2nd = compile_2nd_matrix_iso_cone_free_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi, cone_s1) 
275   
276           
277          self.reduce_and_rot(0.0, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
278   
279           
280          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 281   
282   
284          """Target function for torsionless isotropic cone model optimisation. 
285   
286          This function optimises against alignment tensors. 
287   
288          @param params:  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. 
289          @type params:   list of float 
290          @return:        The chi-squared or SSE value. 
291          @rtype:         float 
292          """ 
293   
294           
295          ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi, cone_theta = params 
296   
297           
298          frame_order_2nd = compile_2nd_matrix_iso_cone_torsionless(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi, cone_theta) 
299   
300           
301          self.reduce_and_rot(0.0, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
302   
303           
304          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 305   
306   
308          """Target function for pseudo-elliptic cone model optimisation. 
309   
310          @param params:  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. 
311          @type params:   list of float 
312          @return:        The chi-squared or SSE value. 
313          @rtype:         float 
314          """ 
315   
316           
317          ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max = params 
318   
319           
320          frame_order_2nd = compile_2nd_matrix_pseudo_ellipse(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y, cone_sigma_max) 
321   
322           
323          self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
324   
325           
326          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 327   
328   
330          """Target function for free_rotor pseudo-elliptic cone model optimisation. 
331   
332          @param params:  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. 
333          @type params:   list of float 
334          @return:        The chi-squared or SSE value. 
335          @rtype:         float 
336          """ 
337   
338           
339          ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y = params 
340   
341           
342          frame_order_2nd = compile_2nd_matrix_pseudo_ellipse_free_rotor(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y) 
343   
344           
345          self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
346   
347           
348          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 349   
350   
352          """Target function for torsionless pseudo-elliptic cone model optimisation. 
353   
354          @param params:  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. 
355          @type params:   list of float 
356          @return:        The chi-squared or SSE value. 
357          @rtype:         float 
358          """ 
359   
360           
361          ave_pos_alpha, ave_pos_beta, ave_pos_gamma, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y = params 
362   
363           
364          frame_order_2nd = compile_2nd_matrix_pseudo_ellipse_torsionless(self.frame_order_2nd, self.rot, eigen_alpha, eigen_beta, eigen_gamma, cone_theta_x, cone_theta_y) 
365   
366           
367          self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
368   
369           
370          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 371   
372   
374          """Target function for rigid model optimisation. 
375   
376          This function optimises against alignment tensors.  The Euler angles for the tensor rotation are the 3 parameters optimised in this model. 
377   
378          @param params:  The vector of parameter values.  These are the tensor rotation angles {alpha, beta, gamma}. 
379          @type params:   list of float 
380          @return:        The chi-squared or SSE value. 
381          @rtype:         float 
382          """ 
383   
384           
385          ave_pos_alpha, ave_pos_beta, ave_pos_gamma = params 
386   
387           
388          self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma) 
389   
390           
391          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 392   
393   
395          """Target function for rotor model optimisation. 
396   
397          This function optimises against alignment tensors.  The Euler angles for the tensor rotation are the first 3 parameters optimised in this model, followed by the polar and azimuthal angles of the cone axis and the torsion angle restriction. 
398   
399          @param params:  The vector of parameter values.  These are the tensor rotation angles {alpha, beta, gamma, theta, phi, sigma_max}. 
400          @type params:   list of float 
401          @return:        The chi-squared or SSE value. 
402          @rtype:         float 
403          """ 
404   
405           
406          ave_pos_alpha, ave_pos_beta, ave_pos_gamma, axis_theta, axis_phi, sigma_max = params 
407   
408           
409          frame_order_2nd = compile_2nd_matrix_rotor(self.frame_order_2nd, self.rot, self.z_axis, self.cone_axis, axis_theta, axis_phi, sigma_max) 
410   
411           
412          self.reduce_and_rot(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, frame_order_2nd) 
413   
414           
415          return chi2(self.red_tensors, self.red_tensors_bc, self.red_errors) 
 416   
417   
418 -    def reduce_and_rot(self, ave_pos_alpha=None, ave_pos_beta=None, ave_pos_gamma=None, daeg=None): 
 419          """Reduce and rotate the alignments tensors using the frame order matrix and Euler angles. 
420   
421          @keyword ave_pos_alpha: The alpha Euler angle describing the average domain position, the tensor rotation. 
422          @type ave_pos_alpha:    float 
423          @keyword ave_pos_beta:  The beta Euler angle describing the average domain position, the tensor rotation. 
424          @type ave_pos_beta:     float 
425          @keyword ave_pos_gamma: The gamma Euler angle describing the average domain position, the tensor rotation. 
426          @type ave_pos_gamma:    float 
427          @keyword daeg:          The 2nd degree frame order matrix. 
428          @type daeg:             rank-2, 9D array 
429          """ 
430   
431           
432          euler_to_R(ave_pos_alpha, ave_pos_beta, ave_pos_gamma, self.rot) 
433   
434           
435          for i in range(self.num_tensors): 
436               
437              index1 = i*5 
438              index2 = i*5+5 
439   
440               
441              if daeg != None: 
442                   
443                  reduce_alignment_tensor(daeg, self.full_tensors[index1:index2], self.red_tensors_bc[index1:index2]) 
444   
445                   
446                  to_tensor(self.tensor_3D, self.red_tensors_bc[index1:index2]) 
447   
448               
449              else: 
450                   
451                  to_tensor(self.tensor_3D, self.full_tensors[index1:index2]) 
452   
453               
454              if self.full_in_ref_frame[i]: 
455                  tensor_3D = dot(self.rot, dot(self.tensor_3D, transpose(self.rot))) 
456   
457               
458              else: 
459                  tensor_3D = dot(transpose(self.rot), dot(self.tensor_3D, self.rot)) 
460   
461               
462              to_5D(self.red_tensors_bc[index1:index2], tensor_3D) 
  463