1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23 """Module containing the internal relax structural object."""
24
25
26 from numpy import array, dot, float64, linalg, zeros
27 import os
28 from os import F_OK, access
29 from re import search
30 from string import digits
31 from warnings import warn
32
33
34 from data.relax_xml import fill_object_contents, xml_to_object
35 from generic_fns import pipes, relax_re
36 from generic_fns.mol_res_spin import spin_loop
37 from generic_fns.mol_res_spin import Selection
38 from generic_fns.structure.api_base import Base_struct_API, ModelList, Displacements
39 from relax_errors import RelaxError, RelaxNoneIntError, RelaxNoPdbError
40 from relax_io import file_root, open_read_file
41 from relax_warnings import RelaxWarning
42
43
44
46 """The internal relax structural data object.
47
48 The structural data object for this class is a container possessing a number of different arrays
49 corresponding to different structural information. These objects are described in the
50 structural container docstring.
51 """
52
53
54 id = 'internal'
55
56
58 """Find the atom named attached_atom directly bonded to the atom located at the index.
59
60 @param attached_atom: The name of the attached atom to return.
61 @type attached_atom: str
62 @param index: The index of the atom which the attached atom is attached to.
63 @type index: int
64 @param mol: The molecule container.
65 @type mol: MolContainer instance
66 @return: A tuple of information about the bonded atom.
67 @rtype: tuple consisting of the atom number (int), atom name (str), element
68 name (str), and atomic position (Numeric array of len 3)
69 """
70
71
72 bonded_found = False
73
74
75 if not mol.bonded[index]:
76
77 if not hasattr(mol, 'type'):
78 self._mol_type(mol)
79
80
81 if mol.type == 'protein':
82 self._protein_connect(mol)
83
84
85 else:
86 self._find_bonded_atoms(index, mol, radius=2)
87
88
89 matching_list = []
90 for bonded_index in mol.bonded[index]:
91 if relax_re.search(mol.atom_name[bonded_index], attached_atom):
92 matching_list.append(bonded_index)
93 num_attached = len(matching_list)
94
95
96 if num_attached > 1:
97
98 matching_names = []
99 for i in matching_list:
100 matching_names.append(mol.atom_name[i])
101
102
103 return None, None, None, None, None, 'More than one attached atom found: ' + repr(matching_names)
104
105
106 if num_attached == 0:
107 if relax_re.search('@*', attached_atom):
108 matching_list = []
109 bonded_num=[]
110 bonded_name=[]
111 element=[]
112 pos=[]
113 for spin, mol_name, res_num, res_name in spin_loop(selection=attached_atom, full_info=True):
114 bonded_num.append(spin.num)
115 bonded_name.append(spin.name)
116 element.append(spin.element)
117 pos.append(spin.pos)
118 if len(bonded_num) == 1:
119 return bonded_num[0], bonded_name[0], element[0], pos[0], attached_atom, None
120 elif len(bonded_num) > 1:
121
122 return None, None, None, None, None, 'More than one attached atom found: ' + repr(matching_names)
123 elif len(bonded_num) > 1:
124
125 return None, None, None, None, None, "No attached atom could be found"
126 else:
127 return None, None, None, None, None, "No attached atom could be found"
128
129
130 index = matching_list[0]
131 bonded_num = mol.atom_num[index]
132 bonded_name = mol.atom_name[index]
133 element = mol.element[index]
134 pos = [mol.x[index], mol.y[index], mol.z[index]]
135 attached_name = mol.atom_name[index]
136
137
138 return bonded_num, bonded_name, element, pos, attached_name, None
139
140
142 """Find all atoms within a sphere and say that they are attached to the central atom.
143
144 The found atoms will be added to the 'bonded' data structure.
145
146
147 @param index: The index of the central atom.
148 @type index: int
149 @param mol: The molecule container.
150 @type mol: MolContainer instance
151 """
152
153
154 centre = array([mol.x[index], mol.y[index], mol.z[index]], float64)
155
156
157 dist_list = []
158 connect_list = {}
159 element_list = {}
160 for i in range(len(mol.atom_num)):
161
162 if mol.element[index] == 'H' and mol.element[i] == 'H':
163 continue
164
165
166 pos = array([mol.x[i], mol.y[i], mol.z[i]], float64)
167
168
169 dist = linalg.norm(centre-pos)
170
171
172 if dist < radius:
173
174 dist_list.append(dist)
175
176
177 connect_list[dist] = i
178
179
180 element_list[dist] = mol.element[i]
181
182
183 max_conn = 1000
184 if mol.element[index] == 'H':
185 max_conn = 1
186 elif mol.element[index] == 'O':
187 max_conn = 2
188 elif mol.element[index] == 'N':
189 max_conn = 3
190 elif mol.element[index] == 'C':
191 max_conn = 4
192
193
194 dist_list.sort()
195
196
197 for i in range(min(max_conn, len(dist_list))):
198 mol.atom_connect(index, connect_list[dist_list[i]])
199
200
202 """Return the chemical name corresponding to the given residue ID.
203
204 The following names are currently returned::
205 ________________________________________________
206 | | |
207 | hetID | Chemical name |
208 |________|_____________________________________|
209 | | |
210 | TNS | Tensor |
211 | COM | Centre of mass |
212 | AXS | Tensor axes |
213 | SIM | Monte Carlo simulation tensor axes |
214 | PIV | Pivot point |
215 | CON | Cone object |
216 | AVE | Average vector |
217 |________|_____________________________________|
218
219 For any other residues, no description is returned.
220
221 @param hetID: The residue ID.
222 @type hetID: str
223 @return: The chemical name.
224 @rtype: str or None
225 """
226
227
228 if hetID == 'TNS':
229 return 'Tensor'
230
231
232 if hetID == 'COM':
233 return 'Centre of mass'
234
235
236 if hetID == 'AXS':
237 return 'Tensor axes'
238
239
240 if hetID == 'SIM':
241 return 'Monte Carlo simulation tensor axes'
242
243
244 if hetID == 'PIV':
245 return 'Pivot point'
246
247
248 if hetID == 'CON':
249 return 'Cone'
250
251
252 if hetID == 'AVE':
253 return 'Average vector'
254
255
257 """Generator function for looping over the models in the PDB file.
258
259 @param file_path: The full path of the PDB file.
260 @type file_path: str
261 @return: The model number and all the records for that model.
262 @rtype: tuple of int and array of str
263 """
264
265
266 file = open_read_file(file_path)
267 lines = file.readlines()
268 file.close()
269
270
271 if lines == []:
272 raise RelaxError("The PDB file is empty.")
273
274
275 model = None
276 records = []
277
278
279 for i in range(len(lines)):
280
281 if search('^MODEL', lines[i]):
282 try:
283 model = int(lines[i].split()[1])
284 except:
285 raise RelaxError("The MODEL record " + repr(lines[i]) + " is corrupt, cannot read the PDB file.")
286
287
288 if not (search('^ATOM', lines[i]) or search('^HETATM', lines[i])) and not len(records):
289 continue
290
291
292 if search('^ENDMDL', lines[i]):
293
294 yield model, records
295
296
297 records = []
298
299
300 continue
301
302
303 records.append(lines[i])
304
305
306 if len(records):
307 yield model, records
308
309
311 """Generator function for looping over the models in the XYZ file.
312
313 @param file_path: The full path of the XYZ file.
314 @type file_path: str
315 @return: The model number and all the records for that model.
316 @rtype: tuple of int and array of str
317 """
318
319
320 file = open_read_file(file_path)
321 lines = file.readlines()
322 file.close()
323
324
325 if lines == []:
326 raise RelaxError("The XYZ file is empty.")
327
328
329 total_atom = 0
330 model = 0
331 records = []
332
333
334 for i in range(len(lines)):
335 num=0
336 word = lines[i].split()
337
338 if (i==0) and (len(word)==1):
339 try:
340 total_atom = int(word[0])
341 num = 1
342 except:
343 raise RelaxError("The MODEL record " + repr(lines[i]) + " is corrupt, cannot read the XYZ file.")
344
345
346 if (len(records) == total_atom):
347
348 yield records
349
350
351 records = []
352
353
354 if (len(word) != 4):
355 continue
356
357
358 records.append(lines[i])
359
360
361 if len(records):
362 yield records
363
364
366 """Generator function for looping over the molecules in the PDB records of a model.
367
368 @param records: The list of PDB records for the model, or if no models exist the entire
369 PDB file.
370 @type records: list of str
371 @return: The molecule number and all the records for that molecule.
372 @rtype: tuple of int and list of str
373 """
374
375
376 if records == []:
377 raise RelaxError("There are no PDB records for this model.")
378
379
380 mol_num = 1
381 mol_records = []
382 end = False
383
384
385 for i in range(len(records)):
386
387 if search('^END', records[i]):
388 break
389
390
391 if search('^ENDMDL', records[i]):
392 end = True
393
394
395 elif i < len(records)-1 and search('^TER', records[i]) and not search('^HETATM', records[i+1]):
396 end = True
397
398
399 elif i < len(records)-1 and search('^HETATM', records[i]) and search('^ATOM', records[i+1]):
400 end = True
401
402
403 if end:
404
405 yield mol_num, mol_records
406
407
408 mol_records = []
409
410
411 mol_num = mol_num + 1
412
413
414 end = False
415
416
417 continue
418
419
420 mol_records.append(records[i])
421
422
423 if len(mol_records):
424 yield mol_num, mol_records
425
426
428 """Check the validity of the data arrays in the given structure object.
429
430 @param struct: The structural object.
431 @type struct: Structure_container instance
432 """
433
434
435 num = len(struct.atom_name)
436
437
438 if len(struct.bonded) != num and len(struct.chain_id) != num and len(struct.element) != num and len(struct.pdb_record) != num and len(struct.res_name) != num and len(struct.res_num) != num and len(struct.seg_id) != num and len(struct.x) != num and len(struct.y) != num and len(struct.z) != num:
439 raise RelaxError("The structural data is invalid.")
440
441
443 """Determine the type of molecule.
444
445 @param mol: The molecule data container.
446 @type mol: MolContainer instance
447 """
448
449
450 aa = ['ALA', 'ARG', 'ASN', 'ASP', 'CYS', 'GLU', 'GLN', 'GLY', 'HIS', 'ILE', 'LEU', 'LYS', 'MET', 'PHE', 'PRO', 'SER', 'THR', 'TRP', 'TYR', 'VAL']
451
452
453 mol.type = 'other'
454
455
456 for res in mol.res_name:
457
458 if res in aa:
459
460 mol.type = 'protein'
461 return
462
463
465 """Set up the connectivities for the protein.
466
467 @param mol: The molecule data container.
468 @type mol: MolContainer instance
469 """
470
471
472 curr_res_num = None
473 res_atoms = []
474
475
476 for i in range(len(mol.atom_num)):
477
478 if mol.res_num[i] != curr_res_num:
479
480 if len(res_atoms):
481 self._protein_intra_connect(mol, res_atoms)
482
483
484 curr_res_num = mol.res_num[i]
485
486
487 res_atoms = []
488
489
490 res_atoms.append(i)
491
492
493 if i == len(mol.atom_num) - 1 and len(res_atoms):
494 self._protein_intra_connect(mol, res_atoms)
495
496
498 """Set up the connectivities for the protein.
499
500 @param mol: The molecule data container.
501 @type mol: MolContainer instance
502 @param res_atoms: The list of atom indices corresponding to the residue.
503 @type res_atoms: list of int
504 """
505
506
507 indices = {
508 'N': None,
509 'C': None,
510 'O': None,
511 'CA': None,
512 'HN': None,
513 'H': None,
514 'HA': None
515 }
516
517
518 for index in res_atoms:
519 if mol.atom_name[index] in indices:
520 indices[mol.atom_name[index]] = index
521
522
523 pairs = [
524 ['N', 'HN'],
525 ['N', 'H'],
526 ['N', 'CA'],
527 ['CA', 'HA'],
528 ['CA', 'C'],
529 ['C', 'O']
530 ]
531
532
533 for pair in pairs:
534 if indices[pair[0]] != None and indices[pair[1]] != None:
535 mol.atom_connect(indices[pair[0]], indices[pair[1]])
536
537
539 """Convert the data into a format for writing to file.
540
541 @param data: The data to convert to the required format.
542 @type data: anything
543 @keyword format: The format to convert to. This can be 'str', 'float', or 'int'.
544 @type format: str
545 @return: The converted version of the data.
546 @rtype: str
547 """
548
549
550 if format == 'str':
551
552 if data == None:
553 data = ''
554
555
556 if not isinstance(data, str):
557 data = repr(data)
558
559
560 if format == 'float':
561
562 if data == None:
563 data = 0.0
564
565
566 if not isinstance(data, float):
567 data = float(data)
568
569
570 return data
571
572
573 - def add_atom(self, mol_name=None, atom_name=None, res_name=None, res_num=None, pos=[None, None, None], element=None, atom_num=None, chain_id=None, segment_id=None, pdb_record=None):
574 """Add a new atom to the structural data object.
575
576 @keyword mol_name: The name of the molecule.
577 @type mol_name: str
578 @keyword atom_name: The atom name, e.g. 'H1'.
579 @type atom_name: str or None
580 @keyword res_name: The residue name.
581 @type res_name: str or None
582 @keyword res_num: The residue number.
583 @type res_num: int or None
584 @keyword pos: The position vector of coordinates.
585 @type pos: list (length = 3)
586 @keyword element: The element symbol.
587 @type element: str or None
588 @keyword atom_num: The atom number.
589 @type atom_num: int or None
590 @keyword chain_id: The chain identifier.
591 @type chain_id: str or None
592 @keyword segment_id: The segment identifier.
593 @type segment_id: str or None
594 @keyword pdb_record: The optional PDB record name, e.g. 'ATOM' or 'HETATM'.
595 @type pdb_record: str or None
596 """
597
598
599 pipes.test()
600
601
602 if len(self.structural_data) == 0:
603 self.add_model()
604
605
606 for model in self.structural_data:
607
608 mol = self.get_molecule(mol_name, model=model.num)
609
610
611 if mol == None:
612 self.add_molecule(name=mol_name)
613 mol = self.get_molecule(mol_name, model=model.num)
614
615
616 mol.atom_add(atom_name=atom_name, res_name=res_name, res_num=res_num, pos=pos, element=element, atom_num=atom_num, chain_id=chain_id, segment_id=segment_id, pdb_record=pdb_record)
617
618
619 - def add_model(self, model=None, coords_from=None):
620 """Add a new model to the store.
621
622 The new model will be constructured with the structural information from the other models currently present. The coords_from argument allows the atomic positions to be taken from a certain model. If this argument is not set, then the atomic positions from the first model will be used.
623
624 @keyword model: The number of the model to create.
625 @type model: int or None
626 @keyword coords_from: The model number to take the coordinates from.
627 @type coords_from: int or None
628 @return: The model container.
629 @rtype: ModelContainer instance
630 """
631
632
633 if model != None:
634 for i in range(len(self.structural_data)):
635 if model == self.structural_data[i].num:
636 raise RelaxError("The model '%s' already exists." % model)
637
638
639 self.structural_data.add_item(model_num=model)
640
641
642 if coords_from == None:
643 coords_from = self.structural_data[0].num
644
645
646 for mol_name, res_num, res_name, atom_num, atom_name, element, pos in self.atom_loop(model_num=coords_from, mol_name_flag=True, res_num_flag=True, res_name_flag=True, atom_num_flag=True, atom_name_flag=True, element_flag=True, pos_flag=True):
647
648 self.add_atom(self, mol_name=mol_name, atom_name=atom_name, res_name=res_name, res_num=res_num, pos=pos, element=element, atom_num=atom_num)
649
650
651 return self.structural_data[-1]
652
653
655 """Add a new molecule to the store.
656
657 @keyword name: The molecule identifier string.
658 @type name: str
659 """
660
661
662 if len(self.structural_data) == 0:
663 self.add_model()
664
665
666 for i in range(len(self.structural_data)):
667
668 self.structural_data[i].mol.add_item(mol_name=name, mol_cont=MolContainer())
669
670
671 - def are_bonded(self, atom_id1=None, atom_id2=None):
672 """Determine if two atoms are directly bonded to each other.
673
674 @keyword atom_id1: The molecule, residue, and atom identifier string of the first atom.
675 @type atom_id1: str
676 @keyword atom_id2: The molecule, residue, and atom identifier string of the second atom.
677 @type atom_id2: str
678 @return: True if the atoms are directly bonded.
679 @rtype: bool
680 """
681
682
683 sel_obj1 = Selection(atom_id1)
684 sel_obj2 = Selection(atom_id2)
685
686
687 for mol in self.structural_data[0].mol:
688 for i in range(len(mol.atom_num)):
689 if not len(mol.bonded[i]):
690 self._find_bonded_atoms(i, mol, radius=2)
691
692
693 for mol in self.structural_data[0].mol:
694
695 if not sel_obj1.contains_mol(mol.mol_name):
696 continue
697 if not sel_obj2.contains_mol(mol.mol_name):
698 continue
699
700
701 index1 = None
702 for i in range(len(mol.atom_num)):
703
704 if sel_obj1.contains_spin(mol.atom_num[i], mol.atom_name[i], mol.res_num[i], mol.res_name[i], mol.mol_name):
705 index1 = i
706 break
707
708
709 index2 = None
710 for i in range(len(mol.atom_num)):
711
712 if sel_obj2.contains_spin(mol.atom_num[i], mol.atom_name[i], mol.res_num[i], mol.res_name[i], mol.mol_name):
713 index2 = i
714 break
715
716
717 if index1 < len(mol.bonded):
718 if index2 in mol.bonded[index1]:
719 return True
720 else:
721 return False
722
723
724 - def atom_loop(self, atom_id=None, str_id=None, model_num=None, model_num_flag=False, mol_name_flag=False, res_num_flag=False, res_name_flag=False, atom_num_flag=False, atom_name_flag=False, element_flag=False, pos_flag=False, ave=False):
725 """Generator function for looping over all atoms in the internal relax structural object.
726
727 @keyword atom_id: The molecule, residue, and atom identifier string. Only atoms matching this selection will be yielded.
728 @type atom_id: str
729 @keyword str_id: The structure identifier. This can be the file name, model number, or structure number. If None, then all structures will be looped over.
730 @type str_id: str, int, or None
731 @keyword model_num: Only loop over a specific model.
732 @type model_num: int or None
733 @keyword model_num_flag: A flag which if True will cause the model number to be yielded.
734 @type model_num_flag: bool
735 @keyword mol_name_flag: A flag which if True will cause the molecule name to be yielded.
736 @type mol_name_flag: bool
737 @keyword res_num_flag: A flag which if True will cause the residue number to be yielded.
738 @type res_num_flag: bool
739 @keyword res_name_flag: A flag which if True will cause the residue name to be yielded.
740 @type res_name_flag: bool
741 @keyword atom_num_flag: A flag which if True will cause the atom number to be yielded.
742 @type atom_num_flag: bool
743 @keyword atom_name_flag: A flag which if True will cause the atom name to be yielded.
744 @type atom_name_flag: bool
745 @keyword element_flag: A flag which if True will cause the element name to be yielded.
746 @type element_flag: bool
747 @keyword pos_flag: A flag which if True will cause the atomic position to be yielded.
748 @type pos_flag: bool
749 @keyword ave: A flag which if True will result in this method returning the average atom properties across all loaded structures.
750 @type ave: bool
751 @return: A tuple of atomic information, as described in the docstring.
752 @rtype: tuple consisting of optional molecule name (str), residue number (int), residue name (str), atom number (int), atom name(str), element name (str), and atomic position (array of len 3).
753 """
754
755
756 if not len(self.structural_data):
757 raise RelaxNoPdbError
758
759
760 sel_obj = None
761 if atom_id:
762 sel_obj = Selection(atom_id)
763
764
765 for model in self.model_loop(model_num):
766
767 for mol_index in range(len(model.mol)):
768 mol = model.mol[mol_index]
769
770
771 if sel_obj and not sel_obj.contains_mol(mol.mol_name):
772 continue
773
774
775 for i in range(len(mol.atom_name)):
776
777 if sel_obj and not sel_obj.contains_spin(mol.atom_num[i], mol.atom_name[i], mol.res_num[i], mol.res_name[i], mol.mol_name):
778 continue
779
780
781 res_num = mol.res_num[i]
782 res_name = mol.res_name[i]
783 atom_num = mol.atom_num[i]
784 atom_name = mol.atom_name[i]
785 element = mol.element[i]
786 pos = zeros(3, float64)
787
788
789 if ave:
790
791 for model in self.model_loop():
792
793 mol = model.mol[mol_index]
794
795
796 if mol.atom_num[i] != atom_num:
797 raise RelaxError("The loaded structures do not contain the same atoms. The average structural properties can not be calculated.")
798
799
800 pos = pos + array([mol.x[i], mol.y[i], mol.z[i]], float64)
801
802
803 pos = pos / len(self.structural_data)
804 else:
805 pos = array([mol.x[i], mol.y[i], mol.z[i]], float64)
806
807
808 mol_name = mol.mol_name
809
810
811 atomic_tuple = ()
812 if model_num_flag:
813 if ave:
814 atomic_tuple = atomic_tuple + (None,)
815 else:
816 atomic_tuple = atomic_tuple + (model.num,)
817 if mol_name_flag:
818 atomic_tuple = atomic_tuple + (mol_name,)
819 if res_num_flag:
820 atomic_tuple = atomic_tuple + (res_num,)
821 if res_name_flag:
822 atomic_tuple = atomic_tuple + (res_name,)
823 if atom_num_flag:
824 atomic_tuple = atomic_tuple + (atom_num,)
825 if atom_name_flag:
826 atomic_tuple = atomic_tuple + (atom_name,)
827 if element_flag:
828 atomic_tuple = atomic_tuple + (element,)
829 if pos_flag:
830 atomic_tuple = atomic_tuple + (pos,)
831
832
833 yield atomic_tuple
834
835
836 if ave:
837 break
838
839
840 - def bond_vectors(self, attached_atom=None, model_num=None, mol_name=None, res_num=None, res_name=None, spin_num=None, spin_name=None, return_name=False, return_warnings=False):
841 """Find the bond vector between the atoms of 'attached_atom' and 'atom_id'.
842
843 @keyword attached_atom: The name of the bonded atom.
844 @type attached_atom: str
845 @keyword model_num: The model of which to return the vectors from. If not supplied and multiple models exist, then vectors from all models will be returned.
846 @type model_num: None or int
847 @keyword mol_name: The name of the molecule that attached_atom belongs to.
848 @type mol_name: str
849 @keyword res_num: The number of the residue that attached_atom belongs to.
850 @type res_num: str
851 @keyword res_name: The name of the residue that attached_atom belongs to.
852 @type res_name: str
853 @keyword spin_num: The number of the spin that attached_atom is attached to.
854 @type spin_num: str
855 @keyword spin_name: The name of the spin that attached_atom is attached to.
856 @type spin_name: str
857 @keyword return_name: A flag which if True will cause the name of the attached atom to be returned together with the bond vectors.
858 @type return_name: bool
859 @keyword return_warnings: A flag which if True will cause warning messages to be returned.
860 @type return_warnings: bool
861 @return: The list of bond vectors for each model.
862 @rtype: list of numpy arrays (or a tuple if return_name or return_warnings are set)
863 """
864
865
866 vectors = []
867 attached_name = None
868 warnings = None
869
870
871 for model in self.model_loop(model_num):
872
873 for mol in model.mol:
874
875 if mol_name and mol_name != mol.mol_name:
876 continue
877
878
879 index = None
880 for i in range(len(mol.atom_name)):
881
882 if (res_num != None and mol.res_num[i] != res_num) or (res_name != None and mol.res_name[i] != res_name):
883 continue
884
885
886 if (spin_num != None and mol.atom_num[i] != spin_num) or (spin_name != None and mol.atom_name[i] != spin_name):
887 continue
888
889
890 index = i
891 break
892
893
894 if index != None:
895
896 bonded_num, bonded_name, element, pos, attached_name, warnings = self._bonded_atom(attached_atom, index, mol)
897
898
899 if (bonded_num, bonded_name, element) == (None, None, None):
900 continue
901
902
903 vector = array(pos, float64) - array([mol.x[index], mol.y[index], mol.z[index]], float64)
904
905
906 vectors.append(vector)
907
908
909 else:
910 warnings = "Cannot find the atom in the structure"
911
912
913 data = (vectors,)
914 if return_name:
915 data = data + (attached_name,)
916 if return_warnings:
917 data = data + (warnings,)
918
919
920 return data
921
922
923 - def connect_atom(self, mol_name=None, index1=None, index2=None):
924 """Connect two atoms in the structural data object.
925
926 @keyword mol_name: The name of the molecule.
927 @type mol_name: str
928 @keyword index1: The global index of the first atom.
929 @type index1: str
930 @keyword index2: The global index of the first atom.
931 @type index2: str
932 """
933
934
935 pipes.test()
936
937
938 if self.get_molecule(mol_name) == None:
939 self.add_molecule(name=mol_name)
940
941
942 for model in self.structural_data:
943
944 mol = self.get_molecule(mol_name)
945
946
947 mol.atom_connect(index1=index1, index2=index2)
948
949
951 """Delete all the structural information."""
952
953
954 print("Deleting the following structural data:\n")
955 print(self.structural_data)
956
957
958 del self.structural_data
959
960
961 self.structural_data = ModelList()
962
963
965 """Return the molecule.
966
967 Only one model can be specified.
968
969
970 @param molecule: The molecule name.
971 @type molecule: int or None
972 @keyword model: The model number.
973 @type model: int or None
974 @raises RelaxError: If the model is not specified and there is more than one model loaded.
975 @return: The MolContainer corresponding to the molecule name and model number.
976 @rtype: MolContainer instance or None
977 """
978
979
980 if model == None and self.num_models() > 1:
981 raise RelaxError("The target molecule cannot be determined as there are %s models already present." % self.num_models())
982
983
984 if not isinstance(model, int) and not model == None:
985 raise RelaxNoneIntError
986
987
988 if not len(self.structural_data):
989 return
990
991
992 for model_cont in self.model_loop(model):
993
994 for mol in model_cont.mol:
995
996 if mol.mol_name == molecule:
997 return mol
998
999
1000 - def load_pdb(self, file_path, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, verbosity=False):
1001 """Method for loading structures from a PDB file.
1002
1003 @param file_path: The full path of the PDB file.
1004 @type file_path: str
1005 @keyword read_mol: The molecule(s) to read from the file, independent of model. The
1006 molecules are determined differently by the different parsers, but
1007 are numbered consecutively from 1. If set to None, then all
1008 molecules will be loaded.
1009 @type read_mol: None, int, or list of int
1010 @keyword set_mol_name: Set the names of the molecules which are loaded. If set to None,
1011 then the molecules will be automatically labelled based on the file
1012 name or other information.
1013 @type set_mol_name: None, str, or list of str
1014 @keyword read_model: The PDB model to extract from the file. If set to None, then all
1015 models will be loaded.
1016 @type read_model: None, int, or list of int
1017 @keyword set_model_num: Set the model number of the loaded molecule. If set to None, then
1018 the PDB model numbers will be preserved, if they exist.
1019 @type set_model_num: None, int, or list of int
1020 @keyword verbosity: A flag which if True will cause messages to be printed.
1021 @type verbosity: bool
1022 @return: The status of the loading of the PDB file.
1023 @rtype: bool
1024 """
1025
1026
1027 if verbosity:
1028 print("\nInternal relax PDB parser.")
1029
1030
1031 if not access(file_path, F_OK):
1032
1033 return False
1034
1035
1036 path, file = os.path.split(file_path)
1037
1038
1039 if read_mol and not isinstance(read_mol, list):
1040 read_mol = [read_mol]
1041 if set_mol_name and not isinstance(set_mol_name, list):
1042 set_mol_name = [set_mol_name]
1043 if read_model and not isinstance(read_model, list):
1044 read_model = [read_model]
1045 if set_model_num and not isinstance(set_model_num, list):
1046 set_model_num = [set_model_num]
1047
1048
1049 model_index = 0
1050 orig_model_num = []
1051 mol_conts = []
1052 for model_num, model_records in self._parse_models_pdb(file_path):
1053
1054 if read_model and model_num not in read_model:
1055 continue
1056
1057
1058 orig_model_num.append(model_num)
1059
1060
1061 mol_conts.append([])
1062 mol_index = 0
1063 orig_mol_num = []
1064 new_mol_name = []
1065 for mol_num, mol_records in self._parse_mols(model_records):
1066
1067 if read_mol and mol_num not in read_mol:
1068 continue
1069
1070
1071 if set_mol_name:
1072 new_mol_name.append(set_mol_name[mol_index])
1073 else:
1074
1075 num_struct = 0
1076 for model in self.structural_data:
1077 if not set_model_num or (model_index <= len(set_model_num) and set_model_num[model_index] == model.num):
1078 num_struct = len(model.mol)
1079
1080
1081 new_mol_name.append(file_root(file) + '_mol' + repr(mol_num+num_struct))
1082
1083
1084 orig_mol_num.append(mol_num)
1085
1086
1087 mol = MolContainer()
1088
1089
1090 mol.fill_object_from_pdb(mol_records)
1091
1092
1093 mol_conts[model_index].append(mol)
1094
1095
1096 mol_index = mol_index + 1
1097
1098
1099 model_index = model_index + 1
1100
1101
1102 if not len(mol_conts):
1103 warn(RelaxWarning("No structural data could be read from the file '%s'." % file_path))
1104 return False
1105
1106
1107 self.pack_structs(mol_conts, orig_model_num=orig_model_num, set_model_num=set_model_num, orig_mol_num=orig_mol_num, set_mol_name=new_mol_name, file_name=file, file_path=path)
1108
1109
1110 return True
1111
1112
1113 - def load_xyz(self, file_path, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, verbosity=False):
1114 """Method for loading structures from a XYZ file.
1115
1116 @param file_path: The full path of the XYZ file.
1117 @type file_path: str
1118 @keyword read_mol: The molecule(s) to read from the file, independent of model. The
1119 molecules are determined differently by the different parsers, but
1120 are numbered consecutively from 1. If set to None, then all
1121 molecules will be loaded.
1122 @type read_mol: None, int, or list of int
1123 @keyword set_mol_name: Set the names of the molecules which are loaded. If set to None,
1124 then the molecules will be automatically labelled based on the file
1125 name or other information.
1126 @type set_mol_name: None, str, or list of str
1127 @keyword read_model: The XYZ model to extract from the file. If set to None, then all
1128 models will be loaded.
1129 @type read_model: None, int, or list of int
1130 @keyword set_model_num: Set the model number of the loaded molecule. If set to None, then
1131 the XYZ model numbers will be preserved, if they exist.
1132 @type set_model_num: None, int, or list of int
1133 @keyword verbosity: A flag which if True will cause messages to be printed.
1134 @type verbosity: bool
1135 @return: The status of the loading of the XYZ file.
1136 @rtype: bool
1137 """
1138
1139
1140 if verbosity:
1141 print("\nInternal relax XYZ parser.")
1142
1143
1144 if not access(file_path, F_OK):
1145
1146 return False
1147
1148
1149 path, file = os.path.split(file_path)
1150
1151
1152 if read_mol and not isinstance(read_mol, list):
1153 read_mol = [read_mol]
1154 if set_mol_name and not isinstance(set_mol_name, list):
1155 set_mol_name = [set_mol_name]
1156 if read_model and not isinstance(read_model, list):
1157 read_model = [read_model]
1158 if set_model_num and not isinstance(set_model_num, list):
1159 set_model_num = [set_model_num]
1160
1161
1162 mol_index=0
1163 model_index = 0
1164 xyz_model_increment = 0
1165 orig_model_num = []
1166 mol_conts = []
1167 orig_mol_num = []
1168 new_mol_name = []
1169 for model_records in self._parse_models_xyz(file_path):
1170
1171 xyz_model_increment = xyz_model_increment +1
1172
1173
1174 if read_model and xyz_model_increment not in read_model:
1175 continue
1176
1177
1178 orig_model_num.append(model_index)
1179
1180
1181 if read_mol and mol_index not in read_mol:
1182 continue
1183
1184
1185 if set_mol_name:
1186 new_mol_name.append(set_mol_name[mol_index])
1187 else:
1188 if mol_index==0:
1189
1190 new_mol_name.append(file_root(file) + '_mol' + repr(mol_index+1))
1191
1192
1193 orig_mol_num.append(mol_index)
1194
1195
1196 mol = MolContainer()
1197
1198
1199 mol.fill_object_from_xyz(model_records)
1200
1201
1202 mol_conts.append([])
1203 mol_conts[model_index].append(mol)
1204
1205
1206 mol_index = mol_index + 1
1207
1208
1209 model_index = model_index + 1
1210
1211 orig_mol_num=[0]
1212
1213 self.pack_structs(mol_conts, orig_model_num=orig_model_num, set_model_num=set_model_num, orig_mol_num=orig_mol_num, set_mol_name=new_mol_name, file_name=file, file_path=path)
1214
1215
1216 return True
1217
1218
1219 - def rotate(self, R=None, origin=None, model=None, atom_id=None):
1220 """Rotate the structural information about the given origin.
1221
1222 @keyword R: The forwards rotation matrix.
1223 @type R: numpy 3D, rank-2 array
1224 @keyword origin: The origin of the rotation.
1225 @type origin: numpy 3D, rank-1 array
1226 @keyword model: The model to rotate. If None, all models will be rotated.
1227 @type model: int
1228 @keyword atom_id: The molecule, residue, and atom identifier string. Only atoms matching this selection will be used.
1229 @type atom_id: str or None
1230 """
1231
1232
1233 sel_obj = None
1234 if atom_id:
1235 sel_obj = Selection(atom_id)
1236
1237
1238 for model_cont in self.model_loop(model):
1239
1240 for mol in model_cont.mol:
1241
1242 if sel_obj and not sel_obj.contains_mol(mol.mol_name):
1243 continue
1244
1245
1246 for i in range(len(mol.atom_num)):
1247
1248 if sel_obj and not sel_obj.contains_spin(mol.atom_num[i], mol.atom_name[i], mol.res_num[i], mol.res_name[i], mol.mol_name):
1249 continue
1250
1251
1252 vect = array([mol.x[i], mol.y[i], mol.z[i]], float64) - origin
1253
1254
1255 rot_vect = dot(R, vect)
1256
1257
1258 pos = rot_vect + origin
1259 mol.x[i] = pos[0]
1260 mol.y[i] = pos[1]
1261 mol.z[i] = pos[2]
1262
1263
1264 - def translate(self, T=None, model=None, atom_id=None):
1265 """Displace the structural information by the given translation vector.
1266
1267 @keyword T: The translation vector.
1268 @type T: numpy 3D, rank-1 array
1269 @keyword model: The model to rotate. If None, all models will be rotated.
1270 @type model: int
1271 @keyword atom_id: The molecule, residue, and atom identifier string. Only atoms matching this selection will be used.
1272 @type atom_id: str or None
1273 """
1274
1275
1276 sel_obj = None
1277 if atom_id:
1278 sel_obj = Selection(atom_id)
1279
1280
1281 for model_cont in self.model_loop(model):
1282
1283 for mol in model_cont.mol:
1284
1285 if sel_obj and not sel_obj.contains_mol(mol.mol_name):
1286 continue
1287
1288
1289 for i in range(len(mol.atom_num)):
1290
1291 if sel_obj and not sel_obj.contains_spin(mol.atom_num[i], mol.atom_name[i], mol.res_num[i], mol.res_name[i], mol.mol_name):
1292 continue
1293
1294
1295 mol.x[i] = mol.x[i] + T[0]
1296 mol.y[i] = mol.y[i] + T[1]
1297 mol.z[i] = mol.z[i] + T[2]
1298
1299
1301 """Check that the models are consistent with each other.
1302
1303 This checks that the primary structure is identical between the models.
1304 """
1305
1306
1307 print("Validating models:")
1308
1309
1310 for i in range(len(self.structural_data)):
1311
1312 if len(self.structural_data[0].mol) != len(self.structural_data[i].mol):
1313 raise RelaxError("The number of molecules, %i, in model %i does not match the %i molecules of the first model." % (len(self.structural_data[i].mol), self.structural_data[i].num, len(self.structural_data[0].mol)))
1314
1315
1316 for j in range(len(self.structural_data[i].mol)):
1317
1318 mol = self.structural_data[i].mol[j]
1319 mol_ref = self.structural_data[0].mol[j]
1320
1321
1322 if mol.mol_name != mol_ref.mol_name:
1323 raise RelaxError("The molecule name '%s' of model %i does not match the name '%s' of the first model." % (mol.mol_name, self.structural_data[i].num, mol_ref.mol_name))
1324
1325
1326 for k in range(len(mol.atom_name)):
1327
1328 atom = "%-6s%5s %4s%1s%3s %1s%4s%1s %8s%8s%8s%6.2f%6.2f %4s%2s%2s" % ('ATOM', mol.atom_num[k], self._translate(mol.atom_name[k]), '', self._translate(mol.res_name[k]), self._translate(mol.chain_id[k]), self._translate(mol.res_num[k]), '', '#', '#', '#', 1.0, 0, self._translate(mol.seg_id[k]), self._translate(mol.element[k]), '')
1329 atom_ref = "%-6s%5s %4s%1s%3s %1s%4s%1s %8s%8s%8s%6.2f%6.2f %4s%2s%2s" % ('ATOM', mol_ref.atom_num[k], self._translate(mol_ref.atom_name[k]), '', self._translate(mol_ref.res_name[k]), self._translate(mol_ref.chain_id[k]), self._translate(mol_ref.res_num[k]), '', '#', '#', '#', 1.0, 0, self._translate(mol_ref.seg_id[k]), self._translate(mol_ref.element[k]), '')
1330
1331
1332 if atom != atom_ref:
1333 print(atom)
1334 print(atom_ref)
1335 raise RelaxError("The atoms of model %i do not match the first model." % self.structural_data[i].num)
1336
1337
1338 print("\tAll models are consistent")
1339
1340
1342 """Method for the creation of a PDB file from the structural data.
1343
1344 A number of PDB records including HET, HETNAM, FORMUL, HETATM, TER, CONECT, MASTER, and END
1345 are created. To create the non-standard residue records HET, HETNAM, and FORMUL, the data
1346 structure 'het_data' is created. It is an array of arrays where the first dimension
1347 corresponds to a different residue and the second dimension has the elements:
1348
1349 0. Residue number.
1350 1. Residue name.
1351 2. Chain ID.
1352 3. Total number of atoms in the residue.
1353 4. Number of H atoms in the residue.
1354 5. Number of C atoms in the residue.
1355
1356
1357 @param file: The PDB file object. This object must be writable.
1358 @type file: file object
1359 @keyword model_num: The model to place into the PDB file. If not supplied, then all
1360 models will be placed into the file.
1361 @type model_num: None or int
1362 """
1363
1364
1365 self.validate()
1366
1367
1368 num_hetatm = 0
1369 num_atom = 0
1370 num_ter = 0
1371 num_conect = 0
1372
1373
1374 print("\nCreating the PDB records\n")
1375
1376
1377 print("REMARK")
1378 file.write("REMARK 4 THIS FILE COMPLIES WITH FORMAT V. 3.1, 1-AUG-2007\n")
1379 file.write("REMARK 40 CREATED BY RELAX (HTTP://NMR-RELAX.COM)\n")
1380 num_remark = 2
1381
1382
1383 model_records = False
1384 for model in self.model_loop():
1385 if hasattr(model, 'num') and model.num != None:
1386 model_records = True
1387
1388
1389
1390
1391
1392
1393
1394 het_data = []
1395 het_data_coll = []
1396
1397
1398 index = 0
1399 for mol in self.structural_data[0].mol:
1400
1401 self._validate_data_arrays(mol)
1402
1403
1404 het_data.append([])
1405
1406
1407 for i in range(len(mol.atom_name)):
1408
1409 if mol.pdb_record[i] != 'HETATM' or mol.res_name[i] == None:
1410 continue
1411
1412
1413
1414 if not het_data[index] or not mol.res_num[i] == het_data[index][-1][0]:
1415 het_data[index].append([mol.res_num[i], mol.res_name[i], mol.chain_id[i], 0, []])
1416
1417
1418 if het_data[index][-1][2] == None:
1419 het_data[index][-1][2] = ''
1420
1421
1422 het_data[index][-1][3] = het_data[index][-1][3] + 1
1423
1424
1425 entry = False
1426 for j in range(len(het_data[index][-1][4])):
1427 if mol.element[i] == het_data[index][-1][4][j][0]:
1428 entry = True
1429
1430
1431 if not entry:
1432 het_data[index][-1][4].append([mol.element[i], 0])
1433
1434
1435 for j in range(len(het_data[index][-1][4])):
1436 if mol.element[i] == het_data[index][-1][4][j][0]:
1437 het_data[index][-1][4][j][1] = het_data[index][-1][4][j][1] + 1
1438
1439
1440 for i in range(len(het_data[index])):
1441
1442 found = False
1443 for j in range(len(het_data_coll)):
1444
1445 if het_data[index][i][0] == het_data_coll[j][0]:
1446
1447 found = True
1448
1449
1450 if het_data_coll[j][1] != het_data[index][i][1]:
1451 raise RelaxError("The " + repr(het_data[index][i][1]) + " residue name of hetrogen " + repr(het_data[index][i][0]) + " " + het_data[index][i][1] + " of structure " + repr(index) + " does not match the " + repr(het_data_coll[j][1]) + " name of the previous structures.")
1452
1453 elif het_data_coll[j][2] != het_data[index][i][2]:
1454 raise RelaxError("The hetrogen chain id " + repr(het_data[index][i][2]) + " does not match " + repr(het_data_coll[j][2]) + " of residue " + repr(het_data_coll[j][0]) + " " + het_data_coll[j][1] + " of the previous structures.")
1455
1456 elif het_data_coll[j][3] != het_data[index][i][3]:
1457 raise RelaxError("The " + repr(het_data[index][i][3]) + " atoms of hetrogen " + repr(het_data_coll[j][0]) + " " + het_data_coll[j][1] + " of structure " + repr(index) + " does not match the " + repr(het_data_coll[j][3]) + " of the previous structures.")
1458
1459 elif het_data_coll[j][4] != het_data[index][i][4]:
1460 raise RelaxError("The atom counts " + repr(het_data[index][i][4]) + " for the hetrogen residue " + repr(het_data_coll[j][0]) + " " + het_data_coll[j][1] + " of structure " + repr(index) + " do not match the counts " + repr(het_data_coll[j][4]) + " of the previous structures.")
1461
1462
1463 if not found:
1464 het_data_coll.append(het_data[index][i])
1465
1466
1467 index = index + 1
1468
1469
1470
1471
1472
1473
1474 print("HET")
1475
1476
1477 for het in het_data_coll:
1478 file.write("%-6s %3s %1s%4s%1s %5s %-40s\n" % ('HET', het[2], het[1], het[0], '', het[3], ''))
1479
1480
1481
1482
1483
1484
1485 print("HETNAM")
1486
1487
1488 residues = []
1489 for het in het_data_coll:
1490
1491 if het[1] in residues:
1492 continue
1493 else:
1494 residues.append(het[1])
1495
1496
1497 chemical_name = self._get_chemical_name(het[1])
1498 if not chemical_name:
1499 chemical_name = 'Unknown'
1500
1501
1502 file.write("%-6s %2s %3s %-55s\n" % ('HETNAM', '', het[1], chemical_name))
1503
1504
1505
1506
1507
1508
1509 print("FORMUL")
1510
1511
1512 residues = []
1513 for het in het_data_coll:
1514
1515 if het[1] in residues:
1516 continue
1517 else:
1518 residues.append(het[1])
1519
1520
1521 formula = ''
1522
1523
1524 for atom_count in het[4]:
1525 formula = formula + atom_count[0] + repr(atom_count[1])
1526
1527
1528 file.write("%-6s %2s %3s %2s%1s%-51s\n" % ('FORMUL', het[0], het[1], '', '', formula))
1529
1530
1531
1532
1533
1534
1535
1536 for model in self.model_loop(model_num):
1537
1538
1539
1540 if model_records:
1541
1542 print("\nMODEL %s" % model.num)
1543
1544
1545 file.write("%-6s %4i\n" % ('MODEL', model.num))
1546
1547
1548
1549
1550
1551
1552 for mol in model.mol:
1553
1554 print("ATOM, HETATM, TER")
1555
1556
1557 atom_record = False
1558 for i in range(len(mol.atom_name)):
1559
1560 if mol.pdb_record[i] in [None, 'ATOM']:
1561 atom_record = True
1562
1563
1564 atom_num = mol.atom_num[i]
1565 if atom_num == None:
1566 atom_num = i + 1
1567
1568
1569 file.write("%-6s%5s %4s%1s%3s %1s%4s%1s %8.3f%8.3f%8.3f%6.2f%6.2f %4s%2s%2s\n" % ('ATOM', atom_num, self._translate(mol.atom_name[i]), '', self._translate(mol.res_name[i]), self._translate(mol.chain_id[i]), self._translate(mol.res_num[i]), '', self._translate(mol.x[i], 'float'), self._translate(mol.y[i], 'float'), self._translate(mol.z[i], 'float'), 1.0, 0, self._translate(mol.seg_id[i]), self._translate(mol.element[i]), ''))
1570 num_atom = num_atom + 1
1571
1572
1573 ter_num = atom_num + 1
1574 ter_name = mol.res_name[i]
1575 ter_chain_id = mol.chain_id[i]
1576 ter_res_num = mol.res_num[i]
1577
1578
1579 if atom_record:
1580 file.write("%-6s%5s %3s %1s%4s%1s\n" % ('TER', ter_num, self._translate(ter_name), self._translate(ter_chain_id), self._translate(ter_res_num), ''))
1581 num_ter = num_ter + 1
1582
1583
1584 count_shift = False
1585 for i in range(len(mol.atom_name)):
1586
1587 if mol.pdb_record[i] == 'HETATM':
1588
1589 atom_num = mol.atom_num[i]
1590 if atom_num == None:
1591 atom_num = i + 1
1592
1593
1594 if atom_record and atom_num == ter_num:
1595 count_shift = True
1596 if atom_record and count_shift:
1597 atom_num += 1
1598
1599
1600 file.write("%-6s%5s %4s%1s%3s %1s%4s%1s %8.3f%8.3f%8.3f%6.2f%6.2f %4s%2s%2s\n" % ('HETATM', atom_num, self._translate(mol.atom_name[i]), '', self._translate(mol.res_name[i]), self._translate(mol.chain_id[i]), self._translate(mol.res_num[i]), '', self._translate(mol.x[i], 'float'), self._translate(mol.y[i], 'float'), self._translate(mol.z[i], 'float'), 1.0, 0, self._translate(mol.seg_id[i]), self._translate(mol.element[i]), ''))
1601 num_hetatm = num_hetatm + 1
1602
1603
1604
1605
1606
1607 if model_records:
1608
1609 print("ENDMDL")
1610
1611
1612 file.write("%-6s\n" % 'ENDMDL')
1613
1614
1615
1616
1617
1618
1619 print("CONECT")
1620
1621
1622 for mol in self.structural_data[0].mol:
1623
1624 for i in range(len(mol.atom_name)):
1625
1626 if not len(mol.bonded[i]):
1627 continue
1628
1629
1630 flush = 0
1631 bonded_index = 0
1632 bonded = ['', '', '', '']
1633
1634
1635 for j in range(len(mol.bonded[i])):
1636
1637 if j == len(mol.bonded[i])-1:
1638 flush = True
1639
1640
1641 if bonded_index == 3:
1642 flush = True
1643
1644
1645 bonded[bonded_index] = mol.bonded[i][j]
1646
1647
1648 bonded_index = bonded_index + 1
1649
1650
1651 if flush:
1652
1653 for k in range(4):
1654 if bonded[k] != '':
1655 if mol.atom_num[bonded[k]] != None:
1656 bonded[k] = mol.atom_num[bonded[k]]
1657 else:
1658 bonded[k] = bonded[k] + 1
1659
1660
1661 file.write("%-6s%5s%5s%5s%5s%5s%5s%5s%5s%5s%5s%5s\n" % ('CONECT', i+1, bonded[0], bonded[1], bonded[2], bonded[3], '', '', '', '', '', ''))
1662
1663
1664 flush = False
1665 bonded_index = 0
1666 bonded = ['', '', '', '']
1667
1668
1669 num_conect = num_conect + 1
1670
1671
1672
1673
1674
1675
1676
1677 print("\nMASTER")
1678
1679
1680 file.write("%-6s %5s%5s%5s%5s%5s%5s%5s%5s%5s%5s%5s%5s\n" % ('MASTER', 0, 0, len(het_data_coll), 0, 0, 0, 0, 0, num_atom+num_hetatm, num_ter, num_conect, 0))
1681
1682
1683
1684
1685
1686
1687 print("END")
1688
1689
1690 file.write("END\n")
1691
1692
1694 """The container for the molecular information.
1695
1696 The structural data object for this class is a container possessing a number of different arrays
1697 corresponding to different structural information. These objects include:
1698
1699 - atom_num: The atom name.
1700 - atom_name: The atom name.
1701 - bonded: Each element an array of bonded atom indices.
1702 - chain_id: The chain ID.
1703 - element: The element symbol.
1704 - pdb_record: The optional PDB record name (one of ATOM, HETATM, or TER).
1705 - res_name: The residue name.
1706 - res_num: The residue number.
1707 - seg_id: The segment ID.
1708 - x: The x coordinate of the atom.
1709 - y: The y coordinate of the atom.
1710 - z: The z coordinate of the atom.
1711
1712 All arrays should be of equal length so that an atom index can retrieve all the corresponding
1713 data. Only the atom identification string is compulsory, all other arrays can contain None.
1714 """
1715
1716
1718 """Initialise the molecular container."""
1719
1720
1721 self.atom_num = []
1722
1723
1724 self.atom_name = []
1725
1726
1727 self.bonded = []
1728
1729
1730 self.chain_id = []
1731
1732
1733 self.element = []
1734
1735
1736 self.pdb_record = []
1737
1738
1739 self.res_name = []
1740
1741
1742 self.res_num = []
1743
1744
1745 self.seg_id = []
1746
1747
1748 self.x = []
1749
1750
1751 self.y = []
1752
1753
1754 self.z = []
1755
1756
1758 """Find the atom index corresponding to the given atom number.
1759
1760 @param atom_num: The atom number to find the index of.
1761 @type atom_num: int
1762 @return: The atom index corresponding to the atom.
1763 @rtype: int
1764 """
1765
1766
1767 for j in range(len(self.atom_num)):
1768
1769 if self.atom_num[j] == atom_num:
1770 return j
1771
1772
1773 warn(RelaxWarning("The atom number " + repr(atom_num) + " from the CONECT record cannot be found within the ATOM and HETATM records."))
1774
1775
1777 """Try to determine the element from the PDB atom name.
1778
1779 @param atom_name: The PDB atom name.
1780 @type atom_name: str
1781 @return: The element name, or None if unsuccessful.
1782 @rtype: str or None
1783 """
1784
1785
1786 element = atom_name.strip("'")
1787
1788
1789 element = element.strip(digits)
1790
1791
1792 table = {'C': ['CA', 'CB', 'CG', 'CD', 'CE', 'CH', 'CZ'],
1793 'N': ['ND', 'NE', 'NH', 'NZ'],
1794 'H': ['HA', 'HB', 'HG', 'HD', 'HE', 'HH', 'HT', 'HZ'],
1795 'O': ['OG', 'OD', 'OE', 'OH', 'OT'],
1796 'S': ['SD', 'SG']
1797 }
1798
1799
1800 for key in list(table.keys()):
1801 if element in table[key]:
1802 element = key
1803 break
1804
1805
1806 elements = ['H', 'C', 'N', 'O', 'F', 'P', 'S']
1807
1808
1809 if element in elements:
1810 return element
1811
1812
1813 warn(RelaxWarning("Cannot determine the element associated with atom '%s'." % atom_name))
1814
1815
1817 """Parse the PDB record string and return an array of the corresponding atomic information.
1818
1819 The format of the ATOM and HETATM records is::
1820 __________________________________________________________________________________________
1821 | | | | |
1822 | Columns | Data type | Field | Definition |
1823 |_________|______________|______________|________________________________________________|
1824 | | | | |
1825 | 1 - 6 | Record name | "ATOM" | |
1826 | 7 - 11 | Integer | serial | Atom serial number. |
1827 | 13 - 16 | Atom | name | Atom name. |
1828 | 17 | Character | altLoc | Alternate location indicator. |
1829 | 18 - 20 | Residue name | resName | Residue name. |
1830 | 22 | Character | chainID | Chain identifier. |
1831 | 23 - 26 | Integer | resSeq | Residue sequence number. |
1832 | 27 | AChar | iCode | Code for insertion of residues. |
1833 | 31 - 38 | Real(8.3) | x | Orthogonal coordinates for X in Angstroms. |
1834 | 39 - 46 | Real(8.3) | y | Orthogonal coordinates for Y in Angstroms. |
1835 | 47 - 54 | Real(8.3) | z | Orthogonal coordinates for Z in Angstroms. |
1836 | 55 - 60 | Real(6.2) | occupancy | Occupancy. |
1837 | 61 - 66 | Real(6.2) | tempFactor | Temperature factor. |
1838 | 73 - 76 | LString(4) | segID | Segment identifier, left-justified. |
1839 | 77 - 78 | LString(2) | element | Element symbol, right-justified. |
1840 | 79 - 80 | LString(2) | charge | Charge on the atom. |
1841 |_________|______________|______________|________________________________________________|
1842
1843
1844 The format of the TER record is::
1845 __________________________________________________________________________________________
1846 | | | | |
1847 | Columns | Data type | Field | Definition |
1848 |_________|______________|______________|________________________________________________|
1849 | | | | |
1850 | 1 - 6 | Record name | "TER " | |
1851 | 7 - 11 | Integer | serial | Serial number. |
1852 | 18 - 20 | Residue name | resName | Residue name. |
1853 | 22 | Character | chainID | Chain identifier. |
1854 | 23 - 26 | Integer | resSeq | Residue sequence number. |
1855 | 27 | AChar | iCode | Insertion code. |
1856 |_________|______________|______________|________________________________________________|
1857
1858
1859 The format of the CONECT record is::
1860 __________________________________________________________________________________________
1861 | | | | |
1862 | Columns | Data type | Field | Definition |
1863 |_________|______________|______________|________________________________________________|
1864 | | | | |
1865 | 1 - 6 | Record name | "CONECT" | |
1866 | 7 - 11 | Integer | serial | Atom serial number |
1867 | 12 - 16 | Integer | serial | Serial number of bonded atom |
1868 | 17 - 21 | Integer | serial | Serial number of bonded atom |
1869 | 22 - 26 | Integer | serial | Serial number of bonded atom |
1870 | 27 - 31 | Integer | serial | Serial number of bonded atom |
1871 | 32 - 36 | Integer | serial | Serial number of hydrogen bonded atom |
1872 | 37 - 41 | Integer | serial | Serial number of hydrogen bonded atom |
1873 | 42 - 46 | Integer | serial | Serial number of salt bridged atom |
1874 | 47 - 51 | Integer | serial | Serial number of hydrogen bonded atom |
1875 | 52 - 56 | Integer | serial | Serial number of hydrogen bonded atom |
1876 | 57 - 61 | Integer | serial | Serial number of salt bridged atom |
1877 |_________|______________|______________|________________________________________________|
1878
1879
1880 @param record: The single line PDB record.
1881 @type record: str
1882 @return: The list of atomic information, each element corresponding to the PDB fields
1883 as defined in "Protein Data Bank Contents Guide: Atomic Coordinate Entry
1884 Format Description" version 2.1 (draft), October 25, 1996.
1885 @rtype: list of str
1886 """
1887
1888
1889 fields = []
1890
1891
1892 if search('^ATOM', record) or search('^HETATM', record):
1893
1894 fields.append(record[0:6])
1895 fields.append(record[6:11])
1896 fields.append(record[12:16])
1897 fields.append(record[16])
1898 fields.append(record[17:20])
1899 fields.append(record[21])
1900 fields.append(record[22:26])
1901 fields.append(record[26])
1902 fields.append(record[30:38])
1903 fields.append(record[38:46])
1904 fields.append(record[46:54])
1905 fields.append(record[54:60])
1906 fields.append(record[60:66])
1907 fields.append(record[72:76])
1908 fields.append(record[76:78])
1909 fields.append(record[78:80])
1910
1911
1912 for i in range(len(fields)):
1913
1914 fields[i] = fields[i].strip()
1915
1916
1917 if fields[i] == '':
1918 fields[i] = None
1919
1920
1921 if fields[1]:
1922 fields[1] = int(fields[1])
1923 if fields[6]:
1924 fields[6] = int(fields[6])
1925 if fields[8]:
1926 fields[8] = float(fields[8])
1927 if fields[9]:
1928 fields[9] = float(fields[9])
1929 if fields[10]:
1930 fields[10] = float(fields[10])
1931 if fields[11]:
1932 fields[11] = float(fields[11])
1933 if fields[12]:
1934 fields[12] = float(fields[12])
1935
1936
1937 if search('^CONECT', record):
1938
1939 fields.append(record[0:6])
1940 fields.append(record[6:11])
1941 fields.append(record[11:16])
1942 fields.append(record[16:21])
1943 fields.append(record[21:26])
1944 fields.append(record[26:31])
1945
1946
1947 for i in range(len(fields)):
1948
1949 fields[i] = fields[i].strip()
1950
1951
1952 if fields[i] == '':
1953 fields[i] = None
1954
1955
1956 if fields[1]:
1957 fields[1] = int(fields[1])
1958 if fields[2]:
1959 fields[2] = int(fields[2])
1960 if fields[3]:
1961 fields[3] = int(fields[3])
1962 if fields[4]:
1963 fields[4] = int(fields[4])
1964 if fields[5]:
1965 fields[5] = int(fields[5])
1966
1967
1968 return fields
1969
1970
1972 """Parse the XYZ record string and return an array of the corresponding atomic information.
1973
1974 The format of the XYZ records is::
1975 __________________________________________________________________________________________
1976 | | | | |
1977 | Columns | Data type | Field | Definition |
1978 |_________|______________|______________|________________________________________________|
1979 | | | | |
1980 | 1 | String | element | |
1981 | 2 | Real | x | Orthogonal coordinates for X in Angstroms |
1982 | 3 | Real | y | Orthogonal coordinates for Y in Angstroms |
1983 | 4 | Real | z | Orthogonal coordinates for Z in Angstroms |
1984 |_________|______________|______________|________________________________________________|
1985
1986
1987 @param record: The single line PDB record.
1988 @type record: str
1989 @return: The list of atomic information
1990 @rtype: list of str
1991 """
1992
1993
1994 fields = []
1995 word = record.split()
1996
1997
1998 if len(word)==4:
1999
2000 fields.append(word[0])
2001 fields.append(word[1])
2002 fields.append(word[2])
2003 fields.append(word[3])
2004
2005
2006 for i in range(len(fields)):
2007
2008 fields[i] = fields[i].strip()
2009
2010
2011 if fields[i] == '':
2012 fields[i] = None
2013
2014
2015 if fields[1]:
2016 fields[1] = float(fields[1])
2017 if fields[2]:
2018 fields[2] = float(fields[2])
2019 if fields[3]:
2020 fields[3] = float(fields[3])
2021
2022
2023 return fields
2024
2025
2026 - def atom_add(self, atom_name=None, res_name=None, res_num=None, pos=[None, None, None], element=None, atom_num=None, chain_id=None, segment_id=None, pdb_record=None):
2027 """Method for adding an atom to the structural data object.
2028
2029 This method will create the key-value pair for the given atom.
2030
2031
2032 @keyword atom_name: The atom name, e.g. 'H1'.
2033 @type atom_name: str or None
2034 @keyword res_name: The residue name.
2035 @type res_name: str or None
2036 @keyword res_num: The residue number.
2037 @type res_num: int or None
2038 @keyword pos: The position vector of coordinates.
2039 @type pos: list (length = 3)
2040 @keyword element: The element symbol.
2041 @type element: str or None
2042 @keyword atom_num: The atom number.
2043 @type atom_num: int or None
2044 @keyword chain_id: The chain identifier.
2045 @type chain_id: str or None
2046 @keyword segment_id: The segment identifier.
2047 @type segment_id: str or None
2048 @keyword pdb_record: The optional PDB record name, e.g. 'ATOM' or 'HETATM'.
2049 @type pdb_record: str or None
2050 """
2051
2052
2053 self.atom_num.append(atom_num)
2054 self.atom_name.append(atom_name)
2055 self.bonded.append([])
2056 self.chain_id.append(chain_id)
2057 self.element.append(element)
2058 self.pdb_record.append(pdb_record)
2059 self.res_name.append(res_name)
2060 self.res_num.append(res_num)
2061 self.seg_id.append(segment_id)
2062 self.x.append(pos[0])
2063 self.y.append(pos[1])
2064 self.z.append(pos[2])
2065
2066
2068 """Method for connecting two atoms within the data structure object.
2069
2070 This method will append index2 to the array at bonded[index1] and vice versa.
2071
2072
2073 @keyword index1: The index of the first atom.
2074 @type index1: int
2075 @keyword index2: The index of the second atom.
2076 @type index2: int
2077 """
2078
2079
2080 if index2 not in self.bonded[index1]:
2081 self.bonded[index1].append(index2)
2082 if index1 not in self.bonded[index2]:
2083 self.bonded[index2].append(index1)
2084
2085
2087 """Method for generating a complete Structure_container object from the given PDB records.
2088
2089 @param records: A list of structural PDB records.
2090 @type records: list of str
2091 """
2092
2093
2094 for record in records:
2095
2096 record = self._parse_pdb_record(record)
2097
2098
2099 if not record:
2100 continue
2101
2102
2103 if record[0] == 'ATOM' or record[0] == 'HETATM':
2104
2105 element = record[14]
2106 if not element:
2107 element = self._det_pdb_element(record[2])
2108
2109
2110 self.atom_add(pdb_record=record[0], atom_num=record[1], atom_name=record[2], res_name=record[4], chain_id=record[5], res_num=record[6], pos=[record[8], record[9], record[10]], segment_id=record[13], element=element)
2111
2112
2113 if record[0] == 'CONECT':
2114
2115 for i in range(len(record)-2):
2116
2117 if record[i+2] == None:
2118 continue
2119
2120
2121 self.atom_connect(index1=self._atom_index(record[1]), index2=self._atom_index(record[i+2]))
2122
2123
2125 """Method for generating a complete Structure_container object from the given xyz records.
2126
2127 @param records: A list of structural xyz records.
2128 @type records: list of str
2129 """
2130
2131
2132 atom_number = 1
2133
2134
2135 for record in records:
2136
2137 record = self._parse_xyz_record(record)
2138
2139
2140 if not record:
2141 continue
2142
2143
2144 if len(record) == 4:
2145
2146 element = record[0]
2147 if not element:
2148 element = self._det_pdb_element(record[2])
2149
2150
2151 self.atom_add(atom_num=atom_number, pos=[record[1], record[2], record[3]], element=element)
2152
2153
2154 atom_number = atom_number + 1
2155
2156
2157 - def from_xml(self, mol_node, file_version=1):
2158 """Recreate the MolContainer from the XML molecule node.
2159
2160 @param mol_node: The molecule XML node.
2161 @type mol_node: xml.dom.minicompat.NodeList instance
2162 @keyword file_version: The relax XML version of the XML file.
2163 @type file_version: int
2164 """
2165
2166
2167 xml_to_object(mol_node, self, file_version=file_version)
2168
2169
2171 """Check if the container is empty."""
2172
2173
2174 if hasattr(self, 'mol_name'): return False
2175 if hasattr(self, 'file_name'): return False
2176 if hasattr(self, 'file_path'): return False
2177 if hasattr(self, 'file_mol_num'): return False
2178 if hasattr(self, 'file_model'): return False
2179
2180
2181 if not self.atom_num == []: return False
2182 if not self.atom_name == []: return False
2183 if not self.bonded == []: return False
2184 if not self.chain_id == []: return False
2185 if not self.element == []: return False
2186 if not self.pdb_record == []: return False
2187 if not self.res_name == []: return False
2188 if not self.res_num == []: return False
2189 if not self.seg_id == []: return False
2190 if not self.x == []: return False
2191 if not self.y == []: return False
2192 if not self.z == []: return False
2193
2194
2195 return True
2196
2197
2198 - def to_xml(self, doc, element):
2199 """Create XML elements for the contents of this molecule container.
2200
2201 @param doc: The XML document object.
2202 @type doc: xml.dom.minidom.Document instance
2203 @param element: The element to add the molecule XML elements to.
2204 @type element: XML element object
2205 """
2206
2207
2208 mol_element = doc.createElement('mol_cont')
2209 element.appendChild(mol_element)
2210
2211
2212 mol_element.setAttribute('desc', 'Molecule container')
2213 mol_element.setAttribute('name', str(self.mol_name))
2214
2215
2216 fill_object_contents(doc, mol_element, object=self, blacklist=list(self.__class__.__dict__.keys()))
2217