pipe_control.interatomic

1 ############################################################################### 2 # # 3 # Copyright (C) 2012-2013 Edward d'Auvergne # 4 # # 5 # This file is part of the program relax (http://www.nmr-relax.com). # 6 # # 7 # This program 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 3 of the License, or # 10 # (at your option) any later version. # 11 # # 12 # This program 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 this program. If not, see <http://www.gnu.org/licenses/>. # 19 # # 20 ############################################################################### 21 22 # Module docstring. 23 """Module for the manipulation of the interatomic data structures in the relax data store.""" 24 25 # Python module imports. 26 from copy import deepcopy 27 from numpy import float64, zeros 28 from numpy.linalg import norm 29 from re import search 30 import sys 31 from warnings import warn 32 33 # relax module imports. 34 from lib.arg_check import is_float 35 from lib.errors import RelaxError, RelaxInteratomInconsistentError, RelaxNoInteratomError, RelaxNoSpinError 36 from lib.io import extract_data, strip, write_data 37 from lib.warnings import RelaxWarning, RelaxZeroVectorWarning 38 from pipe_control import pipes 39 from pipe_control.mol_res_spin import Selection, count_spins, exists_mol_res_spin_data, return_spin, spin_loop 40 41

42 -def copy(pipe_from=None, pipe_to=None, spin_id1=None, spin_id2=None, verbose=True):

43 """Copy the interatomic data from one data pipe to another. 44 45 @keyword pipe_from: The data pipe to copy the interatomic data from. This defaults to the current data pipe. 46 @type pipe_from: str 47 @keyword pipe_to: The data pipe to copy the interatomic data to. This defaults to the current data pipe. 48 @type pipe_to: str 49 @keyword spin_id1: The spin ID string of the first atom. 50 @type spin_id1: str 51 @keyword spin_id2: The spin ID string of the second atom. 52 @type spin_id2: str 53 @keyword verbose: A flag which if True will cause info about each spin pair to be printed out. 54 @type verbose: bool 55 """ 56 57 # Defaults. 58 if pipe_from == None and pipe_to == None: 59 raise RelaxError("The pipe_from and pipe_to arguments cannot both be set to None.") 60 elif pipe_from == None: 61 pipe_from = pipes.cdp_name() 62 elif pipe_to == None: 63 pipe_to = pipes.cdp_name() 64 65 # Test if the pipe_from and pipe_to data pipes exist. 66 pipes.test(pipe_from) 67 pipes.test(pipe_to) 68 69 # Check that the spin IDs exist. 70 if spin_id1: 71 if count_spins(selection=spin_id1, pipe=pipe_from, skip_desel=False) == 0: 72 raise RelaxNoSpinError(spin_id1, pipe_from) 73 if count_spins(selection=spin_id1, pipe=pipe_to, skip_desel=False) == 0: 74 raise RelaxNoSpinError(spin_id1, pipe_to) 75 if spin_id2: 76 if count_spins(selection=spin_id2, pipe=pipe_from, skip_desel=False) == 0: 77 raise RelaxNoSpinError(spin_id2, pipe_from) 78 if count_spins(selection=spin_id2, pipe=pipe_to, skip_desel=False) == 0: 79 raise RelaxNoSpinError(spin_id2, pipe_to) 80 81 # Check for the sequence data in the target pipe if no spin IDs are given. 82 if not spin_id1 and not spin_id2: 83 for spin, spin_id in spin_loop(pipe=pipe_from, return_id=True): 84 if not return_spin(spin_id, pipe=pipe_to): 85 raise RelaxNoSpinError(spin_id, pipe_to) 86 87 # Test if pipe_from contains interatomic data (skipping the rest of the function if it is missing). 88 if not exists_data(pipe_from): 89 return 90 91 # Loop over the interatomic data of the pipe_from data pipe. 92 ids = [] 93 for interatom in interatomic_loop(selection1=spin_id1, selection2=spin_id2, pipe=pipe_from): 94 # Create a new container. 95 new_interatom = create_interatom(spin_id1=interatom.spin_id1, spin_id2=interatom.spin_id2, pipe=pipe_to) 96 97 # Duplicate all the objects of the container. 98 for name in dir(interatom): 99 # Skip special objects. 100 if search('^_', name): 101 continue 102 103 # Skip the spin IDs. 104 if name in ['spin_id1', 'spin_id2']: 105 continue 106 107 # Skip class methods. 108 if name in list(interatom.__class__.__dict__.keys()): 109 continue 110 111 # Duplicate all other objects. 112 obj = deepcopy(getattr(interatom, name)) 113 setattr(new_interatom, name, obj) 114 115 # Store the IDs for the printout. 116 ids.append([repr(interatom.spin_id1), repr(interatom.spin_id2)]) 117 118 # Print out. 119 if verbose: 120 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2"], data=ids)

