generic_fns.structure.main

1 ############################################################################### 2 # # 3 # Copyright (C) 2003-2011 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 # Python module imports. 24 from math import sqrt 25 from minfx.generic import generic_minimise 26 from numpy import array, dot, float64, ndarray, zeros 27 from numpy.linalg import norm 28 from os import F_OK, access 29 from re import search 30 from string import replace 31 import sys 32 from warnings import warn 33 34 # relax module imports. 35 from generic_fns import molmol, relax_re 36 from generic_fns.mol_res_spin import create_spin, exists_mol_res_spin_data, generate_spin_id, linear_ave, return_molecule, return_residue, return_spin, spin_loop 37 from generic_fns import pipes 38 from generic_fns.structure.api_base import Displacements 39 from generic_fns.structure.internal import Internal 40 from generic_fns.structure.scientific import Scientific_data 41 from generic_fns.structure.superimpose import fit_to_first, fit_to_mean, Pivot_finder 42 from relax_errors import RelaxError, RelaxFileError, RelaxNoPdbError, RelaxNoSequenceError 43 from relax_io import get_file_path, open_write_file, write_spin_data 44 from relax_warnings import RelaxWarning, RelaxNoPDBFileWarning, RelaxZeroVectorWarning 45 46

47 -def add_atom(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):

48 """Add a new atom to the structural data object. 49 50 @keyword mol_name: The name of the molecule. 51 @type mol_name: str 52 @keyword atom_name: The atom name, e.g. 'H1'. 53 @type atom_name: str or None 54 @keyword res_name: The residue name. 55 @type res_name: str or None 56 @keyword res_num: The residue number. 57 @type res_num: int or None 58 @keyword pos: The position vector of coordinates. 59 @type pos: list (length = 3) 60 @keyword element: The element symbol. 61 @type element: str or None 62 @keyword atom_num: The atom number. 63 @type atom_num: int or None 64 @keyword chain_id: The chain identifier. 65 @type chain_id: str or None 66 @keyword segment_id: The segment identifier. 67 @type segment_id: str or None 68 @keyword pdb_record: The optional PDB record name, e.g. 'ATOM' or 'HETATM'. 69 @type pdb_record: str or None 70 """ 71 72 # Test if the current data pipe exists. 73 pipes.test() 74 75 # Place the structural object into the relax data store if needed. 76 if not hasattr(cdp, 'structure'): 77 cdp.structure = Internal() 78 79 # Add the atoms. 80 cdp.structure.add_atom(mol_name=mol_name, 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)

81 82

83 -def connect_atom(index1=None, index2=None):

84 """Connect two atoms. 85 86 @keyword index1: The global index of the first atom. 87 @type index1: str 88 @keyword index2: The global index of the first atom. 89 @type index2: str 90 """ 91 92 # Test if the current data pipe exists. 93 pipes.test() 94 95 # Place the structural object into the relax data store if needed. 96 if not hasattr(cdp, 'structure'): 97 cdp.structure = Internal() 98 99 # Add the atoms. 100 cdp.structure.connect_atom(index1=index1, index2=index2)

101 102

103 -def delete():

104 """Simple function for deleting all structural data.""" 105 106 # Run the object method. 107 cdp.structure.delete() 108 109 # Then remove any spin specific structural info. 110 for spin in spin_loop(): 111 # Delete positional information. 112 if hasattr(spin, 'pos'): 113 del spin.pos 114 115 # Delete bond vectors. 116 if hasattr(spin, 'bond_vect'): 117 del spin.bond_vect 118 if hasattr(spin, 'xh_vect'): 119 del spin.xh_vect

120 121

122 -def displacement(model_from=None, model_to=None, atom_id=None, centroid=None):

