pipe_control.interatomic

1 ############################################################################### 2 # # 3 # Copyright (C) 2012-2015,2017 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, generate_spin_id_unique, return_spin, spin_loop 40 from pipe_control.pipes import check_pipe 41 42

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

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

125 126

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

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

156 157

158 -def create_interatom(spin_id1=None, spin_id2=None, spin1=None, spin2=None, pipe=None, verbose=False):

159 """Create and return the interatomic data container for the two spins. 160 161 @keyword spin_id1: The spin ID string of the first atom. 162 @type spin_id1: str 163 @keyword spin_id2: The spin ID string of the second atom. 164 @type spin_id2: str 165 @keyword spin1: The optional spin container for the first atom. This is for speeding up the interatomic data container creation, if the spin containers are already available in the calling function. 166 @type spin1: str 167 @keyword spin2: The optional spin container for the second atom. This is for speeding up the interatomic data container creation, if the spin containers are already available in the calling function. 168 @type spin2: str 169 @keyword pipe: The data pipe to create the interatomic data container for. This defaults to the current data pipe if not supplied. 170 @type pipe: str or None 171 @keyword verbose: A flag which if True will result printouts. 172 @type verbose: bool 173 @return: The newly created interatomic data container. 174 @rtype: data.interatomic.InteratomContainer instance 175 """ 176 177 # Printout. 178 if verbose: 179 print("Creating an interatomic data container between the spins '%s' and '%s'." % (spin_id1, spin_id2)) 180 181 # The data pipe. 182 if pipe == None: 183 pipe = pipes.cdp_name() 184 185 # Get the data pipe. 186 dp = pipes.get_pipe(pipe) 187 188 # Check that the spin IDs exist. 189 if spin1 == None: 190 spin1 = return_spin(spin_id=spin_id1, pipe=pipe) 191 if spin1 == None: 192 raise RelaxNoSpinError(spin_id1) 193 if spin2 == None: 194 spin2 = return_spin(spin_id=spin_id2, pipe=pipe) 195 if spin2 == None: 196 raise RelaxNoSpinError(spin_id2) 197 198 # Check if the two spin IDs have already been added. 199 for i in range(len(dp.interatomic)): 200 if spin1._hash != spin2._hash and spin1._hash in [dp.interatomic[i]._spin_hash1, dp.interatomic[i]._spin_hash2] and spin2._hash in [dp.interatomic[i]._spin_hash1, dp.interatomic[i]._spin_hash2]: 201 raise RelaxError("The spin pair %s and %s have already been added." % (spin_id1, spin_id2)) 202 203 # Add the data. 204 interatom = dp.interatomic.add_item(spin_id1=spin_id1, spin_id2=spin_id2, spin_hash1=spin1._hash, spin_hash2=spin2._hash) 205 206 # Store the interatom hash in the spin containers. 207 spin1._interatomic_hashes.append(interatom._hash) 208 spin2._interatomic_hashes.append(interatom._hash) 209 210 # Return the interatomic data container. 211 return interatom

212 213

214 -def define_dipole_pair(spin_id1=None, spin_id2=None, spin1=None, spin2=None, pipe=None, direct_bond=False, spin_selection=False, verbose=True):

