generic_fns.rdc

1 ############################################################################### 2 # # 3 # Copyright (C) 2003-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 RDC data.""" 24 25 # Python module imports. 26 from copy import deepcopy 27 from math import pi, sqrt 28 from numpy import float64, ones, zeros 29 from numpy.linalg import norm 30 import sys 31 from warnings import warn 32 33 # relax module imports. 34 from check_types import is_float 35 from float import nan 36 from generic_fns import grace, pipes 37 from generic_fns.align_tensor import get_tensor_index 38 from generic_fns.interatomic import consistent_interatomic_data, create_interatom, interatomic_loop, return_interatom 39 from generic_fns.mol_res_spin import exists_mol_res_spin_data, return_spin 40 from maths_fns.rdc import ave_rdc_tensor 41 from physical_constants import dipolar_constant, return_gyromagnetic_ratio 42 from relax_errors import RelaxError, RelaxNoAlignError, RelaxNoRDCError, RelaxNoSequenceError, RelaxSpinTypeError 43 from relax_io import extract_data, open_write_file, write_data 44 from relax_warnings import RelaxWarning 45 46

47 -def back_calc(align_id=None):

48 """Back calculate the RDC from the alignment tensor and unit bond vectors. 49 50 @keyword align_id: The alignment tensor ID string. 51 @type align_id: str 52 """ 53 54 # Check the pipe setup. 55 check_pipe_setup(rdc_id=align_id, sequence=True, N=True, tensors=True) 56 57 # Convert the align IDs to an array, or take all IDs. 58 if align_id: 59 align_ids = [align_id] 60 else: 61 align_ids = cdp.align_ids 62 63 # Add the ID to the RDC IDs, if needed. 64 for align_id in align_ids: 65 # Init. 66 if not hasattr(cdp, 'rdc_ids'): 67 cdp.rdc_ids = [] 68 69 # Add the ID. 70 if align_id not in cdp.rdc_ids: 71 cdp.rdc_ids.append(align_id) 72 73 # The weights. 74 weights = ones(cdp.N, float64) / cdp.N 75 76 # Unit vector data structure init. 77 unit_vect = zeros((cdp.N, 3), float64) 78 79 # Loop over the interatomic data. 80 count = 0 81 for interatom in interatomic_loop(): 82 # Skip containers with no interatomic vectors. 83 if not hasattr(interatom, 'vector'): 84 continue 85 86 # Get the spins. 87 spin1 = return_spin(interatom.spin_id1) 88 spin2 = return_spin(interatom.spin_id2) 89 90 # Checks. 91 if not hasattr(spin1, 'isotope'): 92 raise RelaxSpinTypeError(interatom.spin_id1) 93 if not hasattr(spin2, 'isotope'): 94 raise RelaxSpinTypeError(interatom.spin_id2) 95 96 # Single vector. 97 if is_float(interatom.vector[0]): 98 vectors = [interatom.vector] 99 else: 100 vectors = interatom.vector 101 102 # Gyromagnetic ratios. 103 g1 = return_gyromagnetic_ratio(spin1.isotope) 104 g2 = return_gyromagnetic_ratio(spin2.isotope) 105 106 # Calculate the RDC dipolar constant (in Hertz, and the 3 comes from the alignment tensor), and append it to the list. 107 dj = 3.0/(2.0*pi) * dipolar_constant(g1, g2, interatom.r) 108 109 # Unit vectors. 110 for c in range(cdp.N): 111 unit_vect[c] = vectors[c] / norm(vectors[c]) 112 113 # Initialise if necessary. 114 if not hasattr(interatom, 'rdc_bc'): 115 interatom.rdc_bc = {} 116 117 # Calculate the RDCs. 118 for id in align_ids: 119 # The signed value. 120 interatom.rdc_bc[id] = ave_rdc_tensor(dj, unit_vect, cdp.N, cdp.align_tensors[get_tensor_index(id)].A, weights=weights) 121 122 # The absolute value. 123 if hasattr(interatom, 'absolute_rdc') and id in interatom.absolute_rdc.keys() and interatom.absolute_rdc[id]: 124 interatom.rdc_bc[id] = abs(interatom.rdc_bc[id]) 125 126 # Increment the counter. 127 count += 1 128 129 # No RDCs calculated. 130 if not count: 131 warn(RelaxWarning("No RDCs have been back calculated, probably due to missing bond vector information."))