121 122

123 -def consistent_interatomic_data(pipe1=None, pipe2=None):

124 """Check that the interatomic data is consistent between two data pipes. 125 126 @keyword pipe1: The name of the first data pipe to compare. 127 @type pipe1: str 128 @keyword pipe2: The name of the second data pipe to compare. 129 @type pipe2: str 130 @raises RelaxError: If the data is inconsistent. 131 """ 132 133 # Get the data pipes. 134 dp1 = pipes.get_pipe(pipe1) 135 dp2 = pipes.get_pipe(pipe2) 136 137 # Check the data lengths. 138 if len(dp1.interatomic) != len(dp2.interatomic): 139 raise RelaxInteratomInconsistentError(pipe1, pipe2) 140 141 # Loop over the interatomic data. 142 for i in range(len(dp1.interatomic)): 143 # Alias the containers. 144 interatom1 = dp1.interatomic[i] 145 interatom2 = dp2.interatomic[i] 146 147 # Check the spin IDs. 148 if interatom1.spin_id1 != interatom2.spin_id1: 149 raise RelaxInteratomInconsistentError(pipe1, pipe2) 150 if interatom1.spin_id2 != interatom2.spin_id2: 151 raise RelaxInteratomInconsistentError(pipe1, pipe2)

152 153

154 -def create_interatom(spin_id1=None, spin_id2=None, pipe=None, verbose=False):

155 """Create and return the interatomic data container for the two spins. 156 157 @keyword spin_id1: The spin ID string of the first atom. 158 @type spin_id1: str 159 @keyword spin_id2: The spin ID string of the second atom. 160 @type spin_id2: str 161 @keyword pipe: The data pipe to create the interatomic data container for. This defaults to the current data pipe if not supplied. 162 @type pipe: str or None 163 @keyword verbose: A flag which if True will result printouts. 164 @type verbose: bool 165 @return: The newly created interatomic data container. 166 @rtype: data.interatomic.InteratomContainer instance 167 """ 168 169 # Printout. 170 if verbose: 171 print("Creating an interatomic data container between the spins '%s' and '%s'." % (spin_id1, spin_id2)) 172 173 # The data pipe. 174 if pipe == None: 175 pipe = pipes.cdp_name() 176 177 # Get the data pipe. 178 dp = pipes.get_pipe(pipe) 179 180 # Check that the spin IDs exist. 181 spin = return_spin(spin_id1, pipe) 182 if spin == None: 183 raise RelaxNoSpinError(spin_id1) 184 spin = return_spin(spin_id2, pipe) 185 if spin == None: 186 raise RelaxNoSpinError(spin_id2) 187 188 # Check if the two spin IDs have already been added. 189 for i in range(len(dp.interatomic)): 190 if id_match(spin_id=spin_id1, interatom=dp.interatomic[i], pipe=pipe) and id_match(spin_id=spin_id2, interatom=dp.interatomic[i], pipe=pipe): 191 raise RelaxError("The spin pair %s and %s have already been added." % (spin_id1, spin_id2)) 192 193 # Add and return the data. 194 return dp.interatomic.add_item(spin_id1=spin_id1, spin_id2=spin_id2)

195 196

197 -def define(spin_id1=None, spin_id2=None, pipe=None, direct_bond=False, verbose=True):