215 """Set up the magnetic dipole-dipole interaction. 216 217 @keyword spin_id1: The spin identifier string of the first spin of the pair. 218 @type spin_id1: str 219 @keyword spin_id2: The spin identifier string of the second spin of the pair. 220 @type spin_id2: str 221 @keyword spin1: An optional single spin container for the first atom. This is for speeding up the interatomic data container creation, if the spin containers are already available in the calling function. 222 @type spin1: str 223 @keyword spin2: An optional single spin container for the second atom. This is for speeding up the interatomic data container creation, if the spin containers are already available in the calling function. 224 @type spin2: str 225 @param pipe: The data pipe to operate on. Defaults to the current data pipe. 226 @type pipe: str 227 @keyword direct_bond: A flag specifying if the two spins are directly bonded. 228 @type direct_bond: bool 229 @keyword spin_selection: Define the interatomic data container selection based on the spin selection. If either spin is deselected, the interatomic container will also be deselected. Otherwise the container will be selected. 230 @type spin_selection: bool 231 @keyword verbose: A flag which if True will result in printouts of the created interatomoic data containers. 232 @type verbose: bool 233 """ 234 235 # The data pipe. 236 if pipe == None: 237 pipe = pipes.cdp_name() 238 239 # Get the data pipe. 240 dp = pipes.get_pipe(pipe) 241 242 # Initialise data structures for storing spin data. 243 ids = [] 244 spins = [] 245 spin_selections = [] 246 247 # Pre-supplied spins. 248 if spin1 and spin2: 249 # Store the IDs for the printout. 250 ids.append([spin_id1, spin_id2]) 251 252 # Store the spin data. 253 spins.append([spin1, spin2]) 254 spin_selections.append([spin1.select, spin2.select]) 255 256 # Use the structural data to find connected atoms. 257 elif hasattr(dp, 'structure'): 258 # The selection objects. 259 selection1 = cdp.structure.selection(atom_id=spin_id1) 260 selection2 = cdp.structure.selection(atom_id=spin_id2) 261 262 # Loop over the atoms of the first spin selection. 263 for mol_name1, res_num1, res_name1, atom_num1, atom_name1, mol_index1, atom_index1 in dp.structure.atom_loop(selection=selection1, model_num=1, mol_name_flag=True, res_num_flag=True, res_name_flag=True, atom_num_flag=True, atom_name_flag=True, mol_index_flag=True, index_flag=True): 264 # Generate the first spin ID. 265 id1 = generate_spin_id_unique(pipe_cont=dp, mol_name=mol_name1, res_num=res_num1, res_name=res_name1, spin_num=atom_num1, spin_name=atom_name1) 266 267 # Do the spin exist? 268 spin1 = return_spin(spin_id=id1) 269 if not spin1: 270 continue 271 272 # Loop over the atoms of the second spin selection. 273 for mol_name2, res_num2, res_name2, atom_num2, atom_name2, mol_index2, atom_index2 in dp.structure.atom_loop(selection=selection2, model_num=1, mol_name_flag=True, res_num_flag=True, res_name_flag=True, atom_num_flag=True, atom_name_flag=True, mol_index_flag=True, index_flag=True): 274 # Directly bonded atoms. 275 if direct_bond: 276 # Different molecules. 277 if mol_name1 != mol_name2: 278 continue 279 280 # Skip non-bonded atom pairs. 281 if not dp.structure.are_bonded_index(mol_index1=mol_index1, atom_index1=atom_index1, mol_index2=mol_index2, atom_index2=atom_index2): 282 continue 283 284 # Generate the second spin ID. 285 id2 = generate_spin_id_unique(pipe_cont=dp, mol_name=mol_name2, res_num=res_num2, res_name=res_name2, spin_num=atom_num2, spin_name=atom_name2) 286 287 # Do the spin exist? 288 spin2 = return_spin(spin_id=id2) 289 if not spin2: 290 continue 291 292 # Store the IDs for the printout. 293 ids.append([id1, id2]) 294 295 # Store the spin data. 296 spins.append([spin1, spin2]) 297 spin_selections.append([spin1.select, spin2.select]) 298 299 # No structural data present or the spin IDs are not in the structural data, so use spin loops and some basic rules. 300 if ids == []: 301 for spin1, mol_name1, res_num1, res_name1, id1 in spin_loop(spin_id1, pipe=pipe, full_info=True, return_id=True): 302 for spin2, mol_name2, res_num2, res_name2, id2 in spin_loop(spin_id2, pipe=pipe, full_info=True, return_id=True): 303 # Directly bonded atoms. 304 if direct_bond: 305 # Different molecules. 306 if mol_name1 != mol_name2: 307 continue 308 309 # No element info. 310 if not hasattr(spin1, 'element'): 311 raise RelaxError("The spin '%s' does not have the element type set." % id1) 312 if not hasattr(spin2, 'element'): 313 raise RelaxError("The spin '%s' does not have the element type set." % id2) 314 315 # Backbone NH and CH pairs. 316 pair = False 317 if (spin1.element == 'N' and spin2.element == 'H') or (spin2.element == 'N' and spin1.element == 'H'): 318 pair = True 319 elif (spin1.element == 'C' and spin2.element == 'H') or (spin2.element == 'C' and spin1.element == 'H'): 320 pair = True 321 322 # Same residue, so skip. 323 if pair and res_num1 != None and res_num1 != res_num2: 324 continue 325 elif pair and res_num1 == None and res_name1 != res_name2: 326 continue 327 328 # Store the IDs for the printout. 329 ids.append([id1, id2]) 330 331 # Store the spin data. 332 spins.append([spin1, spin2]) 333 spin_selections.append([spin1.select, spin2.select]) 334 335 # No matches, so fail! 336 if not len(ids): 337 # Find the problem. 338 count1 = 0 339 count2 = 0 340 for spin in spin_loop(spin_id1): 341 count1 += 1 342 for spin in spin_loop(spin_id2): 343 count2 += 1 344 345 # Report the problem. 346 if count1 == 0 and count2 == 0: 347 raise RelaxError("Neither spin IDs '%s' and '%s' match any spins." % (spin_id1, spin_id2)) 348 elif count1 == 0: 349 raise RelaxError("The spin ID '%s' matches no spins." % spin_id1) 350 elif count2 == 0: 351 raise RelaxError("The spin ID '%s' matches no spins." % spin_id2) 352 else: 353 raise RelaxError("Unknown error.") 354 355 # Define the interaction. 356 for i in range(len(ids)): 357 # Unpack. 358 id1, id2 = ids[i] 359 spin1, spin2 = spins[i] 360 361 # Get the interatomic data object, if it exists. 362 interatom = return_interatom(spin_hash1=spin1._hash, spin_hash2=spin2._hash, pipe=pipe) 363 364 # Create the container if needed. 365 if interatom == None: 366 interatom = create_interatom(spin_id1=id1, spin_id2=id2, spin1=spins[i][0], spin2=spins[i][1], pipe=pipe) 367 368 # Check that this has not already been set up. 369 if interatom.dipole_pair: 370 raise RelaxError("The magnetic dipole-dipole interaction already exists between the spins '%s' and '%s'." % (id1, id2)) 371 372 # Set a flag indicating that a dipole-dipole interaction is present. 373 interatom.dipole_pair = True 374 375 # Set the selection. 376 if spin_selection: 377 interatom.select = False 378 if spin_selections[i][0] and spin_selections[i][1]: 379 interatom.select = True 380 381 # Printout. 382 if verbose: 383 # Conversion. 384 for i in range(len(ids)): 385 ids[i][0] = repr(ids[i][0]) 386 ids[i][1] = repr(ids[i][1]) 387 388 # The printout. 389 print("Interatomic interactions are now defined for the following spins:\n") 390 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2"], data=ids)