132 133

134 -def check_pipe_setup(pipe=None, rdc_id=None, sequence=False, N=False, tensors=False, rdc=False):

135 """Check that the current data pipe has been setup sufficiently. 136 137 @keyword pipe: The data pipe to check, defaulting to the current pipe. 138 @type pipe: None or str 139 @keyword rdc_id: The RDC ID string to check for in cdp.rdc_ids. 140 @type rdc_id: None or str 141 @keyword sequence: A flag which when True will invoke the sequence data check. 142 @type sequence: bool 143 @keyword N: A flag which if True will check that cdp.N is set. 144 @type N: bool 145 @keyword tensors: A flag which if True will check that alignment tensors exist. 146 @type tensors: bool 147 @keyword rdc: A flag which if True will check that RDCs exist. 148 @type rdc: bool 149 """ 150 151 # The data pipe. 152 if pipe == None: 153 pipe = pipes.cdp_name() 154 155 # Get the data pipe. 156 dp = pipes.get_pipe(pipe) 157 158 # Test if the current data pipe exists. 159 pipes.test(pipe) 160 161 # Test if sequence data exists. 162 if sequence and not exists_mol_res_spin_data(pipe): 163 raise RelaxNoSequenceError(pipe) 164 165 # Check for dp.N. 166 if N and not hasattr(dp, 'N'): 167 raise RelaxError("The number of states N has not been set.") 168 169 # Check that alignment tensors are present. 170 if tensors and (not hasattr(dp, 'align_tensors') or len(dp.align_tensors) == 0): 171 raise RelaxNoTensorError('alignment') 172 173 # Test for the alignment ID. 174 if rdc_id and (not hasattr(dp, 'align_ids') or rdc_id not in dp.align_ids): 175 raise RelaxNoAlignError(rdc_id, pipe) 176 177 # Test if RDC data exists. 178 if rdc and not hasattr(dp, 'align_ids'): 179 raise RelaxNoAlignError() 180 if rdc and not hasattr(dp, 'rdc_ids'): 181 raise RelaxNoRDCError() 182 elif rdc and rdc_id and rdc_id not in dp.rdc_ids: 183 raise RelaxNoRDCError(rdc_id)

184 185

186 -def convert(value, align_id, to_intern=False):

187 """Convert the RDC based on the 'D' or '2D' data type. 188 189 @param value: The value or error to convert. 190 @type value: float or None 191 @param align_id: The alignment tensor ID string. 192 @type align_id: str 193 @keyword to_intern: A flag which if True will convert to the internal D notation if needed, or if False will convert from the internal D notation to the external D or 2D format. 194 @type to_intern: bool 195 @return: The converted value. 196 @rtype: float or None 197 """ 198 199 # Handle values of None. 200 if value == None: 201 return None 202 203 # The conversion factor. 204 factor = 1.0 205 if hasattr(cdp, 'rdc_data_types') and align_id in cdp.rdc_data_types and cdp.rdc_data_types[align_id] == '2D': 206 # Convert from 2D to D. 207 if to_intern: 208 factor = 2.0 209 210 # Convert from D to 2D. 211 else: 212 factor = 0.5 213 214 # Return the converted value. 215 return value * factor

216 217

218 -def copy(pipe_from=None, pipe_to=None, align_id=None):