198 """Set up the magnetic dipole-dipole interaction. 199 200 @keyword spin_id1: The spin identifier string of the first spin of the pair. 201 @type spin_id1: str 202 @keyword spin_id2: The spin identifier string of the second spin of the pair. 203 @type spin_id2: str 204 @param pipe: The data pipe to operate on. Defaults to the current data pipe. 205 @type pipe: str 206 @keyword direct_bond: A flag specifying if the two spins are directly bonded. 207 @type direct_bond: bool 208 @keyword verbose: A flag which if True will result in printouts of the created interatomoic data containers. 209 @type verbose: bool 210 """ 211 212 # The data pipe. 213 if pipe == None: 214 pipe = pipes.cdp_name() 215 216 # Get the data pipe. 217 dp = pipes.get_pipe(pipe) 218 219 # Loop over both spin selections. 220 ids = [] 221 for spin1, mol_name1, res_num1, res_name1, id1 in spin_loop(spin_id1, pipe=pipe, full_info=True, return_id=True): 222 for spin2, mol_name2, res_num2, res_name2, id2 in spin_loop(spin_id2, pipe=pipe, full_info=True, return_id=True): 223 # Directly bonded atoms. 224 if direct_bond: 225 # Different molecules. 226 if mol_name1 != mol_name2: 227 continue 228 229 # From structural info. 230 if hasattr(dp, 'structure') and dp.structure.get_molecule(mol_name1, model=1): 231 if not dp.structure.are_bonded(atom_id1=id1, atom_id2=id2): 232 continue 233 234 # From the residue info. 235 else: 236 # No element info. 237 if not hasattr(spin1, 'element'): 238 raise RelaxError("The spin '%s' does not have the element type set." % id1) 239 if not hasattr(spin2, 'element'): 240 raise RelaxError("The spin '%s' does not have the element type set." % id2) 241 242 # Backbone NH and CH pairs. 243 pair = False 244 if (spin1.element == 'N' and spin2.element == 'H') or (spin2.element == 'N' and spin1.element == 'H'): 245 pair = True 246 elif (spin1.element == 'C' and spin2.element == 'H') or (spin2.element == 'C' and spin1.element == 'H'): 247 pair = True 248 249 # Same residue, so skip. 250 if pair and res_num1 != None and res_num1 != res_num2: 251 continue 252 elif pair and res_num1 == None and res_name1 != res_name2: 253 continue 254 255 # Get the interatomic data object, if it exists. 256 interatom = return_interatom(id1, id2, pipe=pipe) 257 258 # Create the container if needed. 259 if interatom == None: 260 interatom = create_interatom(spin_id1=id1, spin_id2=id2, pipe=pipe) 261 262 # Check that this has not already been set up. 263 if interatom.dipole_pair: 264 raise RelaxError("The magnetic dipole-dipole interaction already exists between the spins '%s' and '%s'." % (id1, id2)) 265 266 # Set a flag indicating that a dipole-dipole interaction is present. 267 interatom.dipole_pair = True 268 269 # Store the IDs for the printout. 270 ids.append([repr(id1), repr(id2)]) 271 272 # No matches, so fail! 273 if not len(ids): 274 # Find the problem. 275 count1 = 0 276 count2 = 0 277 for spin in spin_loop(spin_id1): 278 count1 += 1 279 for spin in spin_loop(spin_id2): 280 count2 += 1 281 282 # Report the problem. 283 if count1 == 0 and count2 == 0: 284 raise RelaxError("Neither spin IDs '%s' and '%s' match any spins." % (spin_id1, spin_id2)) 285 elif count1 == 0: 286 raise RelaxError("The spin ID '%s' matches no spins." % spin_id1) 287 elif count2 == 0: 288 raise RelaxError("The spin ID '%s' matches no spins." % spin_id2) 289 else: 290 raise RelaxError("Unknown error.") 291 292 # Print out. 293 if verbose: 294 print("Interatomic interactions are now defined for the following spins:\n") 295 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2"], data=ids)

296 297

298 -def exists_data(pipe=None):

299 """Determine if any interatomic data exists. 300 301 @keyword pipe: The data pipe in which the interatomic data will be checked for. 302 @type pipe: str 303 @return: The answer to the question about the existence of data. 304 @rtype: bool 305 """ 306 307 # The current data pipe. 308 if pipe == None: 309 pipe = pipes.cdp_name() 310 311 # Test the data pipe. 312 pipes.test(pipe) 313 314 # Get the data pipe. 315 dp = pipes.get_pipe(pipe) 316 317 # The interatomic data structure is empty. 318 if dp.interatomic.is_empty(): 319 return False 320 321 # Otherwise. 322 return True

323 324

325 -def id_match(spin_id=None, interatom=None, pipe=None):

