1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23 """Module for handling the geometric representation of the frame order motions."""
24
25
26 from copy import deepcopy
27 from math import pi
28 from numpy import cross, dot, eye, float64, zeros
29 from numpy.linalg import norm
30 import sys
31 from warnings import warn
32
33
34 from lib.errors import RelaxFault
35 from lib.frame_order.conversions import create_rotor_axis_alpha, create_rotor_axis_spherical
36 from lib.frame_order.variables import MODEL_DOUBLE_ROTOR, MODEL_FREE_ROTOR, MODEL_ISO_CONE, MODEL_ISO_CONE_FREE_ROTOR, MODEL_ISO_CONE_TORSIONLESS, MODEL_LIST_DOUBLE, MODEL_LIST_FREE_ROTORS, MODEL_LIST_ISO_CONE, MODEL_LIST_PSEUDO_ELLIPSE, MODEL_PSEUDO_ELLIPSE, MODEL_PSEUDO_ELLIPSE_FREE_ROTOR, MODEL_PSEUDO_ELLIPSE_TORSIONLESS, MODEL_ROTOR
37 from lib.geometry.rotations import euler_to_R_zyz
38 from lib.io import open_write_file
39 from lib.structure.cones import Iso_cone, Pseudo_elliptic
40 from lib.structure.geometric import generate_vector_residues
41 from lib.structure.internal.object import Internal
42 from lib.structure.represent.cone import cone
43 from lib.structure.represent.rotor import rotor
44 from lib.text.sectioning import subsection, subsubsection
45 from lib.warnings import RelaxWarning
46 from pipe_control.structure.mass import pipe_centre_of_mass
47 from specific_analyses.frame_order.checks import check_parameters
48 from specific_analyses.frame_order.data import domain_moving, generate_pivot
49
50
51 -def add_axes(structure=None, representation=None, size=None, sims=False):
52 """Add the axis system for the current frame order model to the structural object.
53
54 @keyword structure: The internal structural object to add the rotor objects to.
55 @type structure: lib.structure.internal.object.Internal instance
56 @keyword representation: The representation to create. If this is set to None or 'A', the standard representation will be created. If set to 'B', the axis system will be inverted.
57 @type representation: None or str
58 @keyword size: The size of the geometric object in Angstroms.
59 @type size: float
60 @keyword sims: A flag which if True will add the Monte Carlo simulation rotors to the structural object. There must be one model for each Monte Carlo simulation already present in the structural object.
61 @type sims: bool
62 """
63
64
65 mol_name = 'axes'
66 structure.add_molecule(name=mol_name)
67
68
69 if representation == 'B':
70 T = -eye(3)
71 else:
72 T = eye(3)
73
74
75 model_nums = [None]
76 sim_indices = [None]
77 if sims:
78 model_nums = [i+1 for i in range(cdp.sim_number)]
79 sim_indices = list(range(cdp.sim_number))
80
81
82 for i in range(len(model_nums)):
83
84 mol = structure.get_molecule(mol_name, model=model_nums[i])
85
86
87 pivot1 = generate_pivot(order=1, sim_index=sim_indices[i], pdb_limit=True)
88 pivot2 = generate_pivot(order=2, sim_index=sim_indices[i], pdb_limit=True)
89
90
91 if cdp.model in [MODEL_ISO_CONE_TORSIONLESS]:
92
93 print("\nGenerating the z-axis system.")
94
95
96 if sims:
97 axis = create_rotor_axis_spherical(theta=cdp.axis_theta_sim[sim_indices[i]], phi=cdp.axis_phi_sim[sim_indices[i]])
98 else:
99 axis = create_rotor_axis_spherical(theta=cdp.axis_theta, phi=cdp.axis_phi)
100 print(("Central axis: %s." % axis))
101
102
103 axis = dot(T, axis)
104
105
106 print("\nGenerating the axis vectors.")
107 res_num = generate_vector_residues(mol=mol, vector=axis, atom_name='z-ax', res_name_vect='AXE', res_num=2, origin=pivot1, scale=size)
108
109
110 elif cdp.model in [MODEL_DOUBLE_ROTOR]:
111
112 print("\nGenerating the z-axis linking the two pivot points.")
113
114
115 axis = pivot1 - pivot2
116 print(("Interconnecting axis: %s." % axis))
117
118
119 print("\nGenerating the axis vectors.")
120 res_num = generate_vector_residues(mol=mol, vector=axis, atom_name='z-ax', res_name_vect='AXE', res_num=1, origin=pivot2)
121
122
123 elif cdp.model in MODEL_LIST_PSEUDO_ELLIPSE:
124
125 print("\nGenerating the full axis system.")
126
127
128 mol.atom_add(atom_num=1, atom_name='R', res_name='AXE', res_num=1, pos=pivot1, element='C', pdb_record='HETATM')
129
130
131 axes = generate_axis_system(sim_index=sim_indices[i])
132
133
134 axes = dot(T, axes)
135
136
137 label = ['x', 'y']
138 if cdp.model in [MODEL_PSEUDO_ELLIPSE_TORSIONLESS]:
139 label = ['x', 'y', 'z']
140
141
142 print("\nGenerating the axis vectors.")
143 for j in range(len(label)):
144 res_num = generate_vector_residues(mol=mol, vector=axes[:, j], atom_name='%s-ax'%label[j], res_name_vect='AXE', res_num=2, origin=pivot1, scale=size)
145
146
147 -def add_cones(structure=None, representation=None, size=None, inc=None, sims=False):
148 """Add the cone geometric object for the current frame order model to the structural object.
149
150 @keyword structure: The internal structural object to add the rotor objects to.
151 @type structure: lib.structure.internal.object.Internal instance
152 @keyword representation: The representation to create. If this is set to None or 'A', the standard representation will be created. If set to 'B', the axis system will be inverted.
153 @type representation: str
154 @keyword size: The size of the geometric object in Angstroms.
155 @type size: float
156 @keyword inc: The number of increments for the filling of the cone objects.
157 @type inc: int
158 @keyword sims: A flag which if True will add the Monte Carlo simulation pivots to the structural object. There must be one model for each Monte Carlo simulation already present in the structural object.
159 @type sims: bool
160 """
161
162
163 structure.add_molecule(name='cones')
164
165
166 if representation == 'B':
167 T = -eye(3)
168 else:
169 T = eye(3)
170
171
172 model_nums = [None]
173 sim_indices = [None]
174 if sims:
175 model_nums = [i+1 for i in range(cdp.sim_number)]
176 sim_indices = list(range(cdp.sim_number))
177
178
179 for i in range(len(model_nums)):
180
181 mol = structure.get_molecule('cones', model=model_nums[i])
182
183
184 pivot = generate_pivot(order=1, sim_index=sim_indices[i], pdb_limit=True)
185
186
187 R = generate_axis_system(sim_index=sim_indices[i])
188 print(("Rotation matrix:\n%s" % R))
189
190
191 R = dot(T, R)
192
193
194 if cdp.model in MODEL_LIST_PSEUDO_ELLIPSE:
195 if sims:
196 cone_obj = Pseudo_elliptic(cdp.cone_theta_x_sim[sim_indices[i]], cdp.cone_theta_y_sim[sim_indices[i]])
197 else:
198 cone_obj = Pseudo_elliptic(cdp.cone_theta_x, cdp.cone_theta_y)
199
200
201 else:
202
203 if hasattr(cdp, 'cone_theta'):
204 if sims:
205 cone_theta = cdp.cone_theta_sim[sim_indices[i]]
206 else:
207 cone_theta = cdp.cone_theta
208
209
210 cone_obj = Iso_cone(cone_theta)
211
212
213 cone(mol=mol, cone_obj=cone_obj, start_res=1, apex=pivot, R=R, inc=inc, scale=size, distribution='regular', axis_flag=False)
214
215
217 """Add the pivots for the current frame order model to the structural object.
218
219 @keyword structure: The internal structural object to add the rotor objects to.
220 @type structure: lib.structure.internal.object.Internal instance
221 @keyword sims: A flag which if True will add the Monte Carlo simulation pivots to the structural object. There must be one model for each Monte Carlo simulation already present in the structural object.
222 @type sims: bool
223 """
224
225
226 pivots = []
227 mols = []
228 atom_names = []
229 atom_nums = []
230
231
232 mol_name = 'pivots'
233 structure.add_molecule(name=mol_name)
234
235
236 model_nums = [None]
237 sim_indices = [None]
238 if sims:
239 model_nums = [i+1 for i in range(cdp.sim_number)]
240 sim_indices = list(range(cdp.sim_number))
241
242
243 for i in range(len(model_nums)):
244
245 mol = structure.get_molecule(mol_name, model=model_nums[i])
246
247
248 pivot1 = generate_pivot(order=1, sim_index=sim_indices[i], pdb_limit=True)
249 pivot2 = generate_pivot(order=2, sim_index=sim_indices[i], pdb_limit=True)
250
251
252 if cdp.model in [MODEL_DOUBLE_ROTOR]:
253
254 mols.append(mol)
255 pivots.append(pivot1)
256 atom_names.append('Piv1')
257 atom_nums.append(1)
258
259
260 mols.append(mol)
261 pivots.append(pivot2)
262 atom_names.append('Piv2')
263 atom_nums.append(2)
264
265
266 else:
267 mols.append(mol)
268 pivots.append(pivot1)
269 atom_names.append('Piv')
270 atom_nums.append(1)
271
272
273 for i in range(len(mols)):
274 mols[i].atom_add(atom_num=atom_nums[i], atom_name=atom_names[i], res_name='PIV', res_num=1, pos=pivots[i], element='C', pdb_record='HETATM')
275
276
277 -def add_titles(structure=None, representation=None, displacement=40.0, sims=False):
278 """Add atoms to be used for the titles for the frame order geometric objects.
279
280 @keyword structure: The internal structural object to add the rotor objects to.
281 @type structure: lib.structure.internal.object.Internal instance
282 @keyword representation: The representation to title.
283 @type representation: None or str
284 @keyword displacement: The distance away from the pivot point, in Angstrom, to place the title. The simulation title will be shifted by a few extra Angstrom to avoid clashes.
285 @type displacement: float
286 @keyword sims: A flag which if True will add the Monte Carlo simulation pivots to the structural object. There must be one model for each Monte Carlo simulation already present in the structural object.
287 @type sims: bool
288 """
289
290
291 atom_name = None
292 if representation == None and sims:
293 atom_name = 'mc'
294 elif representation == 'A':
295 if sims:
296 atom_name = 'mc-a'
297 else:
298 atom_name = 'a'
299 elif representation == 'B':
300 if sims:
301 atom_name = 'mc-b'
302 else:
303 atom_name = 'b'
304
305
306 if atom_name == None:
307 return
308
309
310 if sims:
311 displacement += 3
312
313
314 mol_name = 'titles'
315 structure.add_molecule(name=mol_name)
316
317
318 if representation == 'B':
319 T = -eye(3)
320 else:
321 T = eye(3)
322
323
324 pivot1 = generate_pivot(order=1, pdb_limit=True)
325 pivot2 = generate_pivot(order=2, pdb_limit=True)
326
327
328 model_nums = [None]
329 sim_indices = [None]
330 if sims:
331 model_nums = [i+1 for i in range(cdp.sim_number)]
332 sim_indices = list(range(cdp.sim_number))
333
334
335 for i in range(len(model_nums)):
336
337 mol = structure.get_molecule(mol_name, model=model_nums[i])
338
339
340 axes = generate_axis_system(sim_index=sim_indices[i])
341
342
343 axis = dot(T, axes[:, 2])
344
345
346 pos = pivot1 + displacement*axis
347
348
349 mol.atom_add(atom_name=atom_name, res_name='TLE', res_num=1, pos=pos, element='Ti', pdb_record='HETATM')
350
351
352 -def add_rotors(structure=None, representation=None, size=None, sims=False):
353 """Add all rotor objects for the current frame order model to the structural object.
354
355 @keyword structure: The internal structural object to add the rotor objects to.
356 @type structure: lib.structure.internal.object.Internal instance
357 @keyword representation: The representation to create. If this is set to None, the rotor will be dumbbell shaped with propellers at both ends. If set to 'A' or 'B', only half of the rotor will be shown.
358 @type representation: None or str
359 @keyword size: The size of the geometric object in Angstroms.
360 @type size: float
361 @keyword sims: A flag which if True will add the Monte Carlo simulation rotors to the structural object. There must be one model for each Monte Carlo simulation already present in the structural object.
362 @type sims: bool
363 """
364
365
366 axis = []
367 span = []
368 staggered = []
369 pivot = []
370 rotor_angle = []
371 com = []
372 label = []
373 models = []
374
375
376 half_rotor = True
377 if representation == None:
378 half_rotor = False
379
380
381 if representation == 'B':
382 T = -eye(3)
383 else:
384 T = eye(3)
385
386
387 model_nums = [None]
388 sim_indices = [None]
389 if sims:
390 model_nums = [i+1 for i in range(cdp.sim_number)]
391 sim_indices = list(range(cdp.sim_number))
392
393
394 for i in range(len(model_nums)):
395
396 pivot1 = generate_pivot(order=1, sim_index=sim_indices[i], pdb_limit=True)
397 pivot2 = generate_pivot(order=2, sim_index=sim_indices[i], pdb_limit=True)
398
399
400 if cdp.model in [MODEL_ROTOR, MODEL_FREE_ROTOR, MODEL_ISO_CONE, MODEL_ISO_CONE_FREE_ROTOR, MODEL_PSEUDO_ELLIPSE, MODEL_PSEUDO_ELLIPSE_FREE_ROTOR]:
401
402 if cdp.model in MODEL_LIST_FREE_ROTORS:
403 rotor_angle.append(pi)
404 else:
405 if sims:
406 rotor_angle.append(cdp.cone_sigma_max_sim[sim_indices[i]])
407 else:
408 rotor_angle.append(cdp.cone_sigma_max)
409
410
411 if cdp.model in [MODEL_ROTOR, MODEL_FREE_ROTOR]:
412 com.append(pipe_centre_of_mass(verbosity=0, missing_error=False))
413 else:
414 com.append(pivot1)
415
416
417 axes = generate_axis_system(sim_index=sim_indices[i])
418 axis.append(axes[:, 2])
419
420
421 if cdp.model in [MODEL_ROTOR, MODEL_FREE_ROTOR]:
422 span.append(size)
423 else:
424 span.append(size + 5.0)
425
426
427 if cdp.model in MODEL_LIST_FREE_ROTORS:
428 staggered.append(False)
429 else:
430 staggered.append(True)
431
432
433 pivot.append(pivot1)
434
435
436 label.append('z-ax')
437
438
439 if sims:
440 models.append(model_nums[i])
441 else:
442 models.append(None)
443
444
445 elif cdp.model in [MODEL_DOUBLE_ROTOR]:
446
447 if sims:
448 rotor_angle.append(cdp.cone_sigma_max_2_sim[sim_indices[i]])
449 rotor_angle.append(cdp.cone_sigma_max_sim[sim_indices[i]])
450 else:
451 rotor_angle.append(cdp.cone_sigma_max_2)
452 rotor_angle.append(cdp.cone_sigma_max)
453
454
455 com.append(pivot2)
456 com.append(pivot1)
457
458
459 frame = generate_axis_system(sim_index=sim_indices[i])
460
461
462 axis.append(frame[:, 0])
463 axis.append(frame[:, 1])
464
465
466 span.append(size)
467 span.append(size)
468
469
470 staggered.append(True)
471 staggered.append(True)
472
473
474 pivot.append(pivot2)
475 pivot.append(pivot1)
476
477
478 label.append('x-ax')
479 label.append('y-ax')
480
481
482 if sims:
483 models.append(model_nums[i])
484 models.append(model_nums[i])
485 else:
486 models.append(None)
487 models.append(None)
488
489
490 for i in range(len(axis)):
491 rotor(structure=structure, rotor_angle=rotor_angle[i], axis=dot(T, axis[i]), axis_pt=pivot[i], label=label[i], centre=com[i], span=span[i]/1e10, blade_length=5e-10, model_num=models[i], staggered=staggered[i], half_rotor=half_rotor)
492
493
495 """Shift the given structural object to the average position.
496
497 @keyword structure: The structural object to operate on.
498 @type structure: lib.structure.internal.object.Internal instance
499 @keyword models: The list of models to shift.
500 @type models: list of int
501 @keyword sim: A flag which if True will use the Monte Carlo simulation results. In this case, the model list should be set to the simulation indices plus 1 and the structural object should have one model per simulation already set up.
502 @type sim: bool
503 """
504
505
506 selection = structure.selection(atom_id=domain_moving())
507
508
509 for i in range(len(models)):
510
511 R = zeros((3, 3), float64)
512 if hasattr(cdp, 'ave_pos_alpha'):
513 if sim:
514 euler_to_R_zyz(cdp.ave_pos_alpha_sim[i], cdp.ave_pos_beta_sim[i], cdp.ave_pos_gamma_sim[i], R)
515 else:
516 euler_to_R_zyz(cdp.ave_pos_alpha, cdp.ave_pos_beta, cdp.ave_pos_gamma, R)
517 else:
518 if sim:
519 euler_to_R_zyz(0.0, cdp.ave_pos_beta_sim[i], cdp.ave_pos_gamma_sim[i], R)
520 else:
521 euler_to_R_zyz(0.0, cdp.ave_pos_beta, cdp.ave_pos_gamma, R)
522 origin = pipe_centre_of_mass(atom_id=domain_moving(), verbosity=0, missing_error=False)
523 structure.rotate(R=R, origin=origin, model=models[i], selection=selection)
524
525
526 if sim:
527 T = [cdp.ave_pos_x_sim[i], cdp.ave_pos_y_sim[i], cdp.ave_pos_z_sim[i]]
528 else:
529 T = [cdp.ave_pos_x, cdp.ave_pos_y, cdp.ave_pos_z]
530 structure.translate(T=T, model=models[i], selection=selection)
531
532
533 -def create_ave_pos(format='PDB', file=None, dir=None, compress_type=0, model=1, force=False):
534 """Create a PDB file of the molecule with the moving domains shifted to the average position.
535
536 @keyword format: The format for outputting the geometric representation. Currently only the 'PDB' format is supported.
537 @type format: str
538 @keyword file: The name of the file for the average molecule structure.
539 @type file: str
540 @keyword dir: The name of the directory to place the PDB file into.
541 @type dir: str
542 @keyword compress_type: The compression type. The integer values correspond to the compression type: 0, no compression; 1, Bzip2 compression; 2, Gzip compression.
543 @type compress_type: int
544 @keyword model: Only one model from an analysed ensemble can be used for the PDB representation of the Monte Carlo simulations, as these consists of one model per simulation.
545 @type model: int
546 @keyword force: Flag which if set to True will cause any pre-existing file to be overwritten.
547 @type force: bool
548 """
549
550
551 subsection(file=sys.stdout, text="Creating a PDB file with the moving domains shifted to the average position.")
552
553
554 if not hasattr(cdp, 'structure'):
555 warn(RelaxWarning("No structural data is present, cannot create the average position representation."))
556 return
557
558
559 titles = []
560 sims = []
561 file_root = []
562 models = []
563 structures = []
564
565
566 titles.append("real average position")
567 sims.append(False)
568 file_root.append(file)
569 models.append([None])
570
571
572 if hasattr(cdp, 'sim_number'):
573 titles.append("MC simulation representation")
574 sims.append(True)
575 file_root.append("%s_sim" % file)
576 models.append([i+1 for i in range(cdp.sim_number)])
577
578
579 structures.append(deepcopy(cdp.structure))
580 if hasattr(cdp, 'sim_number'):
581 structures.append(deepcopy(cdp.structure))
582
583
584 if hasattr(cdp, 'sim_number') and len(structures[-1].structural_data) > 1:
585
586 to_delete = []
587 for model_cont in structures[-1].model_loop():
588 if model_cont.num != model:
589 to_delete.append(model_cont.num)
590 to_delete.reverse()
591
592
593 for num in to_delete:
594 structures[-1].structural_data.delete_model(model_num=num)
595
596
597 for i in range(len(titles)):
598
599 subsubsection(file=sys.stdout, text="Creating the %s." % titles[i])
600
601
602 for j in range(len(models[i])):
603
604 if models[i][j] == 1:
605 structures[i].set_model(model_new=1)
606 elif models[i][j] != None:
607 structures[i].add_model(model=models[i][j])
608
609
610 average_position(structure=structures[i], models=models[i], sim=sims[i])
611
612
613 if format == 'PDB':
614 pdb_file = open_write_file(file_name=file_root[i]+'.pdb', dir=dir, compress_type=compress_type, force=force)
615 structures[i].write_pdb(file=pdb_file)
616 pdb_file.close()
617
618
619 -def create_geometric_rep(format='PDB', file=None, dir=None, compress_type=0, size=30.0, inc=36, force=False):
620 """Create a PDB file containing a geometric object representing the frame order dynamics.
621
622 @keyword format: The format for outputting the geometric representation. Currently only the 'PDB' format is supported.
623 @type format: str
624 @keyword file: The name of the file of the PDB representation of the frame order dynamics to create.
625 @type file: str
626 @keyword dir: The name of the directory to place the PDB file into.
627 @type dir: str
628 @keyword compress_type: The compression type. The integer values correspond to the compression type: 0, no compression; 1, Bzip2 compression; 2, Gzip compression.
629 @type compress_type: int
630 @keyword size: The size of the geometric object in Angstroms.
631 @type size: float
632 @keyword inc: The number of increments for the filling of the cone objects.
633 @type inc: int
634 @keyword force: Flag which if set to True will cause any pre-existing file to be overwritten.
635 @type force: bool
636 """
637
638
639 subsection(file=sys.stdout, text="Creating a PDB file containing a geometric object representing the frame order dynamics.")
640
641
642 check_parameters(escalate=2)
643
644
645 titles = []
646 structures = []
647 representation = []
648 sims = []
649 file_root = []
650
651
652 sym = True
653 if cdp.model in [MODEL_ROTOR, MODEL_FREE_ROTOR, MODEL_DOUBLE_ROTOR]:
654 sym = False
655
656
657 titles.append("Representation A")
658 structures.append(Internal())
659 if sym:
660 representation.append('A')
661 file_root.append("%s_A" % file)
662 else:
663 representation.append(None)
664 file_root.append(file)
665 sims.append(False)
666
667
668 if sym:
669 titles.append("Representation A")
670 structures.append(Internal())
671 representation.append('B')
672 file_root.append("%s_B" % file)
673 sims.append(False)
674
675
676 if hasattr(cdp, 'sim_number'):
677 titles.append("MC simulation representation A")
678 structures.append(Internal())
679 if sym:
680 representation.append('A')
681 file_root.append("%s_sim_A" % file)
682 else:
683 representation.append(None)
684 file_root.append("%s_sim" % file)
685 sims.append(True)
686
687
688 if hasattr(cdp, 'sim_number') and sym:
689 titles.append("MC simulation representation B")
690 structures.append(Internal())
691 representation.append('B')
692 file_root.append("%s_sim_B" % file)
693 sims.append(True)
694
695
696 for i in range(len(structures)):
697
698 subsubsection(file=sys.stdout, text="Creating the %s." % titles[i])
699
700
701 if sims[i]:
702 for sim_i in range(cdp.sim_number):
703 structures[i].add_model(model=sim_i+1)
704
705
706 add_pivots(structure=structures[i], sims=sims[i])
707
708
709 add_rotors(structure=structures[i], representation=representation[i], size=size, sims=sims[i])
710
711
712 add_axes(structure=structures[i], representation=representation[i], size=size, sims=sims[i])
713
714
715 if cdp.model not in [MODEL_ROTOR, MODEL_FREE_ROTOR, MODEL_DOUBLE_ROTOR]:
716 add_cones(structure=structures[i], representation=representation[i], size=size, inc=inc, sims=sims[i])
717
718
719 add_titles(structure=structures[i], representation=representation[i], displacement=size+10, sims=sims[i])
720
721
722 if format == 'PDB':
723 pdb_file = open_write_file(file_root[i]+'.pdb', dir, compress_type=compress_type, force=force)
724 structures[i].write_pdb(pdb_file)
725 pdb_file.close()
726
727
729 """Build a full 3D frame from the single axis.
730
731 @param axis: The Z-axis of the system.
732 @type axis: numpy rank-1, 3D array
733 @param frame: The empty frame to populate.
734 @type frame: numpy rank-2, 3D array
735 """
736
737
738 frame[:, 2] = axis
739
740
741 frame[0, 0] = 1.0
742
743
744 frame[:, 1] = cross(frame[:, 2], frame[:, 0])
745 frame[:, 1] /= norm(frame[:, 1])
746
747
748 frame[:, 0] = cross(frame[:, 1], frame[:, 2])
749 frame[:, 0] /= norm(frame[:, 0])
750
751
753 """Generate and return the full 3D axis system for the current frame order model.
754
755 @keyword sim_index: The optional MC simulation index to obtain the frame for a given simulation.
756 @type sim_index: None or int
757 @return: The full 3D axis system for the model.
758 @rtype: numpy rank-2, 3D float64 array
759 """
760
761
762 axis = zeros(3, float64)
763 frame = zeros((3, 3), float64)
764
765
766 pivot = generate_pivot(order=1, sim_index=sim_index, pdb_limit=True)
767
768
769 com = pipe_centre_of_mass(verbosity=0, missing_error=False)
770
771
772 if cdp.model in [MODEL_ROTOR, MODEL_FREE_ROTOR]:
773
774 if sim_index == None:
775 axis = create_rotor_axis_alpha(alpha=cdp.axis_alpha, pivot=pivot, point=com)
776 else:
777 axis = create_rotor_axis_alpha(alpha=cdp.axis_alpha_sim[sim_index], pivot=pivot, point=com)
778
779
780 frame_from_axis(axis, frame)
781
782
783 elif cdp.model in MODEL_LIST_ISO_CONE:
784
785 if sim_index == None:
786 axis = create_rotor_axis_spherical(theta=cdp.axis_theta, phi=cdp.axis_phi)
787 else:
788 axis = create_rotor_axis_spherical(theta=cdp.axis_theta_sim[sim_index], phi=cdp.axis_phi_sim[sim_index])
789
790
791 frame_from_axis(axis, frame)
792
793
794 elif cdp.model in MODEL_LIST_PSEUDO_ELLIPSE + MODEL_LIST_DOUBLE:
795
796 if sim_index == None:
797 euler_to_R_zyz(cdp.eigen_alpha, cdp.eigen_beta, cdp.eigen_gamma, frame)
798 else:
799 euler_to_R_zyz(cdp.eigen_alpha_sim[sim_index], cdp.eigen_beta_sim[sim_index], cdp.eigen_gamma_sim[sim_index], frame)
800
801
802 else:
803 raise RelaxFault
804
805
806 return frame
807