219 """Copy the RDC data from one data pipe to another. 220 221 @keyword pipe_from: The data pipe to copy the RDC data from. This defaults to the current data pipe. 222 @type pipe_from: str 223 @keyword pipe_to: The data pipe to copy the RDC data to. This defaults to the current data pipe. 224 @type pipe_to: str 225 @param align_id: The alignment ID string. 226 @type align_id: str 227 """ 228 229 # Defaults. 230 if pipe_from == None and pipe_to == None: 231 raise RelaxError("The pipe_from and pipe_to arguments cannot both be set to None.") 232 elif pipe_from == None: 233 pipe_from = pipes.cdp_name() 234 elif pipe_to == None: 235 pipe_to = pipes.cdp_name() 236 237 # Check the pipe setup. 238 check_pipe_setup(pipe=pipe_from, rdc_id=align_id, sequence=True, rdc=True) 239 check_pipe_setup(pipe=pipe_to, sequence=True) 240 241 # Get the data pipes. 242 dp_from = pipes.get_pipe(pipe_from) 243 dp_to = pipes.get_pipe(pipe_to) 244 245 # Test that the interatomic data is consistent between the two data pipe. 246 consistent_interatomic_data(pipe1=pipe_to, pipe2=pipe_from) 247 248 # The IDs. 249 if align_id == None: 250 align_ids = dp_from.align_ids 251 else: 252 align_ids = [align_id] 253 254 # Init target pipe global structures. 255 if not hasattr(dp_to, 'align_ids'): 256 dp_to.align_ids = [] 257 if not hasattr(dp_to, 'rdc_ids'): 258 dp_to.rdc_ids = [] 259 260 # Loop over the align IDs. 261 for align_id in align_ids: 262 # Copy the global data. 263 if align_id not in dp_to.align_ids and align_id not in dp_to.align_ids: 264 dp_to.align_ids.append(align_id) 265 if align_id in dp_from.rdc_ids and align_id not in dp_to.rdc_ids: 266 dp_to.rdc_ids.append(align_id) 267 268 # Loop over the interatomic data. 269 for i in range(len(dp_from.interatomic)): 270 # Alias the containers. 271 interatom_from = dp_from.interatomic[i] 272 interatom_to = dp_to.interatomic[i] 273 274 # No data or errors. 275 if (not hasattr(interatom_from, 'rdc') or not align_id in interatom_from.rdc) and (not hasattr(interatom_from, 'rdc_err') or not align_id in interatom_from.rdc_err): 276 continue 277 278 # Initialise the data structures if necessary. 279 if hasattr(interatom_from, 'rdc') and not hasattr(interatom_to, 'rdc'): 280 interatom_to.rdc = {} 281 if hasattr(interatom_from, 'rdc_err') and not hasattr(interatom_to, 'rdc_err'): 282 interatom_to.rdc_err = {} 283 284 # Copy the value and error from pipe_from. 285 if hasattr(interatom_from, 'rdc'): 286 interatom_to.rdc[align_id] = interatom_from.rdc[align_id] 287 if hasattr(interatom_from, 'rdc_err'): 288 interatom_to.rdc_err[align_id] = interatom_from.rdc_err[align_id]

289 290

291 -def corr_plot(format=None, file=None, dir=None, force=False):