123 """Calculate the rotational and translational displacement between two structural models. 124 125 @keyword model_from: The optional model number for the starting position of the displacement. 126 @type model_from: int or None 127 @keyword model_to: The optional model number for the ending position of the displacement. 128 @type model_to: int or None 129 @keyword atom_id: The molecule, residue, and atom identifier string. This matches the spin ID string format. 130 @type atom_id: str or None 131 @keyword centroid: An alternative position of the centroid, used for studying pivoted systems. 132 @type centroid: list of float or numpy rank-1, 3D array 133 """ 134 135 # Convert the model_from and model_to args to lists, is supplied. 136 if model_from != None: 137 model_from = [model_from] 138 if model_to != None: 139 model_to = [model_to] 140 141 # Create a list of all models. 142 models = [] 143 for model in cdp.structure.model_loop(): 144 models.append(model.num) 145 146 # Set model_from or model_to to all models if None. 147 if model_from == None: 148 model_from = models 149 if model_to == None: 150 model_to = models 151 152 # Initialise the data structure. 153 if not hasattr(cdp.structure, 'displacements'): 154 cdp.structure.displacements = Displacements() 155 156 # Loop over the starting models. 157 for i in range(len(model_from)): 158 # Assemble the atomic coordinates. 159 coord_from = [] 160 for pos in cdp.structure.atom_loop(atom_id=atom_id, model_num=model_from[i], pos_flag=True): 161 coord_from.append(pos[0]) 162 163 # Loop over the ending models. 164 for j in range(len(model_to)): 165 # Assemble the atomic coordinates. 166 coord_to = [] 167 for pos in cdp.structure.atom_loop(atom_id=atom_id, model_num=model_to[j], pos_flag=True): 168 coord_to.append(pos[0]) 169 170 # Send to the base container for the calculations. 171 cdp.structure.displacements._calculate(model_from=model_from[i], model_to=model_to[j], coord_from=array(coord_from), coord_to=array(coord_to), centroid=centroid)

172 173

174 -def find_pivot(models=None, atom_id=None, init_pos=None):

175 """Superimpose a set of structural models. 176 177 @keyword models: The list of models to use. If set to None, then all models will be used. 178 @type models: list of int or None 179 @keyword atom_id: The molecule, residue, and atom identifier string. This matches the spin ID string format. 180 @type atom_id: str or None 181 @keyword init_pos: The starting pivot position for the pivot point optimisation. 182 @type init_pos: list of float or numpy rank-1, 3D array 183 """ 184 185 # Initialised the starting position if needed. 186 if init_pos == None: 187 init_pos = zeros(3, float64) 188 init_pos = array(init_pos) 189 190 # Validate the models. 191 cdp.structure.validate_models() 192 193 # Create a list of all models. 194 if models == None: 195 models = [] 196 for model in cdp.structure.model_loop(): 197 models.append(model.num) 198 199 # Assemble the atomic coordinates of all models. 200 coord = [] 201 for model in models: 202 coord.append([]) 203 for pos in cdp.structure.atom_loop(atom_id=atom_id, model_num=model, pos_flag=True): 204 coord[-1].append(pos[0]) 205 coord[-1] = array(coord[-1]) 206 coord = array(coord) 207 208 # The target function. 209 finder = Pivot_finder(models, coord) 210 results = generic_minimise(func=finder.func, x0=init_pos, min_algor='simplex', print_flag=1) 211 212 # No result. 213 if results == None: 214 return 215 216 # Store the data. 217 cdp.structure.pivot = results 218 219 # Print out. 220 print("Motional pivot found at: %s" % results)

221 222 223 224

225 -def get_pos(spin_id=None, str_id=None, ave_pos=False):