326 """Test if the spin ID matches one of the two spins of the given container. 327 328 @keyword spin_id: The spin ID string of the first atom. 329 @type spin_id: str 330 @keyword interatom: The interatomic data container. 331 @type interatom: InteratomContainer instance 332 @keyword pipe: The data pipe containing the interatomic data container. Defaults to the current data pipe. 333 @type pipe: str or None 334 @return: True if the spin ID matches one of the two spins, False otherwise. 335 @rtype: bool 336 """ 337 338 # Get the spin containers. 339 spin1 = return_spin(interatom.spin_id1, pipe=pipe) 340 spin2 = return_spin(interatom.spin_id2, pipe=pipe) 341 342 # No spins. 343 if spin1 == None or spin2 == None: 344 return False 345 346 # Check if the ID is in the private metadata list. 347 if spin_id in spin1._spin_ids or spin_id in spin2._spin_ids: 348 return True 349 350 # Nothing found. 351 return False

352 353

354 -def interatomic_loop(selection1=None, selection2=None, pipe=None, skip_desel=True):

355 """Generator function for looping over all the interatomic data containers. 356 357 @keyword selection1: The optional spin ID selection of the first atom. 358 @type selection1: str 359 @keyword selection2: The optional spin ID selection of the second atom. 360 @type selection2: str 361 @keyword pipe: The data pipe containing the spin. Defaults to the current data pipe. 362 @type pipe: str 363 @keyword skip_desel: A flag which if True will cause only selected interatomic data containers to be returned. 364 @type skip_desel: bool 365 """ 366 367 # The data pipe. 368 if pipe == None: 369 pipe = pipes.cdp_name() 370 371 # Get the data pipe. 372 dp = pipes.get_pipe(pipe) 373 374 # Parse the spin ID selection strings. 375 select_obj = None 376 select_obj1 = None 377 select_obj2 = None 378 if selection1 and selection2: 379 select_obj1 = Selection(selection1) 380 select_obj2 = Selection(selection2) 381 elif selection1: 382 select_obj = Selection(selection1) 383 elif selection2: 384 select_obj = Selection(selection2) 385 386 # Loop over the containers, yielding them. 387 for i in range(len(dp.interatomic)): 388 # Skip deselected containers. 389 if skip_desel and not dp.interatomic[i].select: 390 continue 391 392 # Aliases. 393 interatom = dp.interatomic[i] 394 mol_index1, res_index1, spin_index1 = cdp.mol._spin_id_lookup[interatom.spin_id1] 395 mol_index2, res_index2, spin_index2 = cdp.mol._spin_id_lookup[interatom.spin_id2] 396 mol1 = cdp.mol[mol_index1] 397 res1 = cdp.mol[mol_index1].res[res_index1] 398 spin1 = cdp.mol[mol_index1].res[res_index1].spin[spin_index1] 399 mol2 = cdp.mol[mol_index2] 400 res2 = cdp.mol[mol_index2].res[res_index2] 401 spin2 = cdp.mol[mol_index2].res[res_index2].spin[spin_index2] 402 403 # The different selection combinations. 404 if select_obj: 405 sel1 = select_obj.contains_spin(spin_name=spin1.name, spin_num=spin1.num, res_name=res1.name, res_num=res1.num, mol=mol1.name) 406 sel2 = select_obj.contains_spin(spin_name=spin2.name, spin_num=spin2.num, res_name=res2.name, res_num=res2.num, mol=mol2.name) 407 if select_obj1: 408 sel11 = select_obj1.contains_spin(spin_name=spin1.name, spin_num=spin1.num, res_name=res1.name, res_num=res1.num, mol=mol1.name) 409 sel12 = select_obj1.contains_spin(spin_name=spin2.name, spin_num=spin2.num, res_name=res2.name, res_num=res2.num, mol=mol2.name) 410 if select_obj2: 411 sel21 = select_obj2.contains_spin(spin_name=spin1.name, spin_num=spin1.num, res_name=res1.name, res_num=res1.num, mol=mol1.name) 412 sel22 = select_obj2.contains_spin(spin_name=spin2.name, spin_num=spin2.num, res_name=res2.name, res_num=res2.num, mol=mol2.name) 413 414 # Check that the selections are met. 415 if select_obj: 416 if not sel1 and not sel2: 417 continue 418 if select_obj1: 419 if not (sel11 or sel12): 420 continue 421 if select_obj2: 422 if not (sel21 or sel22): 423 continue 424 425 # Return the container. 426 yield interatom