292 """Generate a correlation plot of the measured vs. back-calculated RDCs. 293 294 @keyword format: The format for the plot file. The following values are accepted: 'grace', a Grace plot; None, a plain text file. 295 @type format: str or None 296 @keyword file: The file name or object to write to. 297 @type file: str or file object 298 @keyword dir: The name of the directory to place the file into (defaults to the current directory). 299 @type dir: str 300 @keyword force: A flag which if True will cause any pre-existing file to be overwritten. 301 @type force: bool 302 """ 303 304 # Check the pipe setup. 305 check_pipe_setup(sequence=True) 306 307 # Does RDC data exist? 308 if not hasattr(cdp, 'rdc_ids') or not cdp.rdc_ids: 309 warn(RelaxWarning("No RDC data exists, skipping file creation.")) 310 return 311 312 # Open the file for writing. 313 file = open_write_file(file, dir, force) 314 315 # Init. 316 data = [] 317 318 # The diagonal. 319 data.append([[-100, -100, 0], [100, 100, 0]]) 320 321 # Loop over the RDC data. 322 for align_id in cdp.rdc_ids: 323 # Append a new list for this alignment. 324 data.append([]) 325 326 # Errors present? 327 err_flag = False 328 for interatom in interatomic_loop(): 329 # Error present. 330 if hasattr(interatom, 'rdc_err') and align_id in interatom.rdc_err.keys(): 331 err_flag = True 332 break 333 334 # Loop over the interatomic data. 335 for interatom in interatomic_loop(): 336 # Skip if data is missing. 337 if not hasattr(interatom, 'rdc') or not hasattr(interatom, 'rdc_bc') or not align_id in interatom.rdc.keys() or not align_id in interatom.rdc_bc.keys(): 338 continue 339 340 # Append the data. 341 data[-1].append([convert(interatom.rdc_bc[align_id], align_id), convert(interatom.rdc[align_id], align_id)]) 342 343 # Errors. 344 if err_flag: 345 if hasattr(interatom, 'rdc_err') and align_id in interatom.rdc_err.keys(): 346 data[-1][-1].append(convert(interatom.rdc_err[align_id], align_id)) 347 else: 348 data[-1][-1].append(None) 349 350 # Label. 351 data[-1][-1].append("%s-%s" % (interatom.spin_id1, interatom.spin_id2)) 352 353 # The data size. 354 size = len(data) 355 356 # Only one data set. 357 data = [data] 358 359 # Graph type. 360 if err_flag: 361 graph_type = 'xydy' 362 else: 363 graph_type = 'xy' 364 365 # Grace file. 366 if format == 'grace': 367 # The header. 368 grace.write_xy_header(file=file, title="RDC correlation plot", sets=size, set_names=[None]+cdp.rdc_ids, linestyle=[2]+[0]*size, data_type=['rdc_bc', 'rdc'], axis_min=[-10, -10], axis_max=[10, 10], legend_pos=[1, 0.5]) 369 370 # The main data. 371 grace.write_xy_data(data=data, file=file, graph_type=graph_type)

372 373

374 -def delete(align_id=None):

375 """Delete the RDC data corresponding to the alignment ID. 376 377 @keyword align_id: The alignment tensor ID string. If not specified, all data will be deleted. 378 @type align_id: str or None 379 """ 380 381 # Check the pipe setup. 382 check_pipe_setup(sequence=True, rdc_id=align_id, rdc=True) 383 384 # The IDs. 385 if not align_id: 386 ids = deepcopy(cdp.rdc_ids) 387 else: 388 ids = [align_id] 389 390 # Loop over the alignments, removing all the corresponding data. 391 for id in ids: 392 # The RDC ID. 393 cdp.rdc_ids.pop(cdp.rdc_ids.index(id)) 394 395 # The data type. 396 if hasattr(cdp, 'rdc_data_types') and id in cdp.rdc_data_types: 397 cdp.rdc_data_types.pop(id) 398 399 # The interatomic data. 400 for interatom in interatomic_loop(): 401 # The data. 402 if hasattr(interatom, 'rdc') and id in interatom.rdc: 403 interatom.rdc.pop(id) 404 405 # The error. 406 if hasattr(interatom, 'rdc_err') and id in interatom.rdc_err: 407 interatom.rdc_err.pop(id) 408 409 # Clean the global data. 410 if not hasattr(cdp, 'pcs_ids') or id not in cdp.pcs_ids: 411 cdp.align_ids.pop(cdp.align_ids.index(id))

412 413

414 -def display(align_id=None, bc=False):

415 """Display the RDC data corresponding to the alignment ID. 416 417 @keyword align_id: The alignment tensor ID string. 418 @type align_id: str 419 @keyword bc: The back-calculation flag which if True will cause the back-calculated rather than measured data to be displayed. 420 @type bc: bool 421 """ 422 423 # Check the pipe setup. 424 check_pipe_setup(sequence=True, rdc_id=align_id, rdc=True) 425 426 # Call the write method with sys.stdout as the file. 427 write(align_id=align_id, file=sys.stdout, bc=bc)

428 429

430 -def q_factors(spin_id=None):