226 """Load the spins from the structural object into the relax data store. 227 228 @keyword spin_id: The molecule, residue, and spin identifier string. 229 @type spin_id: str 230 @keyword str_id: The structure identifier. This can be the file name, model number, or structure number. 231 @type str_id: int or str 232 @keyword ave_pos: A flag specifying if the average atom position or the atom position from all loaded structures is loaded into the SpinContainer. 233 @type ave_pos: bool 234 """ 235 236 # Test if the current data pipe exists. 237 pipes.test() 238 239 # Test if the structure exists. 240 if not hasattr(cdp, 'structure') or not cdp.structure.num_models() or not cdp.structure.num_molecules(): 241 raise RelaxNoPdbError 242 243 # Loop over all atoms of the spin_id selection. 244 model_index = -1 245 last_model = None 246 for model_num, mol_name, res_num, res_name, atom_num, atom_name, element, pos in cdp.structure.atom_loop(atom_id=spin_id, str_id=str_id, model_num_flag=True, 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, ave=ave_pos): 247 # Update the model info. 248 if last_model != model_num: 249 model_index = model_index + 1 250 last_model = model_num 251 252 # Remove the '+' regular expression character from the mol, res, and spin names! 253 if mol_name and search('\+', mol_name): 254 mol_name = replace(mol_name, '+', '') 255 if res_name and search('\+', res_name): 256 res_name = replace(res_name, '+', '') 257 if atom_name and search('\+', atom_name): 258 atom_name = replace(atom_name, '+', '') 259 260 # The spin identification string. The residue name and spin num is not included to allow molecules with point mutations to be used as different models. 261 id = generate_spin_id(mol_name=mol_name, res_num=res_num, res_name=None, spin_name=atom_name) 262 263 # Get the spin container. 264 spin_cont = return_spin(id) 265 266 # Skip the spin if it doesn't exist. 267 if spin_cont == None: 268 continue 269 270 # Add the position vector to the spin container. 271 if ave_pos: 272 spin_cont.pos = pos 273 else: 274 if not hasattr(spin_cont, 'pos'): 275 spin_cont.pos = [] 276 spin_cont.pos.append(pos) 277 278 # Update pseudo-atoms. 279 for spin in spin_loop(): 280 if hasattr(spin, 'members'): 281 # Get the spin positions. 282 positions = [] 283 for atom in spin.members: 284 # Get the spin container. 285 subspin = return_spin(atom) 286 287 # Test that the spin exists. 288 if subspin == None: 289 raise RelaxNoSpinError(atom) 290 291 # Test the position. 292 if not hasattr(subspin, 'pos') or not subspin.pos: 293 raise RelaxError("Positional information is not available for the atom '%s'." % atom) 294 295 # Alias the position. 296 pos = subspin.pos 297 298 # Convert to a list of lists if not already. 299 multi_model = True 300 if type(pos[0]) in [float, float64]: 301 multi_model = False 302 pos = [pos] 303 304 # Store the position. 305 positions.append([]) 306 for i in range(len(pos)): 307 positions[-1].append(pos[i].tolist()) 308 309 # The averaging. 310 if spin.averaging == 'linear': 311 # Average pos. 312 ave = linear_ave(positions) 313 314 # Convert to the correct structure. 315 if multi_model: 316 spin.pos = ave 317 else: 318 spin.pos = ave[0]

319 320

321 -def load_spins(spin_id=None, str_id=None, ave_pos=False):