391 392

393 -def exists_data(pipe=None):

394 """Determine if any interatomic data exists. 395 396 @keyword pipe: The data pipe in which the interatomic data will be checked for. 397 @type pipe: str 398 @return: The answer to the question about the existence of data. 399 @rtype: bool 400 """ 401 402 # The current data pipe. 403 if pipe == None: 404 pipe = pipes.cdp_name() 405 406 # Test the data pipe. 407 check_pipe(pipe) 408 409 # Get the data pipe. 410 dp = pipes.get_pipe(pipe) 411 412 # The interatomic data structure is empty. 413 if dp.interatomic.is_empty(): 414 return False 415 416 # Otherwise. 417 return True

418 419

420 -def hash_update(interatom=None, pipe=None):

421 """Recreate the spin hashes for the interatomic data container. 422 423 @keyword interatom: The interatomic data container. 424 @type interatom: InteratomContainer instance 425 @keyword pipe: The data pipe containing the interatomic data container. Defaults to the current data pipe. 426 @type pipe: str or None 427 """ 428 429 # Fetch the spin containers. 430 spin1 = return_spin(spin_hash=interatom._spin_hash1, pipe=pipe) 431 spin2 = return_spin(spin_hash=interatom._spin_hash2, pipe=pipe) 432 433 # Reset the hashes. 434 interatom._spin_hash1 = spin1._hash 435 interatom._spin_hash2 = spin2._hash

436 437

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