431 """Calculate the Q-factors for the RDC data. 432 433 @keyword spin_id: The spin ID string used to restrict the Q-factor calculation to a subset of all spins. 434 @type spin_id: None or str 435 """ 436 437 # Check the pipe setup. 438 check_pipe_setup(sequence=True) 439 440 # No RDCs, so no Q factors can be calculated. 441 if not hasattr(cdp, 'rdc_ids') or not len(cdp.rdc_ids): 442 warn(RelaxWarning("No RDC data exists, Q factors cannot be calculated.")) 443 return 444 445 # Q-factor dictonaries. 446 cdp.q_factors_rdc = {} 447 cdp.q_factors_rdc_norm2 = {} 448 449 # Loop over the alignments. 450 for align_id in cdp.rdc_ids: 451 # Init. 452 D2_sum = 0.0 453 sse = 0.0 454 455 # Interatomic data loop. 456 dj = None 457 N = 0 458 interatom_count = 0 459 rdc_data = False 460 rdc_bc_data = False 461 norm2_flag = True 462 for interatom in interatomic_loop(): 463 # Increment the counter. 464 interatom_count += 1 465 466 # Data checks. 467 if hasattr(interatom, 'rdc') and align_id in interatom.rdc: 468 rdc_data = True 469 if hasattr(interatom, 'rdc_bc') and align_id in interatom.rdc_bc: 470 rdc_bc_data = True 471 472 # Skip containers without RDC data. 473 if not hasattr(interatom, 'rdc') or not hasattr(interatom, 'rdc_bc') or not align_id in interatom.rdc or interatom.rdc[align_id] == None or not align_id in interatom.rdc_bc or interatom.rdc_bc[align_id] == None: 474 continue 475 476 # Get the spins. 477 spin1 = return_spin(interatom.spin_id1) 478 spin2 = return_spin(interatom.spin_id2) 479 480 # Sum of squares. 481 sse = sse + (interatom.rdc[align_id] - interatom.rdc_bc[align_id])**2 482 483 # Sum the RDCs squared (for one type of normalisation). 484 D2_sum = D2_sum + interatom.rdc[align_id]**2 485 486 # Gyromagnetic ratios. 487 g1 = return_gyromagnetic_ratio(spin1.isotope) 488 g2 = return_gyromagnetic_ratio(spin2.isotope) 489 490 # Calculate the RDC dipolar constant (in Hertz, and the 3 comes from the alignment tensor), and append it to the list. 491 dj_new = 3.0/(2.0*pi) * dipolar_constant(g1, g2, interatom.r) 492 if norm2_flag and dj != None and dj_new != dj: 493 warn(RelaxWarning("The dipolar constant is not the same for all RDCs, skipping the Q factor normalised with 2Da^2(4 + 3R)/5.")) 494 norm2_flag = False 495 else: 496 dj = dj_new 497 498 # Increment the number of data sets. 499 N = N + 1 500 501 # Warnings (and then exit). 502 if not interatom_count: 503 warn(RelaxWarning("No interatomic data containers have been used in the calculation, skipping the RDC Q factor calculation.")) 504 return 505 if not rdc_data: 506 warn(RelaxWarning("No RDC data can be found for the alignment ID '%s', skipping the RDC Q factor calculation for this alignment." % align_id)) 507 continue 508 if not rdc_bc_data: 509 warn(RelaxWarning("No back-calculated RDC data can be found for the alignment ID '%s', skipping the RDC Q factor calculation for this alignment." % align_id)) 510 continue 511 512 # Normalisation factor of 2Da^2(4 + 3R)/5. 513 if norm2_flag: 514 D = dj * cdp.align_tensors[cdp.align_ids.index(align_id)].A_diag 515 Da = 1.0/3.0 * (D[2, 2] - (D[0, 0]+D[1, 1])/2.0) 516 Dr = 1.0/3.0 * (D[0, 0] - D[1, 1]) 517 if Da == 0: 518 R = nan 519 else: 520 R = Dr / Da 521 norm = 2.0 * (Da)**2 * (4.0 + 3.0*R**2)/5.0 522 if Da == 0.0: 523 norm = 1e-15 524 525 # The Q-factor for the alignment. 526 cdp.q_factors_rdc[align_id] = sqrt(sse / N / norm) 527 528 else: 529 cdp.q_factors_rdc[align_id] = 0.0 530 531 # The second Q-factor definition. 532 cdp.q_factors_rdc_norm2[align_id] = sqrt(sse / D2_sum) 533 534 # The total Q-factor. 535 cdp.q_rdc = 0.0 536 cdp.q_rdc_norm2 = 0.0 537 for id in cdp.q_factors_rdc: 538 cdp.q_rdc = cdp.q_rdc + cdp.q_factors_rdc[id]**2 539 for id in cdp.q_factors_rdc_norm2: 540 cdp.q_rdc_norm2 = cdp.q_rdc_norm2 + cdp.q_factors_rdc_norm2[id]**2 541 cdp.q_rdc = sqrt(cdp.q_rdc / len(cdp.q_factors_rdc)) 542 cdp.q_rdc_norm2 = sqrt(cdp.q_rdc_norm2 / len(cdp.q_factors_rdc_norm2))