322 """Load the spins from the structural object into the relax data store. 323 324 @keyword spin_id: The molecule, residue, and spin identifier string. 325 @type spin_id: str 326 @keyword str_id: The structure identifier. This can be the file name, model number, or structure number. 327 @type str_id: int or str 328 @keyword ave_pos: A flag specifying if the average atom position or the atom position from all loaded structures is loaded into the SpinContainer. 329 @type ave_pos: bool 330 """ 331 332 # Test if the current data pipe exists. 333 pipes.test() 334 335 # Test if the structure exists. 336 if not hasattr(cdp, 'structure') or not cdp.structure.num_models() or not cdp.structure.num_molecules(): 337 raise RelaxNoPdbError 338 339 # Print out. 340 print("Adding the following spins to the relax data store.\n") 341 342 # Init the data for printing out. 343 mol_names = [] 344 res_nums = [] 345 res_names = [] 346 spin_nums = [] 347 spin_names = [] 348 349 # Initialise data for the atom loop. 350 model_index = -1 351 last_model = 1000000 352 353 # Loop over all atoms of the spin_id selection. 354 for model_num, mol_name, res_num, res_name, atom_num, atom_name, element, pos in cdp.structure.atom_loop(atom_id=spin_id, str_id=str_id, model_num_flag=True, 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, ave=ave_pos): 355 # Update the model info. 356 if last_model != model_num: 357 model_index = model_index + 1 358 last_model = model_num 359 360 # Remove the '+' regular expression character from the mol, res, and spin names! 361 if mol_name and search('\+', mol_name): 362 mol_name = replace(mol_name, '+', '') 363 if res_name and search('\+', res_name): 364 res_name = replace(res_name, '+', '') 365 if atom_name and search('\+', atom_name): 366 atom_name = replace(atom_name, '+', '') 367 368 # Generate a spin ID for the current atom. 369 id = generate_spin_id(mol_name=mol_name, res_num=res_num, res_name=res_name, spin_num=atom_num, spin_name=atom_name) 370 371 # Create the spin. 372 try: 373 spin_cont = create_spin(mol_name=mol_name, res_num=res_num, res_name=res_name, spin_num=atom_num, spin_name=atom_name) 374 except RelaxError: 375 # Throw a warning if still on the first model. 376 if not model_index: 377 warn(RelaxWarning("The spin '%s' already exists." % id)) 378 continue 379 380 # Otherwise, get the spin container. 381 spin_cont = return_spin(id) 382 383 # Append all the spin ID info for the first model for printing later. 384 if model_index == 0: 385 mol_names.append(mol_name) 386 res_nums.append(res_num) 387 res_names.append(res_name) 388 spin_nums.append(atom_num) 389 spin_names.append(atom_name) 390 391 # Convert the position vector to a numpy array. 392 pos = array(pos, float64) 393 394 # Average position vector (already averaged across models in the atom_loop). 395 if ave_pos: 396 spin_cont.pos = pos 397 398 # All positions. 399 else: 400 # Initialise. 401 if not hasattr(spin_cont, 'pos'): 402 spin_cont.pos = [] 403 404 # Add the current model's position. 405 spin_cont.pos.append(pos) 406 407 # Add the element. 408 if not model_index: 409 spin_cont.element = element 410 411 # Catch no data. 412 if len(mol_names) == 0: 413 warn(RelaxWarning("No spins matching the '%s' ID string could be found." % spin_id)) 414 return 415 416 # Print out. 417 write_spin_data(file=sys.stdout, mol_names=mol_names, res_nums=res_nums, res_names=res_names, spin_nums=spin_nums, spin_names=spin_names)

418 419

420 -def read_pdb(file=None, dir=None, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, parser='internal', verbosity=1, fail=True):