439 """Generator function for looping over all the interatomic data containers. 440 441 @keyword selection1: The optional spin ID selection of the first atom. 442 @type selection1: str 443 @keyword selection2: The optional spin ID selection of the second atom. 444 @type selection2: str 445 @keyword pipe: The data pipe containing the spin. Defaults to the current data pipe. 446 @type pipe: str 447 @keyword skip_desel: A flag which if True will cause only selected interatomic data containers to be returned. 448 @type skip_desel: bool 449 """ 450 451 # The data pipe. 452 if pipe == None: 453 pipe = pipes.cdp_name() 454 455 # Get the data pipe. 456 dp = pipes.get_pipe(pipe) 457 458 # Parse the spin ID selection strings. 459 select_obj = None 460 select_obj1 = None 461 select_obj2 = None 462 if selection1 and selection2: 463 select_obj1 = Selection(selection1) 464 select_obj2 = Selection(selection2) 465 elif selection1: 466 select_obj = Selection(selection1) 467 elif selection2: 468 select_obj = Selection(selection2) 469 470 # Loop over the containers, yielding them. 471 for i in range(len(dp.interatomic)): 472 # Skip deselected containers. 473 if skip_desel and not dp.interatomic[i].select: 474 continue 475 476 # Aliases. 477 interatom = dp.interatomic[i] 478 mol_index1, res_index1, spin_index1 = dp.mol._spin_hash_lookup[interatom._spin_hash1] 479 mol_index2, res_index2, spin_index2 = dp.mol._spin_hash_lookup[interatom._spin_hash2] 480 mol1 = dp.mol[mol_index1] 481 res1 = dp.mol[mol_index1].res[res_index1] 482 spin1 = dp.mol[mol_index1].res[res_index1].spin[spin_index1] 483 mol2 = dp.mol[mol_index2] 484 res2 = dp.mol[mol_index2].res[res_index2] 485 spin2 = dp.mol[mol_index2].res[res_index2].spin[spin_index2] 486 487 # The different selection combinations. 488 if select_obj: 489 sel1 = select_obj.contains_spin(spin_name=spin1.name, spin_num=spin1.num, res_name=res1.name, res_num=res1.num, mol=mol1.name) 490 sel2 = select_obj.contains_spin(spin_name=spin2.name, spin_num=spin2.num, res_name=res2.name, res_num=res2.num, mol=mol2.name) 491 if select_obj1: 492 sel11 = select_obj1.contains_spin(spin_name=spin1.name, spin_num=spin1.num, res_name=res1.name, res_num=res1.num, mol=mol1.name) 493 sel12 = select_obj1.contains_spin(spin_name=spin2.name, spin_num=spin2.num, res_name=res2.name, res_num=res2.num, mol=mol2.name) 494 if select_obj2: 495 sel21 = select_obj2.contains_spin(spin_name=spin1.name, spin_num=spin1.num, res_name=res1.name, res_num=res1.num, mol=mol1.name) 496 sel22 = select_obj2.contains_spin(spin_name=spin2.name, spin_num=spin2.num, res_name=res2.name, res_num=res2.num, mol=mol2.name) 497 498 # Check that the selections are met. 499 if select_obj: 500 if not sel1 and not sel2: 501 continue 502 if select_obj1: 503 if not (sel11 or sel12): 504 continue 505 if select_obj2: 506 if not (sel21 or sel22): 507 continue 508 509 # Return the container. 510 yield interatom

511 512

513 -def metadata_update(interatom_index=None, pipe=None):

514 """Update all of the private look up metadata. 515 516 @keyword interatom_index: The index of the interatomic data container to update. If not supplied, all containers will be updated. 517 @type interatom_index: int or None 518 @keyword pipe: The data pipe to update, defaulting to the current data pipe. 519 @type pipe: str or None 520 """ 521 522 # The data pipe. 523 if pipe == None: 524 pipe = pipes.cdp_name() 525 526 # Test the data pipe. 527 check_pipe(pipe) 528 529 # Get the data pipe. 530 dp = pipes.get_pipe(pipe) 531 532 # Loop over the containers. 533 for i in range(len(dp.interatomic)): 534 # Interatom skipping. 535 if interatom_index != None and interatom_index != i: 536 continue 537 538 # Alias. 539 interatom = dp.interatomic[i] 540 541 # Get the matching spin from the IDs, rather than from the hashes. 542 spin1 = return_spin(spin_id=interatom.spin_id1, pipe=pipe) 543 spin2 = return_spin(spin_id=interatom.spin_id2, pipe=pipe) 544 545 # Update the hashes. 546 interatom._spin_hash1 = None 547 interatom._spin_hash2 = None 548 if spin1: 549 interatom._spin_hash1 = spin1._hash 550 if spin2: 551 interatom._spin_hash2 = spin2._hash

