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 functions for the handling of alignment tensors."""
25
26
27 from copy import deepcopy
28 from math import pi, sqrt
29 from numpy import arccos, dot, float64, linalg, zeros
30 from numpy.linalg import norm
31 from re import search
32 import sys
33
34
35 from angles import wrap_angles
36 from data.align_tensor import AlignTensorList
37 from generic_fns import pipes
38 from physical_constants import g1H, h_bar, kB, mu0, return_gyromagnetic_ratio
39 from relax_errors import RelaxError, RelaxNoTensorError, RelaxStrError, RelaxTensorError, RelaxUnknownParamCombError, RelaxUnknownParamError
40
41
43 """Function for determining if alignment data exists in the current data pipe.
44
45 @param tensor: The alignment tensor identification string.
46 @type tensor: str
47 @param pipe: The data pipe to search for data in.
48 @type pipe: str
49 @return: The answer to the question.
50 @rtype: bool
51 """
52
53
54 if pipe == None:
55 pipe = pipes.cdp_name()
56
57
58 pipe = pipes.get_pipe(pipe)
59
60
61 if hasattr(pipe, 'align_tensors'):
62 for data in pipe.align_tensors:
63 if data.name == tensor:
64 return True
65 else:
66 return False
67
68
70 """Determine if all alignment tensors are fixed.
71
72 @return: True if all tensors are fixed, False otherwise.
73 @rtype: bool
74 """
75
76
77 for i in range(len(cdp.align_tensors)):
78
79 if not cdp.align_tensors[i].fixed:
80 return False
81
82
83 return True
84
85
87 """Convert the alignment tensor into the magnetic susceptibility (chi) tensor.
88
89 A can be either the full tensor (3D or 5D), a component Aij of the tensor, Aa, or Ar, anything that can be multiplied by the constants to convert from one to the other.
90
91
92 @param A: The alignment tensor or alignment tensor component.
93 @type A: numpy array or float
94 @param B0: The magnetic field strength in Hz.
95 @type B0: float
96 @param T: The temperature in Kalvin.
97 @type T: float
98 @return: A multiplied by the PCS constant.
99 @rtype: numpy array or float
100 """
101
102
103 B0 = 2.0 * pi * B0 / g1H
104
105
106 conv = 15.0 * mu0 * kB * T / B0**2
107
108
109 return conv * A
110
111
112 -def copy(tensor_from=None, pipe_from=None, tensor_to=None, pipe_to=None):
113 """Function for copying alignment tensor data from one data pipe to another.
114
115 @param tensor_from: The identification string of the alignment tensor to copy the data from.
116 @type tensor_from: str
117 @param pipe_from: The data pipe to copy the alignment tensor data from. This defaults to the
118 current data pipe.
119 @type pipe_from: str
120 @param tensor_to: The identification string of the alignment tensor to copy the data to.
121 @type tensor_to: str
122 @param pipe_to: The data pipe to copy the alignment tensor data to. This defaults to the
123 current data pipe.
124 @type pipe_to: str
125 """
126
127
128 if tensor_from == tensor_to and pipe_from == None and pipe_to == None:
129 raise RelaxError("The pipe_from and pipe_to arguments cannot both be set to None when the tensor names are the same.")
130 elif pipe_from == None:
131 pipe_from = pipes.cdp_name()
132 elif pipe_to == None:
133 pipe_to = pipes.cdp_name()
134
135
136 pipes.test(pipe_from)
137 pipes.test(pipe_to)
138
139
140 dp_from = pipes.get_pipe(pipe_from)
141 dp_to = pipes.get_pipe(pipe_to)
142
143
144 if not align_data_exists(tensor_from, pipe_from):
145 raise RelaxNoTensorError('alignment')
146
147
148 if align_data_exists(tensor_to, pipe_to):
149 raise RelaxTensorError('alignment')
150
151
152 if not hasattr(dp_to, 'align_tensors'):
153 dp_to.align_tensors = AlignTensorList()
154
155
156 index_from = get_tensor_index(tensor_from, pipe_from)
157 index_to = get_tensor_index(tensor_to, pipe_to)
158
159
160 if index_to == None:
161 dp_to.align_tensors.append(deepcopy(dp_from.align_tensors[index_from]))
162 else:
163 dp_to.align_tensors[index_to] = deepcopy(dp_from.align_tensors[index_from])
164
165
167 """Function for returning a list of names of data structures associated with the sequence."""
168
169 names = [ 'align_params' ]
170
171 return names
172
173
175 """Return the default values for the alignment tensor parameters.
176
177 @param param: The name of the parameter.
178 @type param: str
179 @return: The default value, which for all parameters is set to zero.
180 @rtype: float
181 """
182
183
184 return 0.0
185
186
187 __default_value_prompt_doc__ = """
188 Alignment tensor parameter default values
189 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
190
191 ________________________________________________________________________
192 | | | |
193 | Data type | Object name | Value |
194 |________________________|____________________|________________________|
195 | | | |
196 | Axx | 'Axx' | 0.0 |
197 | | | |
198 | Ayy | 'Ayy' | 0.0 |
199 | | | |
200 | Azz | 'Azz' | 0.0 |
201 | | | |
202 | Axxyy | 'Axxyy' | 0.0 |
203 | | | |
204 | Axy | 'Axy' | 0.0 |
205 | | | |
206 | Axz | 'Axz' | 0.0 |
207 | | | |
208 | Ayz | 'Ayz' | 0.0 |
209 | | | |
210 | alpha | 'alpha' | 0.0 |
211 | | | |
212 | beta | 'beta' | 0.0 |
213 | | | |
214 | gamma | 'gamma' | 0.0 |
215 |________________________|____________________|________________________|
216
217 """
218
219
221 """Function for deleting alignment tensor data.
222
223 @param tensor: The alignment tensor identification string.
224 @type tensor: str or None
225 """
226
227
228 pipes.test()
229
230
231 if tensor and not align_data_exists(tensor):
232 raise RelaxNoTensorError('alignment')
233 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'):
234 raise RelaxNoTensorError('alignment')
235
236
237 if tensor:
238 tensors = [tensor]
239 else:
240 tensors = cdp.align_ids
241
242
243 for tensor in tensors:
244
245 print("Removing the '%s' tensor." % tensor)
246
247
248 index = get_tensor_index(tensor)
249
250
251 cdp.align_tensors.pop(index)
252
253
254 if not len(cdp.align_tensors):
255 del(cdp.align_tensors)
256
257
259 """Function for displaying the alignment tensor.
260
261 @param tensor: The alignment tensor identification string.
262 @type tensor: str or None
263 """
264
265
266 pipes.test()
267
268
269 if tensor and not align_data_exists(tensor):
270 raise RelaxNoTensorError('alignment')
271 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'):
272 raise RelaxNoTensorError('alignment')
273
274
275 tensor_list = []
276 if not tensor:
277 for tensor_cont in cdp.align_tensors:
278 tensor_list.append(tensor_cont.name)
279 else:
280 tensor_list.append(tensor)
281
282
283 for tensor in tensor_list:
284
285 if not align_data_exists(tensor):
286 raise RelaxNoTensorError('alignment')
287
288
289 data = get_tensor_object(tensor)
290
291
292 head = "# Tensor: %s #" % tensor
293 print("\n\n\n" + '#' * len(head) + "\n" + head + "\n" + '#' * len(head))
294
295
296
297
298
299 title = "# Saupe order matrix."
300 print("\n\n" + title + '\n' + '#'*len(title) + '\n')
301
302
303 print("# 5D, rank-1 notation {Sxx, Syy, Sxy, Sxz, Syz}:")
304 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Sxx, data.Syy, data.Sxy, data.Sxz, data.Syz)))
305
306
307 print("# 5D, rank-1 notation {Szz, Sxx-yy, Sxy, Sxz, Syz} (the Pales default format).")
308 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Szz, data.Sxxyy, data.Sxy, data.Sxz, data.Syz)))
309
310
311 print("# 3D, rank-2 notation.")
312 print("%s" % (data.S))
313
314
315
316
317
318 title = "# Alignment tensor."
319 print("\n\n" + title + '\n' + '#'*len(title) + '\n')
320
321
322 print("# 5D, rank-1 notation {Axx, Ayy, Axy, Axz, Ayz}:")
323 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Axx, data.Ayy, data.Axy, data.Axz, data.Ayz)))
324
325
326 print("# 5D, rank-1 notation {Azz, Axx-yy, Axy, Axz, Ayz} (the Pales default format).")
327 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Azz, data.Axxyy, data.Axy, data.Axz, data.Ayz)))
328
329
330 print("# 3D, rank-2 notation.")
331 print("%s" % data.A)
332
333
334
335
336
337 title = "# Probability tensor."
338 print("\n\n" + title + '\n' + '#'*len(title) + '\n')
339
340
341 print("# 5D, rank-1 notation {Pxx, Pyy, Pxy, Pxz, Pyz}:")
342 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Pxx, data.Pyy, data.Pxy, data.Pxz, data.Pyz)))
343
344
345 print("# 5D, rank-1 notation {Pzz, Pxx-yy, Pxy, Pxz, Pyz}.")
346 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (data.Pzz, data.Pxxyy, data.Pxy, data.Pxz, data.Pyz)))
347
348
349 print("# 3D, rank-2 notation.")
350 print("%s" % data.P)
351
352
353
354
355
356 title = "# Magnetic susceptibility tensor."
357 print("\n\n" + title + '\n' + '#'*len(title) + '\n')
358 chi_tensor = True
359
360
361 print("# The magnetic field strength (MHz):")
362 if hasattr(cdp, 'frq') and tensor in cdp.frq:
363 print(("%s\n" % (cdp.frq[tensor] / 1e6)))
364 else:
365 print(("Not set.\n"))
366 chi_tensor = False
367
368
369 print("# The temperature (K):")
370 if hasattr(cdp, 'temperature') and tensor in cdp.temperature:
371 print(("%s\n" % cdp.temperature[tensor]))
372 else:
373 print(("Not set.\n"))
374 chi_tensor = False
375
376
377 if not chi_tensor:
378 print("# The chi tensor:\nN/A.\n")
379
380
381 else:
382
383 chi_xx = calc_chi_tensor(data.Axx, cdp.frq[tensor], cdp.temperature[tensor])
384 chi_xy = calc_chi_tensor(data.Axy, cdp.frq[tensor], cdp.temperature[tensor])
385 chi_xz = calc_chi_tensor(data.Axz, cdp.frq[tensor], cdp.temperature[tensor])
386 chi_yy = calc_chi_tensor(data.Ayy, cdp.frq[tensor], cdp.temperature[tensor])
387 chi_yz = calc_chi_tensor(data.Ayz, cdp.frq[tensor], cdp.temperature[tensor])
388 chi_zz = calc_chi_tensor(data.Azz, cdp.frq[tensor], cdp.temperature[tensor])
389 chi_xxyy = calc_chi_tensor(data.Axxyy, cdp.frq[tensor], cdp.temperature[tensor])
390 chi = calc_chi_tensor(data.A, cdp.frq[tensor], cdp.temperature[tensor])
391
392
393 print("# 5D, rank-1 notation {chi_xx, chi_yy, chi_xy, chi_xz, chi_yz}:")
394 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (chi_xx, chi_yy, chi_xy, chi_xz, chi_yz)))
395
396
397 print("# 5D, rank-1 notation {chi_zz, chi_xx-yy, chi_xy, chi_xz, chi_yz}.")
398 print(("[%25.12e, %25.12e, %25.12e, %25.12e, %25.12e]\n" % (chi_zz, chi_xxyy, chi_xy, chi_xz, chi_yz)))
399
400
401 print("# 3D, rank-2 notation.")
402 print("%s" % chi)
403
404
405
406
407
408 title = "# Eigensystem."
409 print("\n\n" + title + '\n' + '#'*len(title) + '\n')
410
411
412 print("# Saupe order matrix eigenvalues {Sxx, Syy, Szz}.")
413 print(("[%25.12e, %25.12e, %25.12e]\n" % (data.S_diag[0, 0], data.S_diag[1, 1], data.S_diag[2, 2])))
414 print("# Alignment tensor eigenvalues {Axx, Ayy, Azz}.")
415 print(("[%25.12e, %25.12e, %25.12e]\n" % (data.A_diag[0, 0], data.A_diag[1, 1], data.A_diag[2, 2])))
416 print("# Probability tensor eigenvalues {Pxx, Pyy, Pzz}.")
417 print(("[%25.12e, %25.12e, %25.12e]\n" % (data.P_diag[0, 0], data.P_diag[1, 1], data.P_diag[2, 2])))
418 if chi_tensor:
419 chi_diag = calc_chi_tensor(data.A_diag, cdp.frq[tensor], cdp.temperature[tensor])
420 print("# Magnetic susceptibility eigenvalues {chi_xx, chi_yy, chi_zz}.")
421 print(("[%25.12e, %25.12e, %25.12e]\n" % (chi_diag[0, 0], chi_diag[1, 1], chi_diag[2, 2])))
422
423
424 print("# Eigenvector x.")
425 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.unit_x[0], data.unit_x[1], data.unit_x[2])))
426 print("# Eigenvector y.")
427 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.unit_y[0], data.unit_y[1], data.unit_y[2])))
428 print("# Eigenvector z.")
429 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.unit_z[0], data.unit_z[1], data.unit_z[2])))
430
431
432 print("# Rotation matrix.")
433 print("%s\n" % data.rotation)
434
435
436 print("# Euler angles in zyz notation {alpha, beta, gamma}.")
437 print(("[%25.12f, %25.12f, %25.12f]\n" % (data.euler[0], data.euler[1], data.euler[2])))
438
439
440
441
442
443 title = "# Geometric description."
444 print("\n\n" + title + '\n' + '#'*len(title) + '\n')
445
446
447 print("# Generalized degree of order (GDO).")
448 print("GDO = %-25.12e\n" % gdo(data.A))
449
450
451 print("# Alignment tensor axial component (Aa = 3/2 * Azz, where Aii are the eigenvalues).")
452 print("Aa = %-25.12e\n" % data.Aa)
453
454
455 print("# Rhombic component (Ar = Axx - Ayy, where Aii are the eigenvalues).")
456 print("Ar = %-25.12e\n" % data.Ar)
457 print("# Rhombicity (R = Ar / Aa).")
458 print("R = %-25.12f\n" % data.R)
459 print("# Asymmetry parameter (eta = (Axx - Ayy) / Azz, where Aii are the eigenvalues).")
460 print("eta = %-25.12f\n" % data.eta)
461
462
463 if chi_tensor:
464
465 print("# Magnetic susceptibility axial parameter (chi_ax = chi_zz - (chi_xx + chi_yy)/2, where chi_ii are the eigenvalues).")
466 print("chi_ax = %-25.12e\n" % (chi_diag[2, 2] - (chi_diag[0, 0] + chi_diag[1, 1])/2.0))
467
468
469 print("# Magnetic susceptibility rhombicity parameter (chi_rh = chi_xx - chi_yy, where chi_ii are the eigenvalues).")
470 print("chi_rh = %-25.12e\n" % (chi_diag[0, 0] - chi_diag[1, 1]))
471
472
473 print("\n\n\n")
474
475
476 -def fix(id=None, fixed=True):
477 """Fix the alignment tensor during optimisation.
478
479 @keyword id: The alignment tensor ID string. If set to None, then all alignment tensors will be fixed.
480 @type id: str or None
481 @keyword fixed: If True, the alignment tensor will be fixed during optimisation. If False, the alignment tensors will be optimised.
482 @type fixed: bool
483 """
484
485
486 pipes.test()
487
488
489 if not hasattr(cdp, 'align_tensors') or not hasattr(cdp, 'align_ids'):
490 raise RelaxNoTensorError('alignment')
491
492
493 for i in range(len(cdp.align_tensors)):
494
495 if id and cdp.align_tensors[i].name == id:
496 cdp.align_tensors[i].fixed = fixed
497
498
499 if id == None:
500 cdp.align_tensors[i].fixed = fixed
501
502
504 """Wrap the Euler angles and remove the glide reflection and translational symmetries.
505
506 Wrap the angles such that::
507
508 0 <= alpha <= 2pi,
509 0 <= beta <= pi,
510 0 <= gamma <= 2pi.
511
512 For the simulated values, the angles are wrapped as::
513
514 alpha - pi <= alpha_sim <= alpha + pi
515 beta - pi/2 <= beta_sim <= beta + pi/2
516 gamma - pi <= gamma_sim <= gamma + pi
517
518
519 @param sim_index: The simulation index. If set to None then the actual values will be folded
520 rather than the simulation values.
521 @type sim_index: int or None
522 """
523
524
525
526
527
528
529 alpha = cdp.align_tensors.alpha
530 beta = cdp.align_tensors.beta
531 gamma = cdp.align_tensors.gamma
532
533
534 if sim_index != None:
535 alpha_sim = cdp.align_tensors.alpha_sim[sim_index]
536 beta_sim = cdp.align_tensors.beta_sim[sim_index]
537 gamma_sim = cdp.align_tensors.gamma_sim[sim_index]
538
539
540 if sim_index == None:
541 cdp.align_tensors.alpha = wrap_angles(alpha, 0.0, 2.0*pi)
542 cdp.align_tensors.beta = wrap_angles(beta, 0.0, 2.0*pi)
543 cdp.align_tensors.gamma = wrap_angles(gamma, 0.0, 2.0*pi)
544
545
546 else:
547 cdp.align_tensors.alpha_sim[sim_index] = wrap_angles(alpha_sim, alpha - pi, alpha + pi)
548 cdp.align_tensors.beta_sim[sim_index] = wrap_angles(beta_sim, beta - pi, beta + pi)
549 cdp.align_tensors.gamma_sim[sim_index] = wrap_angles(gamma_sim, gamma - pi, gamma + pi)
550
551
552
553
554
555
556 if sim_index == None:
557
558 if cdp.align_tensors.beta >= pi:
559 cdp.align_tensors.alpha = pi - cdp.align_tensors.alpha
560 cdp.align_tensors.beta = cdp.align_tensors.beta - pi
561
562
563 else:
564
565 if cdp.align_tensors.beta_sim[sim_index] >= cdp.align_tensors.beta + pi/2.0:
566 cdp.align_tensors.alpha_sim[sim_index] = pi - cdp.align_tensors.alpha_sim[sim_index]
567 cdp.align_tensors.beta_sim[sim_index] = cdp.align_tensors.beta_sim[sim_index] - pi
568 elif cdp.align_tensors.beta_sim[sim_index] <= cdp.align_tensors.beta - pi/2.0:
569 cdp.align_tensors.alpha_sim[sim_index] = pi - cdp.align_tensors.alpha_sim[sim_index]
570 cdp.align_tensors.beta_sim[sim_index] = cdp.align_tensors.beta_sim[sim_index] + pi
571
572
574 """Calculate the generalized degree of order (GDO) for the given alignment tensor.
575
576 @param A: The alignment tensor.
577 @type A: rank-2, 3D numpy array
578 @return: The GDO value.
579 @rtype: float
580 """
581
582
583 gdo = sqrt(3.0/2.0) * norm(A)
584
585
586 return gdo
587
588
590 """Function for returning the index corresponding to the 'tensor' argument.
591
592 @param tensor: The alignment tensor identification string.
593 @type tensor: str
594 @param pipe: The data pipe to search for data in.
595 @type pipe: str
596 @return: The index corresponding to the 'tensor' arg.
597 @rtype: int
598 """
599
600
601 if pipe == None:
602 pipe = pipes.cdp_name()
603
604
605 dp = pipes.get_pipe(pipe)
606
607
608 index = None
609
610
611 for i in xrange(len(dp.align_tensors)):
612 if dp.align_tensors[i].name == tensor:
613 index = i
614
615
616 return index
617
618
620 """Function for returning the AlignTensorData instance corresponding to the 'tensor' argument.
621
622 @param tensor: The alignment tensor identification string.
623 @type tensor: str
624 @param pipe: The data pipe to search for data in.
625 @type pipe: str
626 @return: The alignment tensor object corresponding to the 'tensor' arg.
627 @rtype: AlignTensorData instance
628 """
629
630
631 if pipe == None:
632 pipe = pipes.cdp_name()
633
634
635 data = None
636
637
638 for i in xrange(len(cdp.align_tensors)):
639 if cdp.align_tensors[i].name == tensor:
640 data = cdp.align_tensors[i]
641
642
643 return data
644
645
646 -def init(tensor=None, params=None, scale=1.0, angle_units='deg', param_types=0, errors=False):
647 """Function for initialising the alignment tensor.
648
649 @keyword tensor: The alignment tensor identification string.
650 @type tensor: str
651 @keyword params: The alignment tensor parameters.
652 @type params: float
653 @keyword scale: The alignment tensor eigenvalue scaling value.
654 @type scale: float
655 @keyword angle_units: The units for the angle parameters (either 'deg' or 'rad').
656 @type angle_units: str
657 @keyword param_types: The type of parameters supplied. The flag values correspond to, 0:
658 {Axx, Ayy, Axy, Axz, Ayz}, and 1: {Azz, Axx-yy, Axy, Axz, Ayz}.
659 @type param_types: int
660 @keyword errors: A flag which determines if the alignment tensor data or its errors are
661 being input.
662 @type errors: bool
663 """
664
665
666 pipes.test()
667
668
669 if errors and not align_data_exists(tensor):
670 raise RelaxNoTensorError('alignment')
671
672
673 valid_types = ['deg', 'rad']
674 if not angle_units in valid_types:
675 raise RelaxError("The alignment tensor 'angle_units' argument " + repr(angle_units) + " should be either 'deg' or 'rad'.")
676
677
678 if not hasattr(cdp, 'align_ids'):
679 cdp.align_ids = []
680 if tensor not in cdp.align_ids:
681 cdp.align_ids.append(tensor)
682
683
684 if not errors:
685
686 if not hasattr(cdp, 'align_tensors'):
687 cdp.align_tensors = AlignTensorList()
688
689
690 if tensor not in cdp.align_tensors.names():
691 cdp.align_tensors.add_item(tensor)
692
693
694 tensor_obj = get_tensor_object(tensor)
695
696
697 if param_types == 0:
698
699 Sxx, Syy, Sxy, Sxz, Syz = params
700
701
702 Sxx = Sxx * scale
703 Syy = Syy * scale
704 Sxy = Sxy * scale
705 Sxz = Sxz * scale
706 Syz = Syz * scale
707
708
709 set(tensor=tensor_obj, value=[Sxx, Syy, Sxy, Sxz, Syz], param=['Sxx', 'Syy', 'Sxy', 'Sxz', 'Syz'], errors=errors)
710
711
712 elif param_types == 1:
713
714 Szz, Sxxyy, Sxy, Sxz, Syz = params
715
716
717 Szz = Szz * scale
718 Sxxyy = Sxxyy * scale
719 Sxy = Sxy * scale
720 Sxz = Sxz * scale
721 Syz = Syz * scale
722
723
724 set(tensor=tensor_obj, value=[Szz, Sxxyy, Sxy, Sxz, Syz], param=['Szz', 'Sxxyy', 'Sxy', 'Sxz', 'Syz'], errors=errors)
725
726
727 elif param_types == 2:
728
729 Axx, Ayy, Axy, Axz, Ayz = params
730
731
732 Axx = Axx * scale
733 Ayy = Ayy * scale
734 Axy = Axy * scale
735 Axz = Axz * scale
736 Ayz = Ayz * scale
737
738
739 set(tensor=tensor_obj, value=[Axx, Ayy, Axy, Axz, Ayz], param=['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], errors=errors)
740
741
742 elif param_types == 3:
743
744 Azz, Axxyy, Axy, Axz, Ayz = params
745
746
747 Azz = Azz * scale
748 Axxyy = Axxyy * scale
749 Axy = Axy * scale
750 Axz = Axz * scale
751 Ayz = Ayz * scale
752
753
754 set(tensor=tensor_obj, value=[Azz, Axxyy, Axy, Axz, Ayz], param=['Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz'], errors=errors)
755
756
757 elif param_types == 4:
758
759 Axx, Ayy, Axy, Axz, Ayz = params
760
761
762 r = None
763 for spin in spin_loop():
764
765 if r == None:
766 r = spin.r
767
768
769 if r != spin.r:
770 raise RelaxError("Not all spins have the same bond length.")
771
772
773 scale = scale / kappa() * r**3
774 Axx = Axx * scale
775 Ayy = Ayy * scale
776 Axy = Axy * scale
777 Axz = Axz * scale
778 Ayz = Ayz * scale
779
780
781 set(tensor=tensor_obj, value=[Axx, Ayy, Axy, Axz, Ayz], param=['Axx', 'Ayy', 'Axy', 'Axz', 'Ayz'], errors=errors)
782
783
784 elif param_types == 5:
785
786 Azz, Axxyy, Axy, Axz, Ayz = params
787
788
789 r = None
790 for spin in spin_loop():
791
792 if r == None:
793 r = spin.r
794
795
796 if r != spin.r:
797 raise RelaxError("Not all spins have the same bond length.")
798
799
800 scale = scale / kappa() * r**3
801 Azz = Azz * scale
802 Axxyy = Axxyy * scale
803 Axy = Axy * scale
804 Axz = Axz * scale
805 Ayz = Ayz * scale
806
807
808 set(tensor=tensor_obj, value=[Azz, Axxyy, Axy, Axz, Ayz], param=['Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz'], errors=errors)
809
810
811 elif param_types == 6:
812
813 Pxx, Pyy, Pxy, Pxz, Pyz = params
814
815
816 Pxx = Pxx * scale
817 Pyy = Pyy * scale
818 Pxy = Pxy * scale
819 Pxz = Pxz * scale
820 Pyz = Pyz * scale
821
822
823 set(tensor=tensor_obj, value=[Pxx, Pyy, Pxy, Pxz, Pyz], param=['Pxx', 'Pyy', 'Pxy', 'Pxz', 'Pyz'], errors=errors)
824
825
826 elif param_types == 7:
827
828 Pzz, Pxxyy, Pxy, Pxz, Pyz = params
829
830
831 Pzz = Pzz * scale
832 Pxxyy = Pxxyy * scale
833 Pxy = Pxy * scale
834 Pxz = Pxz * scale
835 Pyz = Pyz * scale
836
837
838 set(tensor=tensor_obj, value=[Pzz, Pxxyy, Pxy, Pxz, Pyz], param=['Pzz', 'Pxxyy', 'Pxy', 'Pxz', 'Pyz'], errors=errors)
839
840
841 else:
842 raise RelaxUnknownParamCombError('param_types', param_types)
843
844
846 """The function for creating bounds for the mapping function."""
847
848
849 if param in ['Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz']:
850 return [-50, 50]
851
852
853 elif param == 'alpha':
854 return [0, 2*pi]
855
856
857 elif param == 'beta':
858 return [0, pi]
859
860
861 elif param == 'gamma':
862 return [0, 2*pi]
863
864
865 -def kappa(nuc1='15N', nuc2='1H'):
866 """Function for calculating the kappa constant.
867
868 The kappa constant is::
869
870 kappa = -3/(8pi^2).gI.gS.mu0.h_bar,
871
872 where gI and gS are the gyromagnetic ratios of the I and S spins, mu0 is the permeability of
873 free space, and h_bar is Planck's constant divided by 2pi.
874
875 @param nuc1: The first nucleus type.
876 @type nuc1: str
877 @param nuc2: The first nucleus type.
878 @type nuc2: str
879 @return: The kappa constant value.
880 @rtype: float
881 """
882
883
884 gI = return_gyromagnetic_ratio(nuc1)
885 gS = return_gyromagnetic_ratio(nuc2)
886
887
888 return -3.0/(8.0*pi**2) * gI * gS * mu0 * h_bar
889
890
892 """Function for creating labels, tick locations, and tick values for an OpenDX map.
893
894 @param index: The index (which isn't used here?!?).
895 @type index: int
896 @param params: The list of parameter names.
897 @type params: list of str
898 @param bounds: The bounds of the map.
899 @type bounds: list of lists (of a float and bin)
900 @param swap: An array for switching axes around.
901 @type swap: list of int
902 @param inc: The number of increments of one dimension in the map.
903 @type inc: list of int
904 """
905
906
907 labels = "{"
908 tick_locations = []
909 tick_values = []
910 n = len(params)
911 axis_incs = 5
912 loc_inc = inc / axis_incs
913
914
915 for i in xrange(n):
916
917 factor = return_conversion_factor(params[swap[i]])
918
919
920 units = return_units(params[swap[i]])
921
922
923 if units:
924 labels = labels + "\"" + params[swap[i]] + " (" + units + ")\""
925 else:
926 labels = labels + "\"" + params[swap[i]] + "\""
927
928
929 vals = bounds[swap[i], 0] / factor
930 val_inc = (bounds[swap[i], 1] - bounds[swap[i], 0]) / (axis_incs * factor)
931
932 if i < n - 1:
933 labels = labels + " "
934 else:
935 labels = labels + "}"
936
937
938 string = "{"
939 val = 0.0
940 for j in xrange(axis_incs + 1):
941 string = string + " " + repr(val)
942 val = val + loc_inc
943 string = string + " }"
944 tick_locations.append(string)
945
946
947 string = "{"
948 for j in xrange(axis_incs + 1):
949 string = string + "\"" + "%.2f" % vals + "\" "
950 vals = vals + val_inc
951 string = string + "}"
952 tick_values.append(string)
953
954 return labels, tick_locations, tick_values
955
956
958 """Function for calculating the 5D angles between the alignment tensors.
959
960 The basis set used for the 5D vector construction changes the angles calculated.
961
962 @param basis_set: The basis set to use for constructing the 5D vectors. If set to 0, the
963 basis set is {Sxx, Syy, Sxy, Sxz, Syz}. If 1, then the basis set is {Szz,
964 Sxxyy, Sxy, Sxz, Syz}.
965 @type basis_set: int
966 @param tensors: An array of tensors to apply SVD to. If None, all tensors will be used.
967 @type tensors: None or array of str
968 """
969
970
971 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'):
972 raise RelaxNoTensorError('alignment')
973
974
975 tensor_num = 0
976 for tensor in cdp.align_tensors:
977 if tensors and tensor.name not in tensors:
978 continue
979 tensor_num = tensor_num + 1
980
981
982 matrix = zeros((tensor_num, 5), float64)
983
984
985 i = 0
986 for tensor in cdp.align_tensors:
987
988 if tensors and tensor.name not in tensors:
989 continue
990
991
992 if basis_set == 0:
993
994 matrix[i, 0] = tensor.Sxx
995 matrix[i, 1] = tensor.Syy
996 matrix[i, 2] = tensor.Sxy
997 matrix[i, 3] = tensor.Sxz
998 matrix[i, 4] = tensor.Syz
999
1000
1001 elif basis_set == 1:
1002
1003 matrix[i, 0] = tensor.Szz
1004 matrix[i, 1] = tensor.Sxxyy
1005 matrix[i, 2] = tensor.Sxy
1006 matrix[i, 3] = tensor.Sxz
1007 matrix[i, 4] = tensor.Syz
1008
1009
1010 norm = linalg.norm(matrix[i])
1011 matrix[i] = matrix[i] / norm
1012
1013
1014 i = i + 1
1015
1016
1017 cdp.align_tensors.angles = zeros((tensor_num, tensor_num), float64)
1018
1019
1020 sys.stdout.write("\nData pipe: " + repr(pipes.cdp_name()) + "\n")
1021 sys.stdout.write("\n5D angles in deg between the vectors ")
1022 if basis_set == 0:
1023 sys.stdout.write("{Sxx, Syy, Sxy, Sxz, Syz}")
1024 elif basis_set == 1:
1025 sys.stdout.write("{Szz, Sxx-yy, Sxy, Sxz, Syz}")
1026 sys.stdout.write(":\n")
1027 sys.stdout.write("%8s" % '')
1028 for i in xrange(tensor_num):
1029 sys.stdout.write("%8i" % i)
1030 sys.stdout.write("\n")
1031
1032
1033 for i in xrange(tensor_num):
1034
1035 sys.stdout.write("%8i" % i)
1036
1037
1038 for j in xrange(tensor_num):
1039
1040 delta = dot(matrix[i], matrix[j])
1041
1042
1043 if delta > 1:
1044 delta = 1
1045
1046
1047 cdp.align_tensors.angles[i, j] = arccos(delta)
1048
1049
1050 sys.stdout.write("%8.1f" % (cdp.align_tensors.angles[i, j]*180.0/pi))
1051
1052
1053 sys.stdout.write("\n")
1054
1055
1057 """Count the number of tensors.
1058
1059 @keyword skip_fixed: If set to True, then only the tensors without the fixed flag will be counted. If set to False, then all tensors will be counted.
1060 @type skip_fixed: bool
1061 @return: The number of tensors (excluding fixed tensors by default).
1062 @rtype: int
1063 """
1064
1065
1066 count = 0
1067
1068
1069 for tensor_cont in cdp.align_tensors:
1070
1071 if skip_fixed and tensor_cont.fixed:
1072 continue
1073
1074
1075 count += 1
1076
1077
1078 return count
1079
1080
1081 -def reduction(full_tensor=None, red_tensor=None):
1082 """Specify which tensor is a reduction of which other tensor.
1083
1084 @param full_tensor: The full alignment tensor.
1085 @type full_tensor: str
1086 @param red_tensor: The reduced alignment tensor.
1087 @type red_tensor: str
1088 """
1089
1090
1091 pipes.test()
1092
1093
1094 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'):
1095 raise RelaxNoTensorError('alignment')
1096
1097
1098 match_full = False
1099 match_red = False
1100 i = 0
1101 for tensor_cont in cdp.align_tensors:
1102
1103 if tensor_cont.name == full_tensor:
1104 match_full = True
1105 index_full = i
1106 if tensor_cont.name == red_tensor:
1107 match_red = True
1108 index_red = i
1109
1110
1111 i = i + 1
1112
1113
1114 if not match_full:
1115 raise RelaxNoTensorError('alignment', full_tensor)
1116 if not match_red:
1117 raise RelaxNoTensorError('alignment', red_tensor)
1118
1119
1120 if not hasattr(cdp.align_tensors, 'reduction'):
1121 cdp.align_tensors.reduction = []
1122 cdp.align_tensors.reduction.append([index_full, index_red])
1123
1124
1126 """Function for returning the factor of conversion between different parameter units.
1127
1128 @param param: The parameter name.
1129 @type param: str
1130 @return: The conversion factor.
1131 @rtype: float
1132 """
1133
1134
1135 object_name = return_data_name(param)
1136
1137
1138 if object_name in ['Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz']:
1139 return 1.0
1140
1141
1142 elif object_name in ['alpha', 'beta', 'gamma']:
1143 return (2.0*pi) / 360.0
1144
1145
1146 else:
1147 return 1.0
1148
1149
1151 """Return the parameter name.
1152
1153 @param name: The name of the parameter to return the name of.
1154 @type name: str
1155 @return: The parameter name.
1156 @rtype: str
1157 """
1158
1159
1160 if not isinstance(name, str):
1161 raise RelaxStrError('name', name)
1162
1163
1164 if search('^[Ss]xx$', name):
1165 return 'Sxx'
1166
1167
1168 if search('^[Ss]yy$', name):
1169 return 'Syy'
1170
1171
1172 if search('^[Ss]zz$', name):
1173 return 'Szz'
1174
1175
1176 if search('^[Ss]xy$', name):
1177 return 'Sxy'
1178
1179
1180 if search('^[Ss]xz$', name):
1181 return 'Sxz'
1182
1183
1184 if search('^[Ss]yz$', name):
1185 return 'Syz'
1186
1187
1188 if search('^[Ss]xxyy$', name):
1189 return 'Sxxyy'
1190
1191
1192 if search('^[Aa]xx$', name):
1193 return 'Axx'
1194
1195
1196 if search('^[Aa]yy$', name):
1197 return 'Ayy'
1198
1199
1200 if search('^[Aa]zz$', name):
1201 return 'Azz'
1202
1203
1204 if search('^[Aa]xy$', name):
1205 return 'Axy'
1206
1207
1208 if search('^[Aa]xz$', name):
1209 return 'Axz'
1210
1211
1212 if search('^[Aa]yz$', name):
1213 return 'Ayz'
1214
1215
1216 if search('^[Aa]xxyy$', name):
1217 return 'Axxyy'
1218
1219
1220 if search('^[Pp]xx$', name):
1221 return 'Pxx'
1222
1223
1224 if search('^[Pp]yy$', name):
1225 return 'Pyy'
1226
1227
1228 if search('^[Pp]zz$', name):
1229 return 'Pzz'
1230
1231
1232 if search('^[Pp]xy$', name):
1233 return 'Pxy'
1234
1235
1236 if search('^[Pp]xz$', name):
1237 return 'Pxz'
1238
1239
1240 if search('^[Pp]yz$', name):
1241 return 'Pyz'
1242
1243
1244 if search('^[Pp]xxyy$', name):
1245 return 'Pxxyy'
1246
1247
1248 if search('^a$', name) or search('alpha', name):
1249 return 'alpha'
1250
1251
1252 if search('^b$', name) or search('beta', name):
1253 return 'beta'
1254
1255
1256 if search('^g$', name) or search('gamma', name):
1257 return 'gamma'
1258
1259
1260 raise RelaxUnknownParamError(name)
1261
1262
1263 __return_data_name_prompt_doc__ = """
1264 Alignment tensor parameter string matching patterns
1265 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1266
1267 ____________________________________________________________________________________________
1268 | | | |
1269 | Data type | Object name | Patterns |
1270 |________________________________________________________|______________|__________________|
1271 | | | |
1272 | The xx component of the Saupe order matrix - Sxx | 'Sxx' | '^[Sa]xx$' |
1273 | | | |
1274 | The yy component of the Saupe order matrix - Syy | 'Syy' | '^[Sa]yy$' |
1275 | | | |
1276 | The zz component of the Saupe order matrix - Szz | 'Szz' | '^[Sa]zz$' |
1277 | | | |
1278 | The xy component of the Saupe order matrix - Sxy | 'Sxy' | '^[Sa]xy$' |
1279 | | | |
1280 | The xz component of the Saupe order matrix - Sxz | 'Sxz' | '^[Sa]xz$' |
1281 | | | |
1282 | The yz component of the Saupe order matrix - Syz | 'Syz' | '^[Sa]yz$' |
1283 | | | |
1284 | The xx-yy component of the Saupe order matrix - Sxx-yy | 'Sxxyy' | '^[Sa]xxyy$' |
1285 | | | |
1286 | The xx component of the alignment tensor - Axx | 'Axx' | '^[Aa]xx$' |
1287 | | | |
1288 | The yy component of the alignment tensor - Ayy | 'Ayy' | '^[Aa]yy$' |
1289 | | | |
1290 | The zz component of the alignment tensor - Azz | 'Azz' | '^[Aa]zz$' |
1291 | | | |
1292 | The xy component of the alignment tensor - Axy | 'Axy' | '^[Aa]xy$' |
1293 | | | |
1294 | The xz component of the alignment tensor - Axz | 'Axz' | '^[Aa]xz$' |
1295 | | | |
1296 | The yz component of the alignment tensor - Ayz | 'Ayz' | '^[Aa]yz$' |
1297 | | | |
1298 | The xx-yy component of the alignment tensor - Axx-yy | 'Axxyy' | '^[Aa]xxyy$' |
1299 | | | |
1300 | The xx component of the probability matrix - Pxx | 'Pxx' | '^[Pa]xx$' |
1301 | | | |
1302 | The yy component of the probability matrix - Pyy | 'Pyy' | '^[Pa]yy$' |
1303 | | | |
1304 | The zz component of the probability matrix - Pzz | 'Pzz' | '^[Pa]zz$' |
1305 | | | |
1306 | The xy component of the probability matrix - Pxy | 'Pxy' | '^[Pa]xy$' |
1307 | | | |
1308 | The xz component of the probability matrix - Pxz | 'Pxz' | '^[Pa]xz$' |
1309 | | | |
1310 | The yz component of the probability matrix - Pyz | 'Pyz' | '^[Pa]yz$' |
1311 | | | |
1312 | The xx-yy component of the probability matrix - Pxx-yy | 'Pxxyy' | '^[Pa]xxyy$' |
1313 | | | |
1314 | The first Euler angle of the alignment tensor - alpha | 'alpha' | '^a$' or 'alpha' |
1315 | | | |
1316 | The second Euler angle of the alignment tensor - beta | 'beta' | '^b$' or 'beta' |
1317 | | | |
1318 | The third Euler angle of the alignment tensor - gamma | 'gamma' | '^g$' or 'gamma' |
1319 |________________________________________________________|______________|__________________|
1320 """
1321
1322
1324 """Return the tensor container for the given index, skipping fixed tensors if required.
1325
1326 @param index: The index of the tensor (if skip_fixed is True, then fixed tensors are not included in the index count).
1327 @type index: int
1328 @keyword skip_fixed: A flag which if True will exclude fixed tensors from the indexation.
1329 @type skip_fixed: bool
1330 @return: The tensor corresponding to the index.
1331 @rtype: data.align_tensor.AlignTensorData instance
1332 """
1333
1334
1335 count = 0
1336
1337
1338 for tensor_cont in cdp.align_tensors:
1339
1340 if skip_fixed and tensor_cont.fixed:
1341 continue
1342
1343
1344 if index == count:
1345 return tensor_cont
1346
1347
1348 count += 1
1349
1350
1351 return False
1352
1353
1355 """Function for returning a string representing the parameters units.
1356
1357 @param param: The parameter name.
1358 @type param: str
1359 @return: The string representation of the units.
1360 @rtype: str
1361 """
1362
1363
1364 object_name = return_data_name(param)
1365
1366
1367 if object_name in ['Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz']:
1368 return 'Hz'
1369
1370
1371 elif object_name in ['alpha', 'beta', 'gamma']:
1372 return 'deg'
1373
1374
1375 -def set(tensor=None, value=None, param=None, errors=False):
1376 """Set the tensor.
1377
1378 @keyword tensor: The alignment tensor object.
1379 @type tensor: AlignTensorData instance
1380 @keyword value: The list of values to set the parameters to.
1381 @type value: list of float
1382 @keyword param: The list of parameter names.
1383 @type param: list of str
1384 @keyword errors: A flag which determines if the alignment tensor data or its errors are being
1385 input.
1386 @type errors: bool
1387 """
1388
1389
1390 geo_params = []
1391 geo_values = []
1392 orient_params = []
1393 orient_values = []
1394
1395
1396 for i in xrange(len(param)):
1397
1398 param[i] = return_data_name(param[i])
1399
1400
1401 if param[i] == None:
1402 raise RelaxUnknownParamError("alignment tensor", 'None')
1403
1404
1405 if value[i] == None:
1406 value[i] = default_value(object_names[i])
1407
1408
1409 if param[i] in ['Sxx', 'Syy', 'Szz', 'Sxxyy', 'Sxy', 'Sxz', 'Syz', 'Axx', 'Ayy', 'Azz', 'Axxyy', 'Axy', 'Axz', 'Ayz', 'Pxx', 'Pyy', 'Pzz', 'Pxxyy', 'Pxy', 'Pxz', 'Pyz']:
1410 geo_params.append(param[i])
1411 geo_values.append(value[i])
1412
1413
1414 if param[i] in ['alpha', 'beta', 'gamma']:
1415 orient_params.append(param[i])
1416 orient_values.append(value[i])
1417
1418
1419
1420
1421
1422
1423 if len(geo_params) == 1:
1424
1425
1426
1427
1428 if geo_params[0] == 'Sxx':
1429 if errors:
1430 tensor.Sxx_err = geo_values[0]
1431 else:
1432 tensor.Sxx = geo_values[0]
1433
1434
1435 elif geo_params[0] == 'Syy':
1436 if errors:
1437 tensor.Syy_err = geo_values[0]
1438 else:
1439 tensor.Syy = geo_values[0]
1440
1441
1442 elif geo_params[0] == 'Sxy':
1443 if errors:
1444 tensor.Sxy_err = geo_values[0]
1445 else:
1446 tensor.Sxy = geo_values[0]
1447
1448
1449 elif geo_params[0] == 'Sxz':
1450 if errors:
1451 tensor.Sxz_err = geo_values[0]
1452 else:
1453 tensor.Sxz = geo_values[0]
1454
1455
1456 elif geo_params[0] == 'Syz':
1457 if errors:
1458 tensor.Syz_err = geo_values[0]
1459 else:
1460 tensor.Syz = geo_values[0]
1461
1462
1463
1464
1465
1466
1467 elif geo_params[0] == 'Axx':
1468 if errors:
1469 tensor.Sxx_err = 3.0/2.0 * geo_values[0]
1470 else:
1471 tensor.Sxx = 3.0/2.0 * geo_values[0]
1472
1473
1474 elif geo_params[0] == 'Ayy':
1475 if errors:
1476 tensor.Syy_err = 3.0/2.0 * geo_values[0]
1477 else:
1478 tensor.Syy = 3.0/2.0 * geo_values[0]
1479
1480
1481 elif geo_params[0] == 'Axy':
1482 if errors:
1483 tensor.Sxy_err = 3.0/2.0 * geo_values[0]
1484 else:
1485 tensor.Sxy = 3.0/2.0 * geo_values[0]
1486
1487
1488 elif geo_params[0] == 'Axz':
1489 if errors:
1490 tensor.Sxz_err = 3.0/2.0 * geo_values[0]
1491 else:
1492 tensor.Sxz = 3.0/2.0 * geo_values[0]
1493
1494
1495 elif geo_params[0] == 'Ayz':
1496 if errors:
1497 tensor.Syz_err = 3.0/2.0 * geo_values[0]
1498 else:
1499 tensor.Syz = 3.0/2.0 * geo_values[0]
1500
1501
1502
1503
1504
1505
1506 elif geo_params[0] == 'Pxx':
1507 if errors:
1508 tensor.Sxx_err = 3.0/2.0 * (geo_values[0] - 1.0/3.0)
1509 else:
1510 tensor.Sxx = 3.0/2.0 * (geo_values[0] - 1.0/3.0)
1511
1512
1513 elif geo_params[0] == 'Pyy':
1514 if errors:
1515 tensor.Syy_err = 3.0/2.0 * (geo_values[0] - 1.0/3.0)
1516 else:
1517 tensor.Syy = 3.0/2.0 * (geo_values[0] - 1.0/3.0)
1518
1519
1520 elif geo_params[0] == 'Pxy':
1521 if errors:
1522 tensor.Sxy_err = 3.0/2.0 * geo_values[0]
1523 else:
1524 tensor.Sxy = 3.0/2.0 * geo_values[0]
1525
1526
1527 elif geo_params[0] == 'Pxz':
1528 if errors:
1529 tensor.Sxz_err = 3.0/2.0 * geo_values[0]
1530 else:
1531 tensor.Sxz = 3.0/2.0 * geo_values[0]
1532
1533
1534 elif geo_params[0] == 'Pyz':
1535 if errors:
1536 tensor.Syz_err = 3.0/2.0 * geo_values[0]
1537 else:
1538 tensor.Syz = 3.0/2.0 * geo_values[0]
1539
1540
1541 else:
1542 raise RelaxError("The geometric alignment parameter " + repr(geo_params[0]) + " cannot be set.")
1543
1544
1545 elif len(geo_params) == 5:
1546
1547 if geo_params.count('Sxx') == 1 and geo_params.count('Syy') == 1 and geo_params.count('Sxy') == 1 and geo_params.count('Sxz') == 1 and geo_params.count('Syz') == 1:
1548
1549 Sxx = geo_values[geo_params.index('Sxx')]
1550 Syy = geo_values[geo_params.index('Syy')]
1551 Sxy = geo_values[geo_params.index('Sxy')]
1552 Sxz = geo_values[geo_params.index('Sxz')]
1553 Syz = geo_values[geo_params.index('Syz')]
1554
1555
1556 if errors:
1557 tensor.Axx_err = 2.0/3.0 * Sxx
1558 tensor.Ayy_err = 2.0/3.0 * Syy
1559 tensor.Axy_err = 2.0/3.0 * Sxy
1560 tensor.Axz_err = 2.0/3.0 * Sxz
1561 tensor.Ayz_err = 2.0/3.0 * Syz
1562 else:
1563 tensor.Axx = 2.0/3.0 * Sxx
1564 tensor.Ayy = 2.0/3.0 * Syy
1565 tensor.Axy = 2.0/3.0 * Sxy
1566 tensor.Axz = 2.0/3.0 * Sxz
1567 tensor.Ayz = 2.0/3.0 * Syz
1568
1569
1570 elif geo_params.count('Szz') == 1 and geo_params.count('Sxxyy') == 1 and geo_params.count('Sxy') == 1 and geo_params.count('Sxz') == 1 and geo_params.count('Syz') == 1:
1571
1572 Szz = geo_values[geo_params.index('Szz')]
1573 Sxxyy = geo_values[geo_params.index('Sxxyy')]
1574 Sxy = geo_values[geo_params.index('Sxy')]
1575 Sxz = geo_values[geo_params.index('Sxz')]
1576 Syz = geo_values[geo_params.index('Syz')]
1577
1578
1579 if errors:
1580 tensor.Axx_err = 2.0/3.0 * -0.5*(Szz-Sxxyy)
1581 tensor.Ayy_err = 2.0/3.0 * -0.5*(Szz+Sxxyy)
1582 tensor.Axy_err = 2.0/3.0 * Sxy
1583 tensor.Axz_err = 2.0/3.0 * Sxz
1584 tensor.Ayz_err = 2.0/3.0 * Syz
1585 else:
1586 tensor.Axx = 2.0/3.0 * -0.5*(Szz-Sxxyy)
1587 tensor.Ayy = 2.0/3.0 * -0.5*(Szz+Sxxyy)
1588 tensor.Axy = 2.0/3.0 * Sxy
1589 tensor.Axz = 2.0/3.0 * Sxz
1590 tensor.Ayz = 2.0/3.0 * Syz
1591
1592
1593 elif geo_params.count('Axx') == 1 and geo_params.count('Ayy') == 1 and geo_params.count('Axy') == 1 and geo_params.count('Axz') == 1 and geo_params.count('Ayz') == 1:
1594
1595 Axx = geo_values[geo_params.index('Axx')]
1596 Ayy = geo_values[geo_params.index('Ayy')]
1597 Axy = geo_values[geo_params.index('Axy')]
1598 Axz = geo_values[geo_params.index('Axz')]
1599 Ayz = geo_values[geo_params.index('Ayz')]
1600
1601
1602 if errors:
1603 tensor.Axx_err = Axx
1604 tensor.Ayy_err = Ayy
1605 tensor.Axy_err = Axy
1606 tensor.Axz_err = Axz
1607 tensor.Ayz_err = Ayz
1608 else:
1609 tensor.Axx = Axx
1610 tensor.Ayy = Ayy
1611 tensor.Axy = Axy
1612 tensor.Axz = Axz
1613 tensor.Ayz = Ayz
1614
1615
1616 elif geo_params.count('Azz') == 1 and geo_params.count('Axxyy') == 1 and geo_params.count('Axy') == 1 and geo_params.count('Axz') == 1 and geo_params.count('Ayz') == 1:
1617
1618 Azz = geo_values[geo_params.index('Azz')]
1619 Axxyy = geo_values[geo_params.index('Axxyy')]
1620 Axy = geo_values[geo_params.index('Axy')]
1621 Axz = geo_values[geo_params.index('Axz')]
1622 Ayz = geo_values[geo_params.index('Ayz')]
1623
1624
1625 if errors:
1626 tensor.Axx_err = -0.5*(Azz-Axxyy)
1627 tensor.Ayy_err = -0.5*(Azz+Axxyy)
1628 tensor.Axy_err = Axy
1629 tensor.Axz_err = Axz
1630 tensor.Ayz_err = Ayz
1631 else:
1632 tensor.Axx = -0.5*(Azz-Axxyy)
1633 tensor.Ayy = -0.5*(Azz+Axxyy)
1634 tensor.Axy = Axy
1635 tensor.Axz = Axz
1636 tensor.Ayz = Ayz
1637
1638
1639 elif geo_params.count('Pxx') == 1 and geo_params.count('Pyy') == 1 and geo_params.count('Pxy') == 1 and geo_params.count('Pxz') == 1 and geo_params.count('Pyz') == 1:
1640
1641 Pxx = geo_values[geo_params.index('Pxx')]
1642 Pyy = geo_values[geo_params.index('Pyy')]
1643 Pxy = geo_values[geo_params.index('Pxy')]
1644 Pxz = geo_values[geo_params.index('Pxz')]
1645 Pyz = geo_values[geo_params.index('Pyz')]
1646
1647
1648 if errors:
1649 tensor.Axx_err = Pxx - 1.0/3.0
1650 tensor.Ayy_err = Pyy - 1.0/3.0
1651 tensor.Axy_err = Pxy
1652 tensor.Axz_err = Pxz
1653 tensor.Ayz_err = Pyz
1654 else:
1655 tensor.Axx = Pxx - 1.0/3.0
1656 tensor.Ayy = Pyy - 1.0/3.0
1657 tensor.Axy = Pxy
1658 tensor.Axz = Pxz
1659 tensor.Ayz = Pyz
1660
1661
1662 elif geo_params.count('Pzz') == 1 and geo_params.count('Pxxyy') == 1 and geo_params.count('Pxy') == 1 and geo_params.count('Pxz') == 1 and geo_params.count('Pyz') == 1:
1663
1664 Pzz = geo_values[geo_params.index('Pzz')]
1665 Pxxyy = geo_values[geo_params.index('Pxxyy')]
1666 Pxy = geo_values[geo_params.index('Pxy')]
1667 Pxz = geo_values[geo_params.index('Pxz')]
1668 Pyz = geo_values[geo_params.index('Pyz')]
1669
1670
1671 if errors:
1672 tensor.Axx_err = -0.5*(Pzz-Pxxyy) - 1.0/3.0
1673 tensor.Ayy_err = -0.5*(Pzz+Pxxyy) - 1.0/3.0
1674 tensor.Axy_err = Pxy
1675 tensor.Axz_err = Pxz
1676 tensor.Ayz_err = Pyz
1677 else:
1678 tensor.Axx = -0.5*(Pzz-Pxxyy) - 1.0/3.0
1679 tensor.Ayy = -0.5*(Pzz+Pxxyy) - 1.0/3.0
1680 tensor.Axy = Pxy
1681 tensor.Axz = Pxz
1682 tensor.Ayz = Pyz
1683
1684
1685 else:
1686 raise RelaxUnknownParamCombError('geometric parameter set', geo_params)
1687
1688
1689
1690 else:
1691 raise RelaxUnknownParamCombError('geometric parameter set', geo_params)
1692
1693
1694
1695
1696
1697
1698 if len(orient_params) == 1:
1699
1700 if orient_params[0] == 'alpha':
1701 if errors:
1702 tensor.alpha_err = orient_values[orient_params.index('alpha')]
1703 else:
1704 tensor.alpha = orient_values[orient_params.index('alpha')]
1705
1706
1707 elif orient_params[0] == 'beta':
1708 if errors:
1709 tensor.beta_err = orient_values[orient_params.index('beta')]
1710 else:
1711 tensor.beta = orient_values[orient_params.index('beta')]
1712
1713
1714 elif orient_params[0] == 'gamma':
1715 if errors:
1716 tensor.gamma_err = orient_values[orient_params.index('gamma')]
1717 else:
1718 tensor.gamma = orient_values[orient_params.index('gamma')]
1719
1720
1721 elif len(orient_params) == 2:
1722
1723 if orient_params.count('alpha') == 1 and orient_params.count('beta') == 1:
1724 if errors:
1725 tensor.alpha_err = orient_values[orient_params.index('alpha')]
1726 tensor.beta_err = orient_values[orient_params.index('beta')]
1727 else:
1728 tensor.alpha = orient_values[orient_params.index('alpha')]
1729 tensor.beta = orient_values[orient_params.index('beta')]
1730
1731
1732 if orient_params.count('alpha') == 1 and orient_params.count('gamma') == 1:
1733 if errors:
1734 tensor.alpha_err = orient_values[orient_params.index('alpha')]
1735 tensor.gamma_err = orient_values[orient_params.index('gamma')]
1736 else:
1737 tensor.alpha = orient_values[orient_params.index('alpha')]
1738 tensor.gamma = orient_values[orient_params.index('gamma')]
1739
1740
1741 if orient_params.count('beta') == 1 and orient_params.count('gamma') == 1:
1742 if errors:
1743 tensor.beta_err = orient_values[orient_params.index('beta')]
1744 tensor.gamma_err = orient_values[orient_params.index('gamma')]
1745 else:
1746 tensor.beta = orient_values[orient_params.index('beta')]
1747 tensor.gamma = orient_values[orient_params.index('gamma')]
1748
1749
1750 else:
1751 raise RelaxUnknownParamCombError('orientational parameter set', orient_params)
1752
1753
1754 elif len(orient_params) == 3:
1755
1756 if orient_params.count('alpha') == 1 and orient_params.count('beta') == 1:
1757 if errors:
1758 tensor.alpha_err = orient_values[orient_params.index('alpha')]
1759 tensor.beta_err = orient_values[orient_params.index('beta')]
1760 tensor.gamma_err = orient_values[orient_params.index('gamma')]
1761 else:
1762 tensor.alpha = orient_values[orient_params.index('alpha')]
1763 tensor.beta = orient_values[orient_params.index('beta')]
1764 tensor.gamma = orient_values[orient_params.index('gamma')]
1765
1766
1767 else:
1768 raise RelaxUnknownParamCombError('orientational parameter set', orient_params)
1769
1770
1771 elif len(orient_params) > 3:
1772 raise RelaxUnknownParamCombError('orientational parameter set', orient_params)
1773
1774
1775
1776
1777
1778 if orient_params:
1779 fold_angles()
1780
1781
1782 __set_prompt_doc__ = """
1783 Alignment tensor set details
1784 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1785
1786 If the alignment tensor has not been setup, use the more powerful function
1787 'alignment_tensor.init' to initialise the tensor parameters.
1788
1789 The alignment tensor parameters can only be set when the data pipe corresponds to model-free
1790 analysis. The units of the parameters are:
1791
1792 Unitless for Sxx, Syy, Szz, Sxxyy, Sxy, Sxz, Syz.
1793 Unitless for Axx, Ayy, Azz, Axxyy, Axy, Axz, Ayz.
1794 Unitless for Pxx, Pyy, Pzz, Pxxyy, Pxy, Pxz, Pyz.
1795 Radians for all angles (alpha, beta, gamma).
1796
1797 If a single geometric parameter is supplied, it must be one of Bxx, Byy, Bxy, Bxz, Byz, where B
1798 is one of S, A, or P. For the parameters Bzz and Bxxyy, it is not possible to determine how to
1799 use the currently set values together with the supplied value to calculate the new internal
1800 parameters. When supplying multiple geometric parameters, the set must belong to one of
1801
1802 {Sxx, Syy, Sxy, Sxz, Syz},
1803 {Szz, Sxxyy, Sxy, Sxz, Syz}.
1804 {Axx, Ayy, Axy, Axz, Ayz},
1805 {Azz, Axxyy, Axy, Axz, Ayz}.
1806 {Pxx, Pyy, Pxy, Pxz, Pyz},
1807 {Pzz, Pxxyy, Pxy, Pxz, Pyz}.
1808 """
1809
1810
1811 -def set_domain(tensor=None, domain=None):
1812 """Set the domain label for the given tensor.
1813
1814 @param tensor: The alignment tensor label.
1815 @type tensor: str
1816 @param domain: The domain label.
1817 @type domain: str
1818 """
1819
1820
1821 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'):
1822 raise RelaxNoTensorError('alignment')
1823
1824
1825 match = False
1826 for tensor_cont in cdp.align_tensors:
1827
1828 if tensor_cont.name == tensor:
1829 tensor_cont.domain = domain
1830 match = True
1831
1832
1833 if not match:
1834 raise RelaxNoTensorError('alignment', tensor)
1835
1836
1837 -def svd(basis_set=0, tensors=None):
1838 """Function for calculating the singular values of all the loaded tensors.
1839
1840 The matrix on which SVD will be performed is::
1841
1842 | Sxx1 Syy1 Sxy1 Sxz1 Syz1 |
1843 | Sxx2 Syy2 Sxy2 Sxz2 Syz2 |
1844 | Sxx3 Syy3 Sxy3 Sxz3 Syz3 |
1845 | . . . . . |
1846 | . . . . . |
1847 | . . . . . |
1848 | SxxN SyyN SxyN SxzN SyzN |
1849
1850 This is the default unitary basis set (selected when basis_set is 0). Alternatively a geometric
1851 basis set consisting of the stretching and skewing parameters Szz and Sxx-yy respectively
1852 replacing Sxx and Syy can be chosen by setting basis_set to 1. The matrix in this case is::
1853
1854 | Szz1 Sxxyy1 Sxy1 Sxz1 Syz1 |
1855 | Szz2 Sxxyy2 Sxy2 Sxz2 Syz2 |
1856 | Szz3 Sxxyy3 Sxy3 Sxz3 Syz3 |
1857 | . . . . . |
1858 | . . . . . |
1859 | . . . . . |
1860 | SzzN SxxyyN SxyN SxzN SyzN |
1861
1862 The relationships between the geometric and unitary basis sets are::
1863
1864 Szz = - Sxx - Syy,
1865 Sxxyy = Sxx - Syy,
1866
1867 The SVD values and condition number are dependendent upon the basis set chosen.
1868
1869
1870 @param basis_set: The basis set to create the 5 by n matrix on which to perform SVD.
1871 @type basis_set: int
1872 @param tensors: An array of tensors to apply SVD to. If None, all tensors will be used.
1873 @type tensors: None or array of str
1874 """
1875
1876
1877 if not hasattr(cdp, 'align_tensors') or len(cdp.align_tensors) == 0 or not hasattr(cdp, 'align_ids'):
1878 raise RelaxNoTensorError('alignment')
1879
1880
1881 tensor_num = 0
1882 for tensor in cdp.align_tensors:
1883 if tensors and tensor.name not in tensors:
1884 continue
1885 tensor_num = tensor_num + 1
1886
1887
1888 matrix = zeros((tensor_num, 5), float64)
1889
1890
1891 i = 0
1892 for tensor in cdp.align_tensors:
1893
1894 if tensors and tensor.name not in tensors:
1895 continue
1896
1897
1898 if basis_set == 0:
1899 matrix[i, 0] = tensor.Sxx
1900 matrix[i, 1] = tensor.Syy
1901 matrix[i, 2] = tensor.Sxy
1902 matrix[i, 3] = tensor.Sxz
1903 matrix[i, 4] = tensor.Syz
1904
1905
1906 elif basis_set == 1:
1907 matrix[i, 0] = tensor.Szz
1908 matrix[i, 1] = tensor.Sxxyy
1909 matrix[i, 2] = tensor.Sxy
1910 matrix[i, 3] = tensor.Sxz
1911 matrix[i, 4] = tensor.Syz
1912
1913
1914 i = i + 1
1915
1916
1917 u, s, vh = linalg.svd(matrix)
1918
1919
1920 cdp.align_tensors.singular_vals = s
1921
1922
1923 cdp.align_tensors.cond_num = s[0] / s[-1]
1924
1925
1926 print(("\nData pipe: " + repr(pipes.cdp_name())))
1927 print("\nSingular values:")
1928 for val in s:
1929 print(("\t" + repr(val)))
1930 print(("\nCondition number: " + repr(cdp.align_tensors.cond_num)))
1931