543 544

545 -def read(align_id=None, file=None, dir=None, file_data=None, data_type='D', spin_id1_col=None, spin_id2_col=None, data_col=None, error_col=None, sep=None, neg_g_corr=False, absolute=False):

546 """Read the RDC data from file. 547 548 @keyword align_id: The alignment tensor ID string. 549 @type align_id: str 550 @keyword file: The name of the file to open. 551 @type file: str 552 @keyword dir: The directory containing the file (defaults to the current directory if None). 553 @type dir: str or None 554 @keyword file_data: An alternative to opening a file, if the data already exists in the correct format. The format is a list of lists where the first index corresponds to the row and the second the column. 555 @type file_data: list of lists 556 @keyword data_type: A string which is set to 'D' means that the splitting in the aligned sample was assumed to be J + D, or if set to '2D' then the splitting was taken as J + 2D. 557 @keyword spin_id1_col: The column containing the spin ID strings of the first spin. 558 @type spin_id1_col: int 559 @keyword spin_id2_col: The column containing the spin ID strings of the second spin. 560 @type spin_id2_col: int 561 @keyword data_col: The column containing the RDC data in Hz. 562 @type data_col: int or None 563 @keyword error_col: The column containing the RDC errors. 564 @type error_col: int or None 565 @keyword sep: The column separator which, if None, defaults to whitespace. 566 @type sep: str or None 567 @keyword neg_g_corr: A flag which is used to correct for the negative gyromagnetic ratio of 15N. If True, a sign inversion will be applied to all RDC values to be loaded. 568 @type neg_g_corr: bool 569 @keyword absolute: A flag which if True indicates that the RDCs to load are signless. All RDCs will then be converted to positive values. 570 @type absolute: bool 571 """ 572 573 # Check the pipe setup. 574 check_pipe_setup(sequence=True) 575 576 # Either the data or error column must be supplied. 577 if data_col == None and error_col == None: 578 raise RelaxError("One of either the data or error column must be supplied.") 579 580 # Store the data type as global data (need for the conversion of RDC data). 581 if not hasattr(cdp, 'rdc_data_types'): 582 cdp.rdc_data_types = {} 583 if not align_id in cdp.rdc_data_types: 584 cdp.rdc_data_types[align_id] = data_type 585 586 # Spin specific data. 587 ##################### 588 589 # Extract the data from the file. 590 file_data = extract_data(file, dir, sep=sep) 591 592 # Loop over the RDC data. 593 data = [] 594 for line in file_data: 595 # Invalid columns. 596 if spin_id1_col > len(line): 597 warn(RelaxWarning("The data %s is invalid, no first spin ID column can be found." % line)) 598 continue 599 if spin_id2_col > len(line): 600 warn(RelaxWarning("The data %s is invalid, no second spin ID column can be found." % line)) 601 continue 602 if data_col and data_col > len(line): 603 warn(RelaxWarning("The data %s is invalid, no data column can be found." % line)) 604 continue 605 if error_col and error_col > len(line): 606 warn(RelaxWarning("The data %s is invalid, no error column can be found." % line)) 607 continue 608 609 # Unpack. 610 spin_id1 = line[spin_id1_col-1] 611 spin_id2 = line[spin_id2_col-1] 612 value = None 613 if data_col: 614 value = line[data_col-1] 615 error = None 616 if error_col: 617 error = line[error_col-1] 618 619 # Convert the spin IDs. 620 if spin_id1[0] in ["\"", "\'"]: 621 spin_id1 = eval(spin_id1) 622 if spin_id2[0] in ["\"", "\'"]: 623 spin_id2 = eval(spin_id2) 624 625 # Convert and check the value. 626 if value == 'None': 627 value = None 628 if value != None: 629 try: 630 value = float(value) 631 except ValueError: 632 warn(RelaxWarning("The RDC value of the line %s is invalid." % line)) 633 continue 634 635 # Convert and check the error. 636 if error == 'None': 637 error = None 638 if error != None: 639 try: 640 error = float(error) 641 except ValueError: 642 warn(RelaxWarning("The error value of the line %s is invalid." % line)) 643 continue 644 645 # Get the spins. 646 spin1 = return_spin(spin_id1) 647 spin2 = return_spin(spin_id2) 648 649 # Check the spin IDs. 650 if not spin1: 651 warn(RelaxWarning("The spin ID '%s' cannot be found in the current data pipe, skipping the data %s." % (spin_id1, line))) 652 continue 653 if not spin2: 654 warn(RelaxWarning("The spin ID '%s' cannot be found in the current data pipe, skipping the data %s." % (spin_id2, line))) 655 continue 656 657 # Test the error value (cannot be 0.0). 658 if error == 0.0: 659 raise RelaxError("An invalid error value of zero has been encountered.") 660 661 # Get the interatomic data container. 662 interatom = return_interatom(spin_id1, spin_id2) 663 664 # Create the container if needed. 665 if interatom == None: 666 interatom = create_interatom(spin_id1=spin_id1, spin_id2=spin_id2) 667 668 # Convert and add the data. 669 if data_col: 670 # Data conversion. 671 value = convert(value, align_id, to_intern=True) 672 673 # Correction for the negative gyromagnetic ratio of 15N. 674 if neg_g_corr and value != None: 675 value = -value 676 677 # Absolute values. 678 if absolute: 679 # Force the value to be positive. 680 value = abs(value) 681 682 # Initialise. 683 if not hasattr(interatom, 'rdc'): 684 interatom.rdc = {} 685 686 # Add the value. 687 interatom.rdc[align_id] = value 688 689 # Store the absolute value flag. 690 if not hasattr(interatom, 'absolute_rdc'): 691 interatom.absolute_rdc = {} 692 interatom.absolute_rdc[align_id] = absolute 693 694 # Convert and add the error. 695 if error_col: 696 # Data conversion. 697 error = convert(error, align_id, to_intern=True) 698 699 # Initialise. 700 if not hasattr(interatom, 'rdc_err'): 701 interatom.rdc_err = {} 702 703 # Append the error. 704 interatom.rdc_err[align_id] = error 705 706 # Append the data for printout. 707 data.append([spin_id1, spin_id2]) 708 if is_float(value): 709 data[-1].append("%20.15f" % value) 710 else: 711 data[-1].append("%20s" % value) 712 if is_float(error): 713 data[-1].append("%20.15f" % error) 714 else: 715 data[-1].append("%20s" % error) 716 717 # No data, so fail hard! 718 if not len(data): 719 raise RelaxError("No RDC data could be extracted.") 720 721 # Print out. 722 print("The following RDCs have been loaded into the relax data store:\n") 723 write_data(out=sys.stdout, headings=["Spin_ID1", "Spin_ID2", "Value", "Error"], data=data) 724 725 # Initialise some global structures. 726 if not hasattr(cdp, 'align_ids'): 727 cdp.align_ids = [] 728 if not hasattr(cdp, 'rdc_ids'): 729 cdp.rdc_ids = [] 730 731 # Add the RDC id string. 732 if align_id not in cdp.align_ids: 733 cdp.align_ids.append(align_id) 734 if align_id not in cdp.rdc_ids: 735 cdp.rdc_ids.append(align_id)