552 553

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

555 """Set up the magnetic dipole-dipole interaction. 556 557 @keyword file: The name of the file to open. 558 @type file: str 559 @keyword dir: The directory containing the file (defaults to the current directory if None). 560 @type dir: str or None 561 @keyword unit: The measurement unit. This can be either 'meter' or 'Angstrom'. 562 @type unit: str 563 @keyword spin_id1_col: The column containing the spin ID strings of the first spin. 564 @type spin_id1_col: int 565 @keyword spin_id2_col: The column containing the spin ID strings of the second spin. 566 @type spin_id2_col: int 567 @keyword data_col: The column containing the averaged distances in meters. 568 @type data_col: int or None 569 @keyword sep: The column separator which, if None, defaults to whitespace. 570 @type sep: str or None 571 """ 572 573 # Check the units. 574 if unit not in ['meter', 'Angstrom']: 575 raise RelaxError("The measurement unit of '%s' must be one of 'meter' or 'Angstrom'." % unit) 576 577 # Test if the current data pipe exists. 578 check_pipe() 579 580 # Test if sequence data exists. 581 if not exists_mol_res_spin_data(): 582 raise RelaxNoSequenceError 583 584 # Extract the data from the file, and clean it up. 585 file_data = extract_data(file, dir, sep=sep) 586 file_data = strip(file_data, comments=True) 587 588 # Loop over the RDC data. 589 data = [] 590 for line in file_data: 591 # Invalid columns. 592 if spin_id1_col > len(line): 593 warn(RelaxWarning("The data %s is invalid, no first spin ID column can be found." % line)) 594 continue 595 if spin_id2_col > len(line): 596 warn(RelaxWarning("The data %s is invalid, no second spin ID column can be found." % line)) 597 continue 598 if data_col and data_col > len(line): 599 warn(RelaxWarning("The data %s is invalid, no data column can be found." % line)) 600 continue 601 602 # Unpack. 603 spin_id1 = line[spin_id1_col-1] 604 spin_id2 = line[spin_id2_col-1] 605 ave_dist = None 606 if data_col: 607 ave_dist = line[data_col-1] 608 609 # Convert and check the value. 610 if ave_dist != None: 611 try: 612 ave_dist = float(ave_dist) 613 except ValueError: 614 warn(RelaxWarning("The averaged distance of '%s' from the line %s is invalid." % (ave_dist, line))) 615 continue 616 617 # Unit conversion. 618 if unit == 'Angstrom': 619 ave_dist = ave_dist * 1e-10 620 621 # Get the interatomic data container. 622 spin1 = return_spin(spin_id=spin_id1) 623 spin2 = return_spin(spin_id=spin_id2) 624 interatom = return_interatom(spin_hash1=spin1._hash, spin_hash2=spin2._hash) 625 626 # No container found, so create it. 627 if interatom == None: 628 interatom = create_interatom(spin_id1=spin_id1, spin_id2=spin_id2, verbose=True) 629 630 # Store the averaged distance. 631 interatom.r = ave_dist 632 633 # Store the data for the printout. 634 data.append([repr(interatom.spin_id1), repr(interatom.spin_id2), repr(ave_dist)]) 635 636 # No data, so fail! 637 if not len(data): 638 raise RelaxError("No data could be extracted from the file.") 639 640 # Print out. 641 print("The following averaged distances have been read:\n") 642 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2", "Ave_distance(meters)"], data=data)

643 644

645 -def return_interatom(spin_hash1=None, spin_hash2=None, pipe=None):