421 """The PDB loading function. 422 423 Parsers 424 ======= 425 426 A number of parsers are available for reading PDB files. These include: 427 428 - 'scientific', the Scientific Python PDB parser. 429 - 'internal', a low quality yet fast PDB parser built into relax. 430 431 432 @keyword file: The name of the PDB file to read. 433 @type file: str 434 @keyword dir: The directory where the PDB file is located. If set to None, then the file will be searched for in the current directory. 435 @type dir: str or None 436 @keyword read_mol: The molecule(s) to read from the file, independent of model. The molecules are determined differently by the different parsers, but are numbered consecutively from 1. If set to None, then all molecules will be loaded. 437 @type read_mol: None, int, or list of int 438 @keyword set_mol_name: Set the names of the molecules which are loaded. If set to None, then the molecules will be automatically labelled based on the file name or other information. 439 @type set_mol_name: None, str, or list of str 440 @keyword read_model: The PDB model to extract from the file. If set to None, then all models will be loaded. 441 @type read_model: None, int, or list of int 442 @keyword set_model_num: Set the model number of the loaded molecule. If set to None, then the PDB model numbers will be preserved, if they exist. 443 @type set_model_num: None, int, or list of int 444 @keyword parser: The parser to be used to read the PDB file. 445 @type parser: str 446 @keyword fail: A flag which, if True, will cause a RelaxError to be raised if the PDB file does not exist. If False, then a RelaxWarning will be trown instead. 447 @type fail: bool 448 @keyword verbosity: The amount of information to print to screen. Zero corresponds to minimal output while higher values increase the amount of output. The default value is 1. 449 @type verbosity: int 450 @raise RelaxFileError: If the fail flag is set, then a RelaxError is raised if the PDB file does not exist. 451 """ 452 453 # Test if the current data pipe exists. 454 pipes.test() 455 456 # The file path. 457 file_path = get_file_path(file, dir) 458 459 # Try adding file extensions to the end of the file path, if the file can't be found. 460 file_path_orig = file_path 461 if not access(file_path, F_OK): 462 # List of possible extensions. 463 for ext in ['.pdb', '.gz', '.pdb.gz', '.bz2', '.pdb.bz2']: 464 # Add the extension if the file can be found. 465 if access(file_path+ext, F_OK): 466 file_path = file_path + ext 467 468 # Test if the file exists. 469 if not access(file_path, F_OK): 470 if fail: 471 raise RelaxFileError('PDB', file_path_orig) 472 else: 473 warn(RelaxNoPDBFileWarning(file_path)) 474 return 475 476 # Check that the parser is the same as the currently loaded PDB files. 477 if hasattr(cdp, 'structure') and cdp.structure.id != parser: 478 raise RelaxError("The " + repr(parser) + " parser does not match the " + repr(cdp.structure.id) + " parser of the PDB loaded into the current pipe.") 479 480 # Place the parser specific structural object into the relax data store. 481 if not hasattr(cdp, 'structure'): 482 if parser == 'scientific': 483 cdp.structure = Scientific_data() 484 elif parser == 'internal': 485 cdp.structure = Internal() 486 487 # Load the structures. 488 cdp.structure.load_pdb(file_path, read_mol=read_mol, set_mol_name=set_mol_name, read_model=read_model, set_model_num=set_model_num, verbosity=verbosity) 489 490 # Load into Molmol (if running). 491 molmol.molmol_obj.open_pdb()

492 493

494 -def read_xyz(file=None, dir=None, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, verbosity=1, fail=True):

495 """The XYZ loading function. 496 497 498 @keyword file: The name of the XYZ file to read. 499 @type file: str 500 @keyword dir: The directory where the XYZ file is located. If set to None, then the 501 file will be searched for in the current directory. 502 @type dir: str or None 503 @keyword read_mol: The molecule(s) to read from the file, independent of model. 504 If set to None, then all molecules will be loaded. 505 @type read_mol: None, int, or list of int 506 @keyword set_mol_name: Set the names of the molecules which are loaded. If set to None, then 507 the molecules will be automatically labelled based on the file name or 508 other information. 509 @type set_mol_name: None, str, or list of str 510 @keyword read_model: The XYZ model to extract from the file. If set to None, then all models 511 will be loaded. 512 @type read_model: None, int, or list of int 513 @keyword set_model_num: Set the model number of the loaded molecule. If set to None, then the 514 XYZ model numbers will be preserved, if they exist. 515 @type set_model_num: None, int, or list of int 516 @keyword fail: A flag which, if True, will cause a RelaxError to be raised if the XYZ 517 file does not exist. If False, then a RelaxWarning will be trown 518 instead. 519 @type fail: bool 520 @keyword verbosity: The amount of information to print to screen. Zero corresponds to 521 minimal output while higher values increase the amount of output. The 522 default value is 1. 523 @type verbosity: int 524 @raise RelaxFileError: If the fail flag is set, then a RelaxError is raised if the XYZ file 525 does not exist. 526 """ 527 528 # Test if the current data pipe exists. 529 pipes.test() 530 531 # The file path. 532 file_path = get_file_path(file, dir) 533 534 # Try adding '.xyz' to the end of the file path, if the file can't be found. 535 if not access(file_path, F_OK): 536 file_path_orig = file_path 537 file_path = file_path + '.xyz' 538 539 # Test if the file exists. 540 if not access(file_path, F_OK): 541 if fail: 542 raise RelaxFileError('XYZ', file_path_orig) 543 else: 544 warn(RelaxNoPDBFileWarning(file_path)) 545 return 546 547 # Place the structural object into the relax data store. 548 if not hasattr(cdp, 'structure'): 549 cdp.structure = Internal() 550 551 # Load the structures. 552 cdp.structure.load_xyz(file_path, read_mol=read_mol, set_mol_name=set_mol_name, read_model=read_model, set_model_num=set_model_num, verbosity=verbosity)