736 737

738 -def set_errors(align_id=None, spin_id1=None, spin_id2=None, sd=None):

739 """Set the RDC errors if not already present. 740 741 @keyword align_id: The optional alignment tensor ID string. 742 @type align_id: str 743 @keyword spin_id1: The optional spin ID string of the first spin. 744 @type spin_id1: None or str 745 @keyword spin_id2: The optional spin ID string of the second spin. 746 @type spin_id2: None or str 747 @keyword sd: The RDC standard deviation in Hz. 748 @type sd: float or int. 749 """ 750 751 # Check the pipe setup. 752 check_pipe_setup(sequence=True, rdc_id=align_id, rdc=True) 753 754 # Convert the align IDs to an array, or take all IDs. 755 if align_id: 756 align_ids = [align_id] 757 else: 758 align_ids = cdp.rdc_ids 759 760 # Loop over the interatomic data. 761 for interatom in interatomic_loop(selection1=spin_id1, selection2=spin_id2): 762 # No data structure. 763 if not hasattr(interatom, 'rdc_err'): 764 interatom.rdc_err = {} 765 766 # Set the error. 767 for id in align_ids: 768 interatom.rdc_err[id] = sd

769 770

771 -def weight(align_id=None, spin_id=None, weight=1.0):

772 """Set optimisation weights on the RDC data. 773 774 @keyword align_id: The alignment tensor ID string. 775 @type align_id: str 776 @keyword spin_id: The spin ID string. 777 @type spin_id: None or str 778 @keyword weight: The optimisation weight. The higher the value, the more importance the RDC will have. 779 @type weight: float or int. 780 """ 781 782 # Check the pipe setup. 783 check_pipe_setup(sequence=True, rdc_id=align_id, rdc=True) 784 785 # Loop over the interatomic data. 786 for interatom in interatomic_loop(): 787 # No data structure. 788 if not hasattr(interatom, 'rdc_weight'): 789 interatom.rdc_weight = {} 790 791 # Set the weight. 792 interatom.rdc_weight[align_id] = weight