646 """Return the list of interatomic data containers for the two spins. 647 648 @keyword spin_hash1: The unique spin hash for the first atom. 649 @type spin_hash1: str 650 @keyword spin_hash2: The unique spin hash for the second atom. 651 @type spin_hash2: str 652 @keyword pipe: The data pipe holding the container. Defaults to the current data pipe. 653 @type pipe: str or None 654 @return: The matching interatomic data container, if it exists. 655 @rtype: data.interatomic.InteratomContainer instance or None 656 """ 657 658 # Checks. 659 if spin_hash1 == None: 660 raise RelaxError("The first spin hash must be supplied.") 661 if spin_hash2 == None: 662 raise RelaxError("The second spin hash must be supplied.") 663 664 # The data pipe. 665 if pipe == None: 666 pipe = pipes.cdp_name() 667 668 # Get the data pipe. 669 dp = pipes.get_pipe(pipe) 670 671 # Return the matching container. 672 for i in range(len(dp.interatomic)): 673 if spin_hash1 != spin_hash2 and spin_hash1 in [dp.interatomic[i]._spin_hash1, dp.interatomic[i]._spin_hash2] and spin_hash2 in [dp.interatomic[i]._spin_hash1, dp.interatomic[i]._spin_hash2]: 674 return dp.interatomic[i] 675 676 # No matchs. 677 return None

678 679

680 -def return_interatom_list(spin_hash=None, pipe=None):

681 """Return the list of interatomic data containers for the given spin. 682 683 @keyword spin_hash: The unique spin hash. 684 @type spin_hash: str 685 @keyword pipe: The data pipe holding the container. This defaults to the current data pipe. 686 @type pipe: str or None 687 @return: The list of matching interatomic data containers, if any exist. 688 @rtype: list of data.interatomic.InteratomContainer instances 689 """ 690 691 # Check. 692 if spin_hash == None: 693 raise RelaxError("The unique spin hash must be supplied.") 694 695 # The data pipe. 696 if pipe == None: 697 pipe = pipes.cdp_name() 698 699 # Get the data pipe. 700 dp = pipes.get_pipe(pipe) 701 702 # Initialise. 703 interatoms = [] 704 705 # Find and append all containers. 706 for i in range(len(dp.interatomic)): 707 if spin_hash in [dp.interatomic[i]._spin_hash1, dp.interatomic[i]._spin_hash2]: 708 interatoms.append(dp.interatomic[i]) 709 710 # Return the list of containers. 711 return interatoms

712 713

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

715 """Set up the magnetic dipole-dipole interaction. 716 717 @keyword spin_id1: The spin identifier string of the first spin of the pair. 718 @type spin_id1: str 719 @keyword spin_id2: The spin identifier string of the second spin of the pair. 720 @type spin_id2: str 721 @keyword ave_dist: The r^-3 averaged interatomic distance. 722 @type ave_dist: float 723 @keyword unit: The measurement unit. This can be either 'meter' or 'Angstrom'. 724 @type unit: str 725 """ 726 727 # Check the units. 728 if unit not in ['meter', 'Angstrom']: 729 raise RelaxError("The measurement unit of '%s' must be one of 'meter' or 'Angstrom'." % unit) 730 731 # Unit conversion. 732 if unit == 'Angstrom': 733 ave_dist = ave_dist * 1e-10 734 735 # Generate the selection objects. 736 sel_obj1 = Selection(spin_id1) 737 sel_obj2 = Selection(spin_id2) 738 739 # Loop over the interatomic containers. 740 data = [] 741 for interatom in interatomic_loop(): 742 # Get the spin info. 743 mol_name1, res_num1, res_name1, spin1 = return_spin(spin_hash=interatom._spin_hash1, full_info=True) 744 mol_name2, res_num2, res_name2, spin2 = return_spin(spin_hash=interatom._spin_hash2, full_info=True) 745 746 # No match, either way. 747 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)): 748 continue 749 750 # Store the averaged distance. 751 interatom.r = ave_dist 752 753 # Store the data for the printout. 754 data.append([repr(interatom.spin_id1), repr(interatom.spin_id2), repr(ave_dist)]) 755 756 # No data, so fail! 757 if not len(data): 758 raise RelaxError("No data could be set.") 759 760 # Print out. 761 print("The following averaged distances have been set:\n") 762 write_data(out=sys.stdout, headings=["Spin_ID_1", "Spin_ID_2", "Ave_distance(meters)"], data=data)

763 764

765 -def unit_vectors(ave=True):

