Package generic_fns :: Package structure :: Module internal
[hide private]
[frames] | no frames]

Source Code for Module generic_fns.structure.internal

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