553 554

555 -def rotate(R=None, origin=None, model=None, atom_id=None):

556 """Rotate the structural data about the origin by the specified forwards rotation. 557 558 @keyword R: The forwards rotation matrix. 559 @type R: numpy 3D, rank-2 array or a 3x3 list of floats 560 @keyword origin: The origin of the rotation. If not supplied, the origin will be set to [0, 0, 0]. 561 @type origin: numpy 3D, rank-1 array or list of len 3 or None 562 @keyword model: The model to rotate. If None, all models will be rotated. 563 @type model: int 564 @keyword atom_id: The molecule, residue, and atom identifier string. Only atoms matching this selection will be used. 565 @type atom_id: str or None 566 """ 567 568 # Test if the current data pipe exists. 569 pipes.test() 570 571 # Test if the structure exists. 572 if not hasattr(cdp, 'structure') or not cdp.structure.num_models() or not cdp.structure.num_molecules(): 573 raise RelaxNoPdbError 574 575 # Set the origin if not supplied. 576 if origin == None: 577 origin = [0, 0, 0] 578 579 # Convert the args to numpy float data structures. 580 R = array(R, float64) 581 origin = array(origin, float64) 582 583 # Call the specific code. 584 cdp.structure.rotate(R=R, origin=origin, model=model, atom_id=atom_id)

585 586

587 -def set_vector(spin=None, xh_vect=None):

588 """Place the XH unit vector into the spin container object. 589 590 @keyword spin: The spin container object. 591 @type spin: SpinContainer instance 592 @keyword xh_vect: The unit vector parallel to the XH bond. 593 @type xh_vect: array of len 3 594 """ 595 596 # Place the XH unit vector into the container. 597 spin.xh_vect = xh_vect

598 599

600 -def superimpose(models=None, method='fit to mean', atom_id=None, centroid=None):

601 """Superimpose a set of structural models. 602 603 @keyword models: The list of models to superimpose. If set to None, then all models will be used. 604 @type models: list of int or None 605 @keyword method: The superimposition method. It must be one of 'fit to mean' or 'fit to first'. 606 @type method: str 607 @keyword atom_id: The molecule, residue, and atom identifier string. This matches the spin ID string format. 608 @type atom_id: str or None 609 @keyword centroid: An alternative position of the centroid to allow for different superpositions, for example of pivot point motions. 610 @type centroid: list of float or numpy rank-1, 3D array 611 """ 612 613 # Check the method. 614 allowed = ['fit to mean', 'fit to first'] 615 if method not in allowed: 616 raise RelaxError("The superimposition method '%s' is unknown. It must be one of %s." % (method, allowed)) 617 618 # Validate the models. 619 cdp.structure.validate_models() 620 621 # Create a list of all models. 622 if models == None: 623 models = [] 624 for model in cdp.structure.model_loop(): 625 models.append(model.num) 626 627 # Assemble the atomic coordinates of all models. 628 coord = [] 629 for model in models: 630 coord.append([]) 631 for pos in cdp.structure.atom_loop(atom_id=atom_id, model_num=model, pos_flag=True): 632 coord[-1].append(pos[0]) 633 coord[-1] = array(coord[-1]) 634 635 # The different algorithms. 636 if method == 'fit to mean': 637 T, R, pivot = fit_to_mean(models=models, coord=coord, centroid=centroid) 638 elif method == 'fit to first': 639 T, R, pivot = fit_to_first(models=models, coord=coord, centroid=centroid) 640 641 642 # Update to the new coordinates. 643 for i in range(len(models)): 644 # Translate the molecule first (the rotational pivot is defined in the first model). 645 translate(T=T[i], model=models[i]) 646 647 # Rotate the molecule. 648 rotate(R=R[i], origin=pivot[i], model=models[i])

