Package specific_analyses :: Package consistency_tests
[hide private]
[frames] | no frames]

Source Code for Package specific_analyses.consistency_tests

  1  ############################################################################### 
  2  #                                                                             # 
  3  # Copyright (C) 2004-2013 Edward d'Auvergne                                   # 
  4  # Copyright (C) 2007-2009 Sebastien Morin                                     # 
  5  #                                                                             # 
  6  # This file is part of the program relax (http://www.nmr-relax.com).          # 
  7  #                                                                             # 
  8  # This program is free software: you can redistribute it and/or modify        # 
  9  # it under the terms of the GNU General Public License as published by        # 
 10  # the Free Software Foundation, either version 3 of the License, or           # 
 11  # (at your option) any later version.                                         # 
 12  #                                                                             # 
 13  # This program is distributed in the hope that it will be useful,             # 
 14  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
 15  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
 16  # GNU General Public License for more details.                                # 
 17  #                                                                             # 
 18  # You should have received a copy of the GNU General Public License           # 
 19  # along with this program.  If not, see <http://www.gnu.org/licenses/>.       # 
 20  #                                                                             # 
 21  ############################################################################### 
 22   
 23  # Package docstring. 
 24  """The consistency testing analysis of multi-field relaxation data.""" 
 25   
 26  # The available modules. 
 27  __all__ = [ 
 28  ] 
 29   
 30  # Python module imports. 
 31  from re import search 
 32  from warnings import warn 
 33   
 34  # relax module imports. 
 35  from lib.errors import RelaxError, RelaxFuncSetupError, RelaxNoSequenceError, RelaxNoValueError, RelaxSpinTypeError 
 36  from lib.float import isInf 
 37  from lib.physical_constants import N15_CSA, NH_BOND_LENGTH, h_bar, mu0, return_gyromagnetic_ratio 
 38  from lib.warnings import RelaxDeselectWarning 
 39  from pipe_control import pipes 
 40  from pipe_control.interatomic import return_interatom_list 
 41  from pipe_control.mol_res_spin import exists_mol_res_spin_data, return_spin, spin_loop 
 42  import specific_analyses 
 43  from specific_analyses.api_base import API_base 
 44  from specific_analyses.api_common import API_common 
 45  from target_functions.consistency_tests import Consistency 
 46  from user_functions.data import Uf_tables; uf_tables = Uf_tables() 
 47  from user_functions.objects import Desc_container 
 48   
 49   