766 """Extract the bond vectors from the loaded structures and store them in the spin container. 767 768 @keyword ave: A flag which if True will cause the average of all vectors to be calculated. 769 @type ave: bool 770 """ 771 772 # Test if the current data pipe exists. 773 check_pipe() 774 775 # Test if interatomic data exists. 776 if not exists_data(): 777 raise RelaxNoInteratomError 778 779 # Print out. 780 if ave: 781 print("Averaging all vectors.") 782 else: 783 print("No averaging of the vectors.") 784 785 # Loop over the interatomic data containers. 786 no_vectors = True 787 pos_info = False 788 for interatom in interatomic_loop(skip_desel=False): 789 # Get the spin info. 790 spin1 = return_spin(spin_hash=interatom._spin_hash1) 791 spin2 = return_spin(spin_hash=interatom._spin_hash2) 792 793 # No positional information. 794 if not hasattr(spin1, 'pos'): 795 continue 796 if not hasattr(spin2, 'pos'): 797 continue 798 799 # Positional information flag. 800 pos_info = True 801 802 # Both single positions. 803 if is_float(spin1.pos[0], raise_error=False) and is_float(spin2.pos[0], raise_error=False): 804 # The vector. 805 vector_list = [spin2.pos - spin1.pos] 806 807 # A single and multiple position pair. 808 elif is_float(spin1.pos[0], raise_error=False) or is_float(spin2.pos[0], raise_error=False): 809 # The first spin has multiple positions. 810 if is_float(spin2.pos[0], raise_error=False): 811 vector_list = [] 812 for i in range(len(spin1.pos)): 813 vector_list.append(spin2.pos - spin1.pos[i]) 814 815 # The second spin has multiple positions. 816 else: 817 vector_list = [] 818 for i in range(len(spin2.pos)): 819 vector_list.append(spin2.pos[i] - spin1.pos) 820 821 # Both spins have multiple positions. 822 else: 823 # Non-matching number of positions. 824 if len(spin1.pos) != len(spin2.pos): 825 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))) 826 827 # Calculate all vectors. 828 vector_list = [] 829 for i in range(len(spin1.pos)): 830 # No structural information. 831 if spin1.pos[i] is None or spin2.pos[i] is None: 832 warn(RelaxWarning("No structural information for state %i can be found between spins '%s' and '%s'." % (i, interatom.spin_id1, interatom.spin_id2))) 833 vector_list.append(None) 834 835 # All data is present. 836 else: 837 vector_list.append(spin2.pos[i] - spin1.pos[i]) 838 839 # Unit vectors. 840 for i in range(len(vector_list)): 841 # No vector. 842 if vector_list[i] is None: 843 continue 844 845 # Normalisation factor. 846 norm_factor = norm(vector_list[i]) 847 848 # Test for zero length. 849 if norm_factor == 0.0: 850 warn(RelaxZeroVectorWarning(spin_id1=interatom.spin_id1, spin_id2=interatom.spin_id2)) 851 852 # Calculate the normalised vector. 853 else: 854 vector_list[i] = vector_list[i] / norm_factor 855 856 # Average. 857 if ave: 858 ave_vector = zeros(3, float64) 859 count = 0 860 for i in range(len(vector_list)): 861 if vector_list[i] is not None: 862 ave_vector = ave_vector + vector_list[i] 863 count += 1 864 vector_list = [ave_vector / count] 865 866 # Convert to a single vector if needed. 867 if len(vector_list) == 1: 868 vector_list = vector_list[0] 869 870 # Store the unit vector(s). 871 setattr(interatom, 'vector', vector_list) 872 873 # We have a vector! 874 no_vectors = False 875 876 # Printout. 877 num = 1 878 if not is_float(vector_list[0], raise_error=False): 879 num = len(vector_list) 880 plural = 's' 881 if num == 1: 882 plural = '' 883 if spin1.name: 884 spin1_str = spin1.name 885 else: 886 spin1_str = spin1.num 887 if spin2.name: 888 spin2_str = spin2.name 889 else: 890 spin2_str = spin2.num 891 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)) 892 893 # Catch the problem of no positional information being present. 894 if not pos_info: 895 raise RelaxError("Positional information could not be found for any spins.") 896 897 # Right, catch the problem of missing vectors to prevent massive user confusion! 898 if no_vectors: 899 raise RelaxError("No vectors could be extracted.")

900

Source Code for Module pipe_control.interatomic