427 428

429 -def read_dist(file=None, dir=None, unit='meter', spin_id1_col=None, spin_id2_col=None, data_col=None, sep=None):

430 """Set up the magnetic dipole-dipole interaction. 431 432 @keyword file: The name of the file to open. 433 @type file: str 434 @keyword dir: The directory containing the file (defaults to the current directory if None). 435 @type dir: str or None 436 @keyword unit: The measurement unit. This can be either 'meter' or 'Angstrom'. 437 @type unit: str 438 @keyword spin_id1_col: The column containing the spin ID strings of the first spin. 439 @type spin_id1_col: int 440 @keyword spin_id2_col: The column containing the spin ID strings of the second spin. 441 @type spin_id2_col: int 442 @keyword data_col: The column containing the averaged distances in meters. 443 @type data_col: int or None 444 @keyword sep: The column separator which, if None, defaults to whitespace. 445 @type sep: str or None 446 """ 447 448 # Check the units. 449 if unit not in ['meter', 'Angstrom']: 450 raise RelaxError("The measurement unit of '%s' must be one of 'meter' or 'Angstrom'." % unit) 451 452 # Test if the current data pipe exists. 453 pipes.test() 454 455 # Test if sequence data exists. 456 if not exists_mol_res_spin_data(): 457 raise RelaxNoSequenceError 458 459 # Extract the data from the file, and clean it up. 460 file_data = extract_data(file, dir, sep=sep) 461 file_data = strip(file_data, comments=True) 462 463 # Loop over the RDC data. 464 data = [] 465 for line in file_data: 466 # Invalid columns. 467 if spin_id1_col > len(line): 468 warn(RelaxWarning("The data %s is invalid, no first spin ID column can be found." % line)) 469 continue 470 if spin_id2_col > len(line): 471 warn(RelaxWarning("The data %s is invalid, no second spin ID column can be found." % line)) 472 continue 473 if data_col and data_col > len(line): 474 warn(RelaxWarning("The data %s is invalid, no data column can be found." % line)) 475 continue 476 477 # Unpack. 478 spin_id1 = line[spin_id1_col-1] 479 spin_id2 = line[spin_id2_col-1] 480 ave_dist = None 481 if data_col: 482 ave_dist = line[data_col-1] 483 484 # Convert and check the value. 485 if ave_dist != None: 486 try: 487 ave_dist = float(ave_dist) 488 except ValueError: 489 warn(RelaxWarning("The averaged distance of '%s' from the line %s is invalid." % (ave_dist, line))) 490 continue 491 492 # Unit conversion. 493 if unit == 'Angstrom': 494 ave_dist = ave_dist * 1e-10 495 496 # Get the interatomic data container. 497 interatom = return_interatom(spin_id1, spin_id2) 498 499 # No container found, so create it. 500 if interatom == None: 501 interatom = create_interatom(spin_id1=spin_id1, spin_id2=spin_id2, verbose=True) 502 503 # Store the averaged distance. 504 interatom.r = ave_dist 505 506 # Store the data for the printout. 507 data.append([repr(interatom.spin_id1), repr(interatom.spin_id2), repr(ave_dist)]) 508 509 # No data, so fail! 510 if not len(data): 511 raise RelaxError("No data could be extracted from the file.") 512 513 # Print out. 514 print("The following averaged distances have been read:\n") 515 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2", "Ave_distance(meters)"], data=data)

516 517

518 -def return_interatom(spin_id1=None, spin_id2=None, pipe=None):

519 """Return the list of interatomic data containers for the two spins. 520 521 @keyword spin_id1: The spin ID string of the first atom. 522 @type spin_id1: str 523 @keyword spin_id2: The spin ID string of the second atom. 524 @type spin_id2: str 525 @keyword pipe: The data pipe holding the container. Defaults to the current data pipe. 526 @type pipe: str or None 527 @return: The matching interatomic data container, if it exists. 528 @rtype: data.interatomic.InteratomContainer instance or None 529 """ 530 531 # Checks. 532 if spin_id1 == None: 533 raise RelaxError("The first spin ID must be supplied.") 534 if spin_id2 == None: 535 raise RelaxError("The second spin ID must be supplied.") 536 537 # The data pipe. 538 if pipe == None: 539 pipe = pipes.cdp_name() 540 541 # Get the data pipe. 542 dp = pipes.get_pipe(pipe) 543 544 # Return the matching container. 545 for i in range(len(dp.interatomic)): 546 if id_match(spin_id=spin_id1, interatom=dp.interatomic[i], pipe=pipe) and id_match(spin_id=spin_id2, interatom=dp.interatomic[i], pipe=pipe): 547 return dp.interatomic[i] 548 549 # No matchs. 550 return None