50 -class Consistency_tests(API_base, API_common):
51 """Class containing functions specific to consistency testing.""" 52
53 - def __init__(self):
54 """Initialise the class by placing API_common methods into the API.""" 55 56 # Execute the base class __init__ method. 57 super(Consistency_tests, self).__init__() 58 59 # Place methods into the API. 60 self.base_data_loop = self._base_data_loop_spin 61 self.create_mc_data = self._create_mc_relax_data 62 self.model_loop = self._model_loop_spin 63 self.return_conversion_factor = self._return_no_conversion_factor 64 self.return_error = self._return_error_relax_data 65 self.return_value = self._return_value_general 66 self.set_param_values = self._set_param_values_spin 67 self.set_selected_sim = self._set_selected_sim_spin 68 self.sim_pack_data = self._sim_pack_relax_data 69 70 # Set up the spin parameters. 71 self.PARAMS.add('j0', scope='spin', desc='Spectral density value at 0 MHz (from Farrow et al. (1995) JBNMR, 6: 153-162)', py_type=float, grace_string='\\qJ(0)\\Q', err=True, sim=True) 72 self.PARAMS.add('f_eta', scope='spin', desc='Eta-test (from Fushman et al. (1998) JACS, 120: 10947-10952)', py_type=float, grace_string='\\qF\\s\\xh\\Q', err=True, sim=True) 73 self.PARAMS.add('f_r2', scope='spin', desc='R2-test (from Fushman et al. (1998) JACS, 120: 10947-10952)', py_type=float, grace_string='\\qF\\sR2\\Q', err=True, sim=True) 74 self.PARAMS.add('csa', scope='spin', default=N15_CSA, units='ppm', desc='CSA value', py_type=float, grace_string='\\qCSA\\Q') 75 self.PARAMS.add('orientation', scope='spin', default=15.7, units='degrees', desc="Angle between the 15N-1H vector and the principal axis of the 15N chemical shift tensor", py_type=float, grace_string='\\q\\xq\\Q') 76 self.PARAMS.add('tc', scope='spin', default=13 * 1e-9, units='ns', desc="Correlation time", py_type=float, grace_string='\\q\\xt\\f{}c\\Q')
77 78
79 - def _set_frq(self, frq=None):
80 """Function for selecting which relaxation data to use in the consistency tests.""" 81 82 # Test if the current pipe exists. 83 pipes.test() 84 85 # Test if the pipe type is set to 'ct'. 86 function_type = cdp.pipe_type 87 if function_type != 'ct': 88 raise RelaxFuncSetupError(specific_analyses.setup.get_string(function_type)) 89 90 # Test if the frequency has been set. 91 if hasattr(cdp, 'ct_frq'): 92 raise RelaxError("The frequency for the run has already been set.") 93 94 # Create the data structure if it doesn't exist. 95 if not hasattr(cdp, 'ct_frq'): 96 cdp.ct_frq = {} 97 98 # Set the frequency. 99 cdp.ct_frq = frq
100 101
102 - def calculate(self, spin_id=None, verbosity=1, sim_index=None):
103 """Calculation of the consistency functions. 104 105 @keyword spin_id: The spin identification string. 106 @type spin_id: None or str 107 @keyword verbosity: The amount of information to print. The higher the value, the greater the verbosity. 108 @type verbosity: int 109 @keyword sim_index: The optional MC simulation index. 110 @type sim_index: None or int 111 """ 112 113 # Test if the frequency has been set. 114 if not hasattr(cdp, 'ct_frq') or not isinstance(cdp.ct_frq, float): 115 raise RelaxError("The frequency has not been set up.") 116 117 # Test if the sequence data is loaded. 118 if not exists_mol_res_spin_data(): 119 raise RelaxNoSequenceError 120 121 # Test if the spin data has been set. 122 for spin, id in spin_loop(spin_id, return_id=True): 123 # Skip deselected spins. 124 if not spin.select: 125 continue 126 127 # Test if the nuclear isotope type has been set. 128 if not hasattr(spin, 'isotope'): 129 raise RelaxSpinTypeError 130 131 # Test if the CSA value has been set. 132 if not hasattr(spin, 'csa') or spin.csa == None: 133 raise RelaxNoValueError("CSA") 134 135 # Test if the angle Theta has been set. 136 if not hasattr(spin, 'orientation') or spin.orientation == None: 137 raise RelaxNoValueError("angle Theta") 138 139 # Test if the correlation time has been set. 140 if not hasattr(spin, 'tc') or spin.tc == None: 141 raise RelaxNoValueError("correlation time") 142 143 # Test the interatomic data. 144 interatoms = return_interatom_list(id) 145 for interatom in interatoms: 146 # No relaxation mechanism. 147 if not interatom.dipole_pair: 148 continue 149 150 # The interacting spin. 151 if id != interatom.spin_id1: 152 spin_id2 = interatom.spin_id1 153 else: 154 spin_id2 = interatom.spin_id2 155 spin2 = return_spin(spin_id2) 156 157 # Test if the nuclear isotope type has been set. 158 if not hasattr(spin2, 'isotope'): 159 raise RelaxSpinTypeError 160 161 # Test if the interatomic distance has been set. 162 if not hasattr(interatom, 'r') or interatom.r == None: 163 raise RelaxNoValueError("interatomic distance", spin_id=spin_id, spin_id2=spin_id2) 164 165 # Frequency index. 166 if cdp.ct_frq not in cdp.spectrometer_frq.values(): 167 raise RelaxError("No relaxation data corresponding to the frequency %s has been loaded." % cdp.ct_frq) 168 169 # Consistency testing. 170 for spin, id in spin_loop(spin_id, return_id=True): 171 # Skip deselected spins. 172 if not spin.select: 173 continue 174 175 # Set the r1, r2, and NOE to None. 176 r1 = None 177 r2 = None 178 noe = None 179 180 # Get the R1, R2, and NOE values corresponding to the set frequency. 181 for ri_id in cdp.ri_ids: 182 # The frequency does not match. 183 if cdp.spectrometer_frq[ri_id] != cdp.ct_frq: 184 continue 185 186 # R1. 187 if cdp.ri_type[ri_id] == 'R1': 188 if sim_index == None: 189 r1 = spin.ri_data[ri_id] 190 else: 191 r1 = spin.ri_data_sim[ri_id][sim_index] 192 193 # R2. 194 if cdp.ri_type[ri_id] == 'R2': 195 if sim_index == None: 196 r2 = spin.ri_data[ri_id] 197 else: 198 r2 = spin.ri_data_sim[ri_id][sim_index] 199 200 # NOE. 201 if cdp.ri_type[ri_id] == 'NOE': 202 if sim_index == None: 203 noe = spin.ri_data[ri_id] 204 else: 205 noe = spin.ri_data_sim[ri_id][sim_index] 206 207 # Skip the spin if not all of the three value exist. 208 if r1 == None or r2 == None or noe == None: 209 continue 210 211 # Loop over the interatomic data. 212 interatoms = return_interatom_list(id) 213 for i in range(len(interatoms)): 214 # No relaxation mechanism. 215 if not interatoms[i].dipole_pair: 216 continue 217 218 # The surrounding spins. 219 if id != interatoms[i].spin_id1: 220 spin_id2 = interatoms[i].spin_id1 221 else: 222 spin_id2 = interatoms[i].spin_id2 223 spin2 = return_spin(spin_id2) 224 225 # Gyromagnetic ratios. 226 gx = return_gyromagnetic_ratio(spin.isotope) 227 gh = return_gyromagnetic_ratio(spin2.isotope) 228 229 # The interatomic distance. 230 r = interatoms[i].r 231 232 # Initialise the function to calculate. 233 self.ct = Consistency(frq=cdp.ct_frq, gx=gx, gh=gh, mu0=mu0, h_bar=h_bar) 234 235 # Calculate the consistency tests values. 236 j0, f_eta, f_r2 = self.ct.func(orientation=spin.orientation, tc=spin.tc, r=r, csa=spin.csa, r1=r1, r2=r2, noe=noe) 237 238 # Consistency tests values. 239 if sim_index == None: 240 spin.j0 = j0 241 spin.f_eta = f_eta 242 spin.f_r2 = f_r2 243 244 # Monte Carlo simulated consistency tests values. 245 else: 246 # Initialise the simulation data structures. 247 self.data_init(spin, sim=1) 248 if spin.j0_sim == None: 249 spin.j0_sim = [] 250 spin.f_eta_sim = [] 251 spin.f_r2_sim = [] 252 253 # Consistency tests values. 254 spin.j0_sim.append(j0) 255 spin.f_eta_sim.append(f_eta) 256 spin.f_r2_sim.append(f_r2)
257 258
259 - def data_init(self, data_cont, sim=False):
260 """Initialise the data structures. 261 262 @param data_cont: The data container. 263 @type data_cont: instance 264 @keyword sim: The Monte Carlo simulation flag, which if true will initialise the simulation data structure. 265 @type sim: bool 266 """ 267 268 # Get the data names. 269 data_names = self.data_names() 270 271 # Loop over the data structure names. 272 for name in data_names: 273 # Simulation data structures. 274 if sim: 275 # Add '_sim' to the names. 276 name = name + '_sim' 277 278 # If the name is not in 'data_cont', add it. 279 if not hasattr(data_cont, name): 280 # Set the attribute. 281 setattr(data_cont, name, None)
282 283 284 default_value_doc = Desc_container("Consistency testing default values") 285 default_value_doc.add_paragraph("These default values are found in the file 'physical_constants.py'.") 286 _table = uf_tables.add_table(label="table: consistency testing default values", caption="Consistency testing default values.") 287 _table.add_headings(["Data type", "Object name", "Value"]) 288 _table.add_row(["Bond length", "'r'", "1.02 * 1e-10"]) 289 _table.add_row(["CSA", "'csa'", "-172 * 1e-6"]) 290 _table.add_row(["Heteronucleus type", "'heteronuc_type'", "'15N'"]) 291 _table.add_row(["Angle theta", "'proton_type'", "'1H'"]) 292 _table.add_row(["Proton type", "'orientation'", "15.7"]) 293 _table.add_row(["Correlation time", "'tc'", "13 * 1e-9"]) 294 default_value_doc.add_table(_table.label) 295 296
297 - def overfit_deselect(self, data_check=True, verbose=True):
298 """Deselect spins which have insufficient data to support calculation. 299 300 @keyword data_check: A flag to signal if the presence of base data is to be checked for. 301 @type data_check: bool 302 @keyword verbose: A flag which if True will allow printouts. 303 @type verbose: bool 304 """ 305 306 # Print out. 307 if verbose: 308 print("\nOver-fit spin deselection:") 309 310 # Test the sequence data exists. 311 if not exists_mol_res_spin_data(): 312 raise RelaxNoSequenceError 313 314 # Loop over spin data. 315 deselect_flag = False 316 spin_count = 0 317 for spin, spin_id in spin_loop(return_id=True): 318 # Skip deselected spins. 319 if not spin.select: 320 continue 321 322 # The interatomic data. 323 interatoms = return_interatom_list(spin_id) 324 325 # Loop over the interatomic data. 326 dipole_relax = False 327 for i in range(len(interatoms)): 328 # No dipolar relaxation mechanism. 329 if not interatoms[i].dipole_pair: 330 continue 331 332 # The surrounding spins. 333 if spin_id != interatoms[i].spin_id1: 334 spin_id2 = interatoms[i].spin_id1 335 else: 336 spin_id2 = interatoms[i].spin_id2 337 spin2 = return_spin(spin_id2) 338 339 # Dipolar relaxation flag. 340 dipole_relax = True 341 342 # No relaxation mechanism. 343 if not dipole_relax or not hasattr(spin, 'csa') or spin.csa == None: 344 warn(RelaxDeselectWarning(spin_id, 'an absence of relaxation mechanisms')) 345 spin.select = False 346 deselect_flag = True 347 continue 348 349 # Data checks. 350 if data_check: 351 # The number of relaxation data points (and for infinite data). 352 data_points = 0 353 inf_data = False 354 if hasattr(cdp, 'ri_ids') and hasattr(spin, 'ri_data'): 355 for id in cdp.ri_ids: 356 if id in spin.ri_data and spin.ri_data[id] != None: 357 data_points += 1 358 359 # Infinite data! 360 if isInf(spin.ri_data[id]): 361 inf_data = True 362 363 # Infinite data. 364 if inf_data: 365 warn(RelaxDeselectWarning(spin_id, 'infinite relaxation data')) 366 spin.select = False 367 deselect_flag = True 368 continue 369 370 # Relaxation data must exist! 371 if not hasattr(spin, 'ri_data'): 372 warn(RelaxDeselectWarning(spin_id, 'missing relaxation data')) 373 spin.select = False 374 deselect_flag = True 375 continue 376 377 # Require 3 or more relaxation data points. 378 if data_points < 3: 379 warn(RelaxDeselectWarning(spin_id, 'insufficient relaxation data, 3 or more data points are required')) 380 spin.select = False 381 deselect_flag = True 382 continue 383 384 # Increment the spin number. 385 spin_count += 1 386 387 # No spins selected, so fail hard to prevent the user from going any further. 388 if spin_count == 0: 389 warn(RelaxWarning("No spins are selected therefore the optimisation or calculation cannot proceed.")) 390 391 # Final printout. 392 if verbose and not deselect_flag: 393 print("No spins have been deselected.")
394 395 396 return_data_name_doc = Desc_container("Consistency testing data type string matching patterns") 397 _table = uf_tables.add_table(label="table: Consistency testing data types", caption="Consistency testing data type string matching patterns.") 398 _table.add_headings(["Data type", "Object name"]) 399 _table.add_row(["J(0)", "'j0'"]) 400 _table.add_row(["F_eta", "'f_eta'"]) 401 _table.add_row(["F_R2", "'f_r2'"]) 402 _table.add_row(["Bond length", "'r'"]) 403 _table.add_row(["CSA", "'csa'"]) 404 _table.add_row(["Heteronucleus type", "'heteronuc_type'"]) 405 _table.add_row(["Proton type", "'proton_type'"]) 406 _table.add_row(["Angle theta", "'orientation'"]) 407 _table.add_row(["Correlation time", "'tc'"]) 408 return_data_name_doc.add_table(_table.label) 409 410 411 set_doc = Desc_container("Consistency testing set details") 412 set_doc.add_paragraph("In consistency testing, only four values can be set, the bond length, CSA, angle Theta ('orientation') and correlation time values. These must be set prior to the calculation of consistency functions.") 413 414
415 - def set_error(self, model_info, index, error):
416 """Set the parameter errors. 417 418 @param model_info: The spin container originating from model_loop(). 419 @type model_info: SpinContainer instance 420 @param index: The index of the parameter to set the errors for. 421 @type index: int 422 @param error: The error value. 423 @type error: float 424 """ 425 426 # Alias. 427 spin = model_info 428 429 # Return J(0) sim data. 430 if index == 0: 431 spin.j0_err = error 432 433 # Return F_eta sim data. 434 if index == 1: 435 spin.f_eta_err = error 436 437 # Return F_R2 sim data. 438 if index == 2: 439 spin.f_r2_err = error
440 441
442 - def sim_return_param(self, model_info, index):
443 """Return the array of simulation parameter values. 444 445 @param model_info: The spin container originating from model_loop(). 446 @type model_info: SpinContainer instance 447 @param index: The index of the parameter to return the array of values for. 448 @type index: int 449 @return: The array of simulation parameter values. 450 @rtype: list of float 451 """ 452 453 # Alias. 454 spin = model_info 455 456 # Skip deselected spins. 457 if not spin.select: 458 return 459 460 # Return J(0) sim data. 461 if index == 0: 462 return spin.j0_sim 463 464 # Return F_eta sim data. 465 if index == 1: 466 return spin.f_eta_sim 467 468 # Return F_R2 sim data. 469 if index == 2: 470 return spin.f_r2_sim
471 472
473 - def sim_return_selected(self, model_info):
474 """Return the array of selected simulation flags. 475 476 @param model_info: The spin container originating from model_loop(). 477 @type model_info: SpinContainer instance 478 @return: The array of selected simulation flags. 479 @rtype: list of int 480 """ 481 482 # Alias. 483 spin = model_info 484 485 # Multiple spins. 486 return spin.select_sim
487