649 650

651 -def translate(T=None, model=None, atom_id=None):

652 """Shift the structural data by the specified translation vector. 653 654 @keyword T: The translation vector. 655 @type T: numpy rank-1, 3D array or list of float 656 @keyword model: The model to translate. If None, all models will be rotated. 657 @type model: int or None 658 @keyword atom_id: The molecule, residue, and atom identifier string. Only atoms matching this selection will be used. 659 @type atom_id: str or None 660 """ 661 662 # Test if the current data pipe exists. 663 pipes.test() 664 665 # Test if the structure exists. 666 if not hasattr(cdp, 'structure') or not cdp.structure.num_models() or not cdp.structure.num_molecules(): 667 raise RelaxNoPdbError 668 669 # Convert the args to numpy float data structures. 670 T = array(T, float64) 671 672 # Call the specific code. 673 cdp.structure.translate(T=T, model=model, atom_id=atom_id)

674 675

676 -def vectors(attached=None, spin_id=None, model=None, verbosity=1, ave=True, unit=True):

677 """Extract the bond vectors from the loaded structures and store them in the spin container. 678 679 @keyword attached: The name of the atom attached to the spin, as given in the structural file. Regular expression can be used, for example 'H*'. This uses relax rather than Python regular expression (i.e. shell like syntax). 680 @type attached: str 681 @keyword spin_id: The spin identifier string. 682 @type spin_id: str 683 @keyword model: The model to extract the vector from. If None, all vectors will be extracted. 684 @type model: str 685 @keyword verbosity: The higher the value, the more information is printed to screen. 686 @type verbosity: int 687 @keyword ave: A flag which if True will cause the average of all vectors to be extracted. 688 @type ave: bool 689 @keyword unit: A flag which if True will cause the function to calculate the unit vectors. 690 @type unit: bool 691 """ 692 693 # Test if the PDB file has been loaded. 694 if not hasattr(cdp, 'structure'): 695 raise RelaxNoPdbError 696 697 # Test if sequence data is loaded. 698 if not exists_mol_res_spin_data(): 699 raise RelaxNoSequenceError 700 701 # Print out. 702 if verbosity: 703 # Number of models. 704 num_models = cdp.structure.num_models() 705 706 # Multiple models loaded. 707 if num_models > 1: 708 if model: 709 print(("Extracting vectors for model '%s'." % model)) 710 else: 711 print(("Extracting vectors for all %s models." % num_models)) 712 if ave: 713 print("Averaging all vectors.") 714 715 # Single model loaded. 716 else: 717 print("Extracting vectors from the single model.") 718 719 # Unit vectors. 720 if unit: 721 print("Calculating the unit vectors.") 722 723 # Determine if the attached atom is a proton. 724 proton = False 725 if relax_re.search('.*H.*', attached) or relax_re.search(attached, 'H'): 726 proton = True 727 if verbosity: 728 if proton: 729 print("The attached atom is a proton.") 730 else: 731 print("The attached atom is not a proton.") 732 print('') 733 734 # Set the variable name in which the vectors will be stored. 735 if proton: 736 object_name = 'xh_vect' 737 else: 738 object_name = 'bond_vect' 739 740 # Loop over the spins. 741 no_vectors = True 742 for spin, mol_name, res_num, res_name in spin_loop(selection=spin_id, full_info=True): 743 # Skip deselected spins. 744 if not spin.select: 745 continue 746 747 # The spin identification string. The residue name and spin num is not included to allow molecules with point mutations to be used as different models. 748 id = generate_spin_id(res_num=res_num, res_name=None, spin_name=spin.name, spin_num=spin.num) 749 750 # Test that the spin number or name are set (one or both are essential for the identification of the atom). 751 if spin.num == None and spin.name == None: 752 warn(RelaxWarning("Either the spin number or name must be set for the spin " + repr(id) + " to identify the corresponding atom in the molecule.")) 753 continue 754 755 # The bond vector already exists. 756 if hasattr(spin, object_name): 757 obj = getattr(spin, object_name) 758 if obj != None: 759 warn(RelaxWarning("The bond vector for the spin " + repr(id) + " already exists.")) 760 continue 761 762 # Get the bond info. 763 bond_vectors, attached_name, warnings = cdp.structure.bond_vectors(attached_atom=attached, model_num=model, res_num=res_num, spin_name=spin.name, spin_num=spin.num, return_name=True, return_warnings=True) 764 765 # No attached atom. 766 if not bond_vectors: 767 # Warning messages. 768 if warnings: 769 warn(RelaxWarning(warnings + " (atom ID " + repr(id) + ").")) 770 771 # Skip the spin. 772 continue 773 774 # Set the attached atom name. 775 if not hasattr(spin, 'attached_atom'): 776 spin.attached_atom = attached_name 777 elif spin.attached_atom != attached_name: 778 raise RelaxError("The " + repr(spin.attached_atom) + " atom already attached to the spin does not match the attached atom " + repr(attached_name) + ".") 779 780 # Initialise the average vector. 781 if ave: 782 ave_vector = zeros(3, float64) 783 784 # Loop over the individual vectors. 785 for i in xrange(len(bond_vectors)): 786 # Unit vector. 787 if unit: 788 # Normalisation factor. 789 norm_factor = norm(bond_vectors[i]) 790 791 # Test for zero length. 792 if norm_factor == 0.0: 793 warn(RelaxZeroVectorWarning(id)) 794 795 # Calculate the normalised vector. 796 else: 797 bond_vectors[i] = bond_vectors[i] / norm_factor 798 799 # Sum the vectors. 800 if ave: 801 ave_vector = ave_vector + bond_vectors[i] 802 803 # Average. 804 if ave: 805 vector = ave_vector / float(len(bond_vectors)) 806 else: 807 vector = bond_vectors 808 809 # Set the vector. 810 setattr(spin, object_name, vector) 811 812 # We have a vector! 813 no_vectors = False 814 815 # Print out of modified spins. 816 if verbosity: 817 if spin.name: 818 print(("Extracted " + spin.name + "-" + str(attached_name) + " vectors for " + repr(id) + '.')) 819 else: 820 print(("Extracted " + str(spin.num) + "-" + str(attached_name) + " vectors for " + repr(id) + '.')) 821 822 # Right, catch the problem of missing vectors to prevent massive user confusion! 823 if no_vectors: 824 raise RelaxError("No vectors could be extracted.")