551 552

553 -def return_interatom_list(spin_id=None, pipe=None):

554 """Return the list of interatomic data containers for the given spin. 555 556 @keyword spin_id: The spin ID string. 557 @type spin_id: str 558 @keyword pipe: The data pipe holding the container. This defaults to the current data pipe. 559 @type pipe: str or None 560 @return: The list of matching interatomic data containers, if any exist. 561 @rtype: list of data.interatomic.InteratomContainer instances 562 """ 563 564 # Check. 565 if spin_id == None: 566 raise RelaxError("The spin ID must be supplied.") 567 568 # The data pipe. 569 if pipe == None: 570 pipe = pipes.cdp_name() 571 572 # Get the data pipe. 573 dp = pipes.get_pipe(pipe) 574 575 # Initialise. 576 interatoms = [] 577 578 # Find and append all containers. 579 for i in range(len(dp.interatomic)): 580 if id_match(spin_id=spin_id, interatom=dp.interatomic[i], pipe=pipe) or id_match(spin_id=spin_id, interatom=dp.interatomic[i], pipe=pipe): 581 interatoms.append(dp.interatomic[i]) 582 583 # Return the list of containers. 584 return interatoms

585 586

587 -def set_dist(spin_id1=None, spin_id2=None, ave_dist=None, unit='meter'):

588 """Set up the magnetic dipole-dipole interaction. 589 590 @keyword spin_id1: The spin identifier string of the first spin of the pair. 591 @type spin_id1: str 592 @keyword spin_id2: The spin identifier string of the second spin of the pair. 593 @type spin_id2: str 594 @keyword ave_dist: The r^-3 averaged interatomic distance. 595 @type ave_dist: float 596 @keyword unit: The measurement unit. This can be either 'meter' or 'Angstrom'. 597 @type unit: str 598 """ 599 600 # Check the units. 601 if unit not in ['meter', 'Angstrom']: 602 raise RelaxError("The measurement unit of '%s' must be one of 'meter' or 'Angstrom'." % unit) 603 604 # Unit conversion. 605 if unit == 'Angstrom': 606 ave_dist = ave_dist * 1e-10 607 608 # Generate the selection objects. 609 sel_obj1 = Selection(spin_id1) 610 sel_obj2 = Selection(spin_id2) 611 612 # Loop over the interatomic containers. 613 data = [] 614 for interatom in interatomic_loop(): 615 # Get the spin info. 616 mol_name1, res_num1, res_name1, spin1 = return_spin(interatom.spin_id1, full_info=True) 617 mol_name2, res_num2, res_name2, spin2 = return_spin(interatom.spin_id2, full_info=True) 618 619 # No match, either way. 620 if not (sel_obj1.contains_spin(spin_num=spin1.num, spin_name=spin1.name, res_num=res_num1, res_name=res_name1, mol=mol_name1) and sel_obj2.contains_spin(spin_num=spin2.num, spin_name=spin2.name, res_num=res_num2, res_name=res_name2, mol=mol_name2)) and not (sel_obj2.contains_spin(spin_num=spin1.num, spin_name=spin1.name, res_num=res_num1, res_name=res_name1, mol=mol_name1) and sel_obj1.contains_spin(spin_num=spin2.num, spin_name=spin2.name, res_num=res_num2, res_name=res_name2, mol=mol_name2)): 621 continue 622 623 # Store the averaged distance. 624 interatom.r = ave_dist 625 626 # Store the data for the printout. 627 data.append([repr(interatom.spin_id1), repr(interatom.spin_id2), repr(ave_dist)]) 628 629 # No data, so fail! 630 if not len(data): 631 raise RelaxError("No data could be set.") 632 633 # Print out. 634 print("The following averaged distances have been set:\n") 635 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2", "Ave_distance(meters)"], data=data)

636 637

638 -def unit_vectors(ave=True):