793 794

795 -def write(align_id=None, file=None, dir=None, bc=False, force=False):

796 """Display the RDC data corresponding to the alignment ID. 797 798 @keyword align_id: The alignment tensor ID string. 799 @type align_id: str 800 @keyword file: The file name or object to write to. 801 @type file: str or file object 802 @keyword dir: The name of the directory to place the file into (defaults to the current directory). 803 @type dir: str 804 @keyword bc: The back-calculation flag which if True will cause the back-calculated rather than measured data to be written. 805 @type bc: bool 806 @keyword force: A flag which if True will cause any pre-existing file to be overwritten. 807 @type force: bool 808 """ 809 810 # Check the pipe setup. 811 check_pipe_setup(sequence=True, rdc_id=align_id, rdc=True) 812 813 # Open the file for writing. 814 file = open_write_file(file, dir, force) 815 816 # Loop over the interatomic data containers and collect the data. 817 data = [] 818 for interatom in interatomic_loop(): 819 # Skip deselected containers. 820 if not interatom.select: 821 continue 822 823 # Skip containers with no RDCs. 824 if not bc and (not hasattr(interatom, 'rdc') or align_id not in interatom.rdc.keys()): 825 continue 826 elif bc and (not hasattr(interatom, 'rdc_bc') or align_id not in interatom.rdc_bc.keys()): 827 continue 828 829 # Append the spin data. 830 data.append([]) 831 data[-1].append(interatom.spin_id1) 832 data[-1].append(interatom.spin_id2) 833 834 # The value. 835 if bc: 836 data[-1].append(repr(convert(interatom.rdc_bc[align_id], align_id))) 837 else: 838 data[-1].append(repr(convert(interatom.rdc[align_id], align_id))) 839 840 # The error. 841 if hasattr(interatom, 'rdc_err') and align_id in interatom.rdc_err.keys(): 842 data[-1].append(repr(convert(interatom.rdc_err[align_id], align_id))) 843 else: 844 data[-1].append(repr(None)) 845 846 # Write out. 847 write_data(out=file, headings=["Spin_ID1", "Spin_ID2", "RDCs", "RDC_error"], data=data)

848

Source Code for Module generic_fns.rdc