825 826

827 -def write_pdb(file=None, dir=None, model_num=None, compress_type=0, force=False):

828 """The PDB writing function. 829 830 @keyword file: The name of the PDB file to write. 831 @type file: str 832 @keyword dir: The directory where the PDB file will be placed. If set to None, then the file will be placed in the current directory. 833 @type dir: str or None 834 @keyword model_num: The model to place into the PDB file. If not supplied, then all models will be placed into the file. 835 @type model_num: None or int 836 @keyword compress_type: The compression type. The integer values correspond to the compression type: 0, no compression; 1, Bzip2 compression; 2, Gzip compression. 837 @type compress_type: int 838 @keyword force: The force flag which if True will cause the file to be overwritten. 839 @type force: bool 840 """ 841 842 # Test if the current data pipe exists. 843 pipes.test() 844 845 # Check if the structural object exists. 846 if not hasattr(cdp, 'structure'): 847 raise RelaxError("No structural data is present in the current data pipe.") 848 849 # Check if the structural object is writable. 850 if cdp.structure.id in ['scientific']: 851 raise RelaxError("The structures from the " + cdp.structure.id + " parser are not writable.") 852 853 # The file path. 854 file_path = get_file_path(file, dir) 855 856 # Add '.pdb' to the end of the file path if it isn't there yet. 857 if not search(".pdb$", file_path): 858 file_path = file_path + '.pdb' 859 860 # Open the file for writing. 861 file = open_write_file(file_path, compress_type=compress_type, force=force) 862 863 # Write the structures. 864 cdp.structure.write_pdb(file, model_num=model_num)

865

Source Code for Module generic_fns.structure.main