639 """Extract the bond vectors from the loaded structures and store them in the spin container. 640 641 @keyword ave: A flag which if True will cause the average of all vectors to be calculated. 642 @type ave: bool 643 """ 644 645 # Test if the current data pipe exists. 646 pipes.test() 647 648 # Test if interatomic data exists. 649 if not exists_data(): 650 raise RelaxNoInteratomError 651 652 # Print out. 653 if ave: 654 print("Averaging all vectors.") 655 else: 656 print("No averaging of the vectors.") 657 658 # Loop over the interatomic data containers. 659 no_vectors = True 660 pos_info = False 661 for interatom in interatomic_loop(skip_desel=False): 662 # Get the spin info. 663 spin1 = return_spin(interatom.spin_id1) 664 spin2 = return_spin(interatom.spin_id2) 665 666 # No positional information. 667 if not hasattr(spin1, 'pos'): 668 continue 669 if not hasattr(spin2, 'pos'): 670 continue 671 672 # Positional information flag. 673 pos_info = True 674 675 # Both single positions. 676 if is_float(spin1.pos[0], raise_error=False) and is_float(spin2.pos[0], raise_error=False): 677 # The vector. 678 vector_list = [spin2.pos - spin1.pos] 679 680 # A single and multiple position pair. 681 elif is_float(spin1.pos[0], raise_error=False) or is_float(spin2.pos[0], raise_error=False): 682 # The first spin has multiple positions. 683 if is_float(spin2.pos[0], raise_error=False): 684 vector_list = [] 685 for i in range(len(spin1.pos)): 686 vector_list.append(spin2.pos - spin1.pos[i]) 687 688 # The second spin has multiple positions. 689 else: 690 vector_list = [] 691 for i in range(len(spin2.pos)): 692 vector_list.append(spin2.pos[i] - spin1.pos) 693 694 # Both spins have multiple positions. 695 else: 696 # Non-matching number of positions. 697 if len(spin1.pos) != len(spin2.pos): 698 raise RelaxError("The spin '%s' consists of %s positions whereas the spin '%s' consists of %s - these numbers must match." % (interatom.spin_id1, len(spin1.pos), interatom.spin_id1, len(spin1.pos))) 699 700 # Calculate all vectors. 701 vector_list = [] 702 for i in range(len(spin1.pos)): 703 vector_list.append(spin2.pos[i] - spin1.pos[i]) 704 705 # Unit vectors. 706 for i in range(len(vector_list)): 707 # Normalisation factor. 708 norm_factor = norm(vector_list[i]) 709 710 # Test for zero length. 711 if norm_factor == 0.0: 712 warn(RelaxZeroVectorWarning(spin_id1=interatom.spin_id1, spin_id2=interatom.spin_id2)) 713 714 # Calculate the normalised vector. 715 else: 716 vector_list[i] = vector_list[i] / norm_factor 717 718 # Average. 719 if ave: 720 ave_vector = zeros(3, float64) 721 for i in range(len(vector_list)): 722 ave_vector = ave_vector + vector_list[i] 723 vector_list = [ave_vector / len(vector_list)] 724 725 # Convert to a single vector if needed. 726 if len(vector_list) == 1: 727 vector_list = vector_list[0] 728 729 # Store the unit vector(s). 730 setattr(interatom, 'vector', vector_list) 731 732 # We have a vector! 733 no_vectors = False 734 735 # Print out. 736 num = 1 737 if not is_float(vector_list[0], raise_error=False): 738 num = len(vector_list) 739 plural = 's' 740 if num == 1: 741 plural = '' 742 if spin1.name: 743 spin1_str = spin1.name 744 else: 745 spin1_str = spin1.num 746 if spin2.name: 747 spin2_str = spin2.name 748 else: 749 spin2_str = spin2.num 750 print("Calculated %s %s-%s unit vector%s between the spins '%s' and '%s'." % (num, spin1_str, spin2_str, plural, interatom.spin_id1, interatom.spin_id2)) 751 752 # Catch the problem of no positional information being present. 753 if not pos_info: 754 raise RelaxError("Positional information could not be found for any spins.") 755 756 # Right, catch the problem of missing vectors to prevent massive user confusion! 757 if no_vectors: 758 raise RelaxError("No vectors could be extracted.")

759

Source Code for Module pipe_control.interatomic