Package specific_analyses :: Package model_free :: Module parameters
[hide private]
[frames] | no frames]

Source Code for Module specific_analyses.model_free.parameters

   1  ############################################################################### 
   2  #                                                                             # 
   3  # Copyright (C) 2003-2014 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  """The model-free analysis parameter functions.""" 
  24   
  25  # Python module imports. 
  26  from math import pi 
  27  from numpy import array, float64, identity, int8, zeros 
  28  from re import match, search 
  29   
  30  # relax module imports. 
  31  from lib.errors import RelaxError 
  32  from pipe_control import diffusion_tensor 
  33  from pipe_control.mol_res_spin import spin_loop 
  34  from specific_analyses.model_free.model import determine_model_type 
  35   
  36   
37 -def are_mf_params_set(spin):
38 """Test if the model-free parameter values are set. 39 40 @param spin: The spin container object. 41 @type spin: SpinContainer instance 42 @return: The name of the first parameter in the parameter list in which the 43 corresponding parameter value is None. If all parameters are set, then None 44 is returned. 45 @rtype: str or None 46 """ 47 48 # Deselected residue. 49 if spin.select == 0: 50 return 51 52 # Loop over the model-free parameters. 53 for j in range(len(spin.params)): 54 # Local tm. 55 if spin.params[j] == 'local_tm' and spin.local_tm == None: 56 return spin.params[j] 57 58 # S2. 59 elif spin.params[j] == 's2' and spin.s2 == None: 60 return spin.params[j] 61 62 # S2f. 63 elif spin.params[j] == 's2f' and spin.s2f == None: 64 return spin.params[j] 65 66 # S2s. 67 elif spin.params[j] == 's2s' and spin.s2s == None: 68 return spin.params[j] 69 70 # te. 71 elif spin.params[j] == 'te' and spin.te == None: 72 return spin.params[j] 73 74 # tf. 75 elif spin.params[j] == 'tf' and spin.tf == None: 76 return spin.params[j] 77 78 # ts. 79 elif spin.params[j] == 'ts' and spin.ts == None: 80 return spin.params[j] 81 82 # Rex. 83 elif spin.params[j] == 'rex' and spin.rex == None: 84 return spin.params[j] 85 86 # r. 87 elif spin.params[j] == 'r' and spin.r == None: 88 return spin.params[j] 89 90 # CSA. 91 elif spin.params[j] == 'csa' and spin.csa == None: 92 return spin.params[j]
93 94
95 -def assemble_param_names(model_type, spin_id=None):
96 """Function for assembling a list of all the model parameter names. 97 98 @param model_type: The model-free model type. This must be one of 'mf', 'local_tm', 99 'diff', or 'all'. 100 @type model_type: str 101 @param spin_id: The spin identification string. 102 @type spin_id: str 103 @return: A list containing all the parameters of the model-free model. 104 @rtype: list of str 105 """ 106 107 # Initialise. 108 param_names = [] 109 110 # Diffusion tensor parameters. 111 if model_type == 'diff' or model_type == 'all': 112 # Spherical diffusion. 113 if cdp.diff_tensor.type == 'sphere': 114 param_names.append('tm') 115 116 # Spheroidal diffusion. 117 elif cdp.diff_tensor.type == 'spheroid': 118 param_names.append('tm') 119 param_names.append('Da') 120 param_names.append('theta') 121 param_names.append('phi') 122 123 # Ellipsoidal diffusion. 124 elif cdp.diff_tensor.type == 'ellipsoid': 125 param_names.append('tm') 126 param_names.append('Da') 127 param_names.append('Dr') 128 param_names.append('alpha') 129 param_names.append('beta') 130 param_names.append('gamma') 131 132 # Model-free parameters (spin specific parameters). 133 if model_type != 'diff': 134 # Loop over the spins. 135 for spin in spin_loop(spin_id): 136 # Skip deselected spins. 137 if not spin.select: 138 continue 139 140 # Add the spin specific model-free parameters. 141 param_names = param_names + spin.params 142 143 # Return the parameter names. 144 return param_names
145 146
147 -def assemble_param_vector(spin=None, spin_id=None, sim_index=None, model_type=None):
148 """Assemble the model-free parameter vector (as numpy array). 149 150 If the spin argument is supplied, then the spin_id argument will be ignored. 151 152 @keyword spin: The spin data container. 153 @type spin: SpinContainer instance 154 @keyword spin_id: The spin identification string. 155 @type spin_id: str 156 @keyword sim_index: The optional MC simulation index. 157 @type sim_index: int 158 @keyword model_type: The optional model type, one of 'all', 'diff', 'mf', or 'local_tm'. 159 @type model_type: str or None 160 @return: An array of the parameter values of the model-free model. 161 @rtype: numpy array 162 """ 163 164 # Initialise. 165 param_vector = [] 166 167 # Determine the model type. 168 if not model_type: 169 model_type = determine_model_type() 170 171 # Diffusion tensor parameters. 172 if model_type == 'diff' or model_type == 'all': 173 # Normal parameters. 174 if sim_index == None: 175 # Spherical diffusion. 176 if cdp.diff_tensor.type == 'sphere': 177 param_vector.append(cdp.diff_tensor.tm) 178 179 # Spheroidal diffusion. 180 elif cdp.diff_tensor.type == 'spheroid': 181 param_vector.append(cdp.diff_tensor.tm) 182 param_vector.append(cdp.diff_tensor.Da) 183 param_vector.append(cdp.diff_tensor.theta) 184 param_vector.append(cdp.diff_tensor.phi) 185 186 # Ellipsoidal diffusion. 187 elif cdp.diff_tensor.type == 'ellipsoid': 188 param_vector.append(cdp.diff_tensor.tm) 189 param_vector.append(cdp.diff_tensor.Da) 190 param_vector.append(cdp.diff_tensor.Dr) 191 param_vector.append(cdp.diff_tensor.alpha) 192 param_vector.append(cdp.diff_tensor.beta) 193 param_vector.append(cdp.diff_tensor.gamma) 194 195 # Monte Carlo diffusion tensor parameters. 196 else: 197 # Spherical diffusion. 198 if cdp.diff_tensor.type == 'sphere': 199 param_vector.append(cdp.diff_tensor.tm_sim[sim_index]) 200 201 # Spheroidal diffusion. 202 elif cdp.diff_tensor.type == 'spheroid': 203 param_vector.append(cdp.diff_tensor.tm_sim[sim_index]) 204 param_vector.append(cdp.diff_tensor.Da_sim[sim_index]) 205 param_vector.append(cdp.diff_tensor.theta_sim[sim_index]) 206 param_vector.append(cdp.diff_tensor.phi_sim[sim_index]) 207 208 # Ellipsoidal diffusion. 209 elif cdp.diff_tensor.type == 'ellipsoid': 210 param_vector.append(cdp.diff_tensor.tm_sim[sim_index]) 211 param_vector.append(cdp.diff_tensor.Da_sim[sim_index]) 212 param_vector.append(cdp.diff_tensor.Dr_sim[sim_index]) 213 param_vector.append(cdp.diff_tensor.alpha_sim[sim_index]) 214 param_vector.append(cdp.diff_tensor.beta_sim[sim_index]) 215 param_vector.append(cdp.diff_tensor.gamma_sim[sim_index]) 216 217 # Model-free parameters (spin specific parameters). 218 if model_type != 'diff': 219 # The loop. 220 if spin: 221 loop = [spin] 222 else: 223 loop = spin_loop(spin_id) 224 225 # Loop over the spins. 226 for spin in loop: 227 # Skip deselected spins. 228 if not spin.select: 229 continue 230 231 # Skip spins with no parameters. 232 if not hasattr(spin, 'params'): 233 continue 234 235 # Loop over the model-free parameters. 236 for i in range(len(spin.params)): 237 # local tm. 238 if spin.params[i] == 'local_tm': 239 if sim_index == None: 240 param_vector.append(spin.local_tm) 241 else: 242 param_vector.append(spin.local_tm_sim[sim_index]) 243 244 # S2. 245 elif spin.params[i] == 's2': 246 if sim_index == None: 247 param_vector.append(spin.s2) 248 else: 249 param_vector.append(spin.s2_sim[sim_index]) 250 251 # S2f. 252 elif spin.params[i] == 's2f': 253 if sim_index == None: 254 param_vector.append(spin.s2f) 255 else: 256 param_vector.append(spin.s2f_sim[sim_index]) 257 258 # S2s. 259 elif spin.params[i] == 's2s': 260 if sim_index == None: 261 param_vector.append(spin.s2s) 262 else: 263 param_vector.append(spin.s2s_sim[sim_index]) 264 265 # te. 266 elif spin.params[i] == 'te': 267 if sim_index == None: 268 param_vector.append(spin.te) 269 else: 270 param_vector.append(spin.te_sim[sim_index]) 271 272 # tf. 273 elif spin.params[i] == 'tf': 274 if sim_index == None: 275 param_vector.append(spin.tf) 276 else: 277 param_vector.append(spin.tf_sim[sim_index]) 278 279 # ts. 280 elif spin.params[i] == 'ts': 281 if sim_index == None: 282 param_vector.append(spin.ts) 283 else: 284 param_vector.append(spin.ts_sim[sim_index]) 285 286 # Rex. 287 elif spin.params[i] == 'rex': 288 if sim_index == None: 289 param_vector.append(spin.rex) 290 else: 291 param_vector.append(spin.rex_sim[sim_index]) 292 293 # r. 294 elif spin.params[i] == 'r': 295 if sim_index == None: 296 param_vector.append(spin.r) 297 else: 298 param_vector.append(spin.r_sim[sim_index]) 299 300 # CSA. 301 elif spin.params[i] == 'csa': 302 if sim_index == None: 303 param_vector.append(spin.csa) 304 else: 305 param_vector.append(spin.csa_sim[sim_index]) 306 307 # Unknown parameter. 308 else: 309 raise RelaxError("Unknown parameter.") 310 311 # Replace all instances of None with 0.0 to allow the list to be converted to a numpy array. 312 for i in range(len(param_vector)): 313 if param_vector[i] == None: 314 param_vector[i] = 0.0 315 316 # Return a numpy array. 317 return array(param_vector, float64)
318 319
320 -def assemble_scaling_matrix(num_params, model_type=None, spin=None, spin_id=None, scaling=True):
321 """Create and return the scaling matrix. 322 323 If the spin argument is supplied, then the spin_id argument will be ignored. 324 325 @param num_params: The number of parameters in the model. 326 @type num_params: int 327 @keyword model_type: The model type, one of 'all', 'diff', 'mf', or 'local_tm'. 328 @type model_type: str 329 @keyword spin: The spin data container. 330 @type spin: SpinContainer instance 331 @keyword spin_id: The spin identification string. 332 @type spin_id: str 333 @return: The diagonal and square scaling matrix. 334 @rtype: numpy diagonal matrix 335 """ 336 337 # Initialise. 338 if num_params == 0: 339 scaling_matrix = zeros((0, 0), float64) 340 else: 341 scaling_matrix = identity(num_params, float64) 342 i = 0 343 344 # No diagonal scaling, so return the identity matrix. 345 if not scaling: 346 return scaling_matrix 347 348 # tm, te, tf, and ts (must all be the same for diagonal scaling!). 349 ti_scaling = 1e-12 350 351 # Diffusion tensor parameters. 352 if model_type == 'diff' or model_type == 'all': 353 # Spherical diffusion. 354 if cdp.diff_tensor.type == 'sphere': 355 # tm. 356 scaling_matrix[i, i] = ti_scaling 357 358 # Increment i. 359 i = i + 1 360 361 # Spheroidal diffusion. 362 elif cdp.diff_tensor.type == 'spheroid': 363 # tm, Da, theta, phi 364 scaling_matrix[i, i] = ti_scaling 365 scaling_matrix[i+1, i+1] = 1e7 366 scaling_matrix[i+2, i+2] = 1.0 367 scaling_matrix[i+3, i+3] = 1.0 368 369 # Increment i. 370 i = i + 4 371 372 # Ellipsoidal diffusion. 373 elif cdp.diff_tensor.type == 'ellipsoid': 374 # tm, Da, Dr, alpha, beta, gamma. 375 scaling_matrix[i, i] = ti_scaling 376 scaling_matrix[i+1, i+1] = 1e7 377 scaling_matrix[i+2, i+2] = 1.0 378 scaling_matrix[i+3, i+3] = 1.0 379 scaling_matrix[i+4, i+4] = 1.0 380 scaling_matrix[i+5, i+5] = 1.0 381 382 # Increment i. 383 i = i + 6 384 385 # Model-free parameters. 386 if model_type != 'diff': 387 # The loop. 388 if spin: 389 loop = [spin] 390 else: 391 loop = spin_loop(spin_id) 392 393 # Loop over the spins. 394 for spin in loop: 395 # Skip deselected spins. 396 if not spin.select: 397 continue 398 399 # Loop over the model-free parameters. 400 for k in range(len(spin.params)): 401 # Local tm, te, tf, and ts (must all be the same for diagonal scaling!). 402 if spin.params[k] == 'local_tm' or search('^t', spin.params[k]): 403 scaling_matrix[i, i] = ti_scaling 404 405 # Rex. 406 elif spin.params[k] == 'rex': 407 scaling_matrix[i, i] = 1.0 / (2.0 * pi * cdp.spectrometer_frq[cdp.ri_ids[0]]) ** 2 408 409 # Interatomic distances. 410 elif spin.params[k] == 'r': 411 scaling_matrix[i, i] = 1e-10 412 413 # CSA. 414 elif spin.params[k] == 'csa': 415 scaling_matrix[i, i] = 1e-4 416 417 # Increment i. 418 i = i + 1 419 420 # Return the scaling matrix. 421 return scaling_matrix
422 423
424 -def conv_factor_rex():
425 """Calculate and return the Rex conversion factor. 426 427 @return: The Rex conversion factor. 428 @rtype: float 429 """ 430 431 # No frequency info. 432 if not hasattr(cdp, 'spectrometer_frq'): 433 raise RelaxError("No spectrometer frequency information is present in the current data pipe.") 434 435 # The 1st spectrometer frequency. 436 if hasattr(cdp, 'ri_ids'): 437 frq = cdp.spectrometer_frq[cdp.ri_ids[0]] 438 439 # Take the highest frequency, if all else fails. 440 else: 441 frqs = sorted(cdp.spectrometer_frq.values()) 442 frq = frqs[-1] 443 444 # The factor. 445 return 1.0 / (2.0 * pi * frq)**2
446 447
448 -def disassemble_param_vector(model_type, param_vector=None, spin=None, spin_id=None, sim_index=None):
449 """Disassemble the model-free parameter vector. 450 451 @param model_type: The model-free model type. This must be one of 'mf', 'local_tm', 452 'diff', or 'all'. 453 @type model_type: str 454 @keyword param_vector: The model-free parameter vector. 455 @type param_vector: numpy array 456 @keyword spin: The spin data container. If this argument is supplied, then the spin_id 457 argument will be ignored. 458 @type spin: SpinContainer instance 459 @keyword spin_id: The spin identification string. 460 @type spin_id: str 461 @keyword sim_index: The optional MC simulation index. 462 @type sim_index: int 463 """ 464 465 # Initialise. 466 param_index = 0 467 468 # Diffusion tensor parameters of the Monte Carlo simulations. 469 if sim_index != None and (model_type == 'diff' or model_type == 'all'): 470 # Spherical diffusion. 471 if cdp.diff_tensor.type == 'sphere': 472 # Sim values. 473 cdp.diff_tensor.set(param='tm', value=param_vector[0], category='sim', sim_index=sim_index) 474 475 # Parameter index. 476 param_index = param_index + 1 477 478 # Spheroidal diffusion. 479 elif cdp.diff_tensor.type == 'spheroid': 480 # Sim values. 481 cdp.diff_tensor.set(param='tm', value=param_vector[0], category='sim', sim_index=sim_index) 482 cdp.diff_tensor.set(param='Da', value=param_vector[1], category='sim', sim_index=sim_index) 483 cdp.diff_tensor.set(param='theta', value=param_vector[2], category='sim', sim_index=sim_index) 484 cdp.diff_tensor.set(param='phi', value=param_vector[3], category='sim', sim_index=sim_index) 485 diffusion_tensor.fold_angles(sim_index=sim_index) 486 487 # Parameter index. 488 param_index = param_index + 4 489 490 # Ellipsoidal diffusion. 491 elif cdp.diff_tensor.type == 'ellipsoid': 492 # Sim values. 493 cdp.diff_tensor.set(param='tm', value=param_vector[0], category='sim', sim_index=sim_index) 494 cdp.diff_tensor.set(param='Da', value=param_vector[1], category='sim', sim_index=sim_index) 495 cdp.diff_tensor.set(param='Dr', value=param_vector[2], category='sim', sim_index=sim_index) 496 cdp.diff_tensor.set(param='alpha', value=param_vector[3], category='sim', sim_index=sim_index) 497 cdp.diff_tensor.set(param='beta', value=param_vector[4], category='sim', sim_index=sim_index) 498 cdp.diff_tensor.set(param='gamma', value=param_vector[5], category='sim', sim_index=sim_index) 499 diffusion_tensor.fold_angles(sim_index=sim_index) 500 501 # Parameter index. 502 param_index = param_index + 6 503 504 # Diffusion tensor parameters. 505 elif model_type == 'diff' or model_type == 'all': 506 # Spherical diffusion. 507 if cdp.diff_tensor.type == 'sphere': 508 # Values. 509 cdp.diff_tensor.set(param='tm', value=param_vector[0]) 510 511 # Parameter index. 512 param_index = param_index + 1 513 514 # Spheroidal diffusion. 515 elif cdp.diff_tensor.type == 'spheroid': 516 # Values. 517 cdp.diff_tensor.set(param='tm', value=param_vector[0]) 518 cdp.diff_tensor.set(param='Da', value=param_vector[1]) 519 cdp.diff_tensor.set(param='theta', value=param_vector[2]) 520 cdp.diff_tensor.set(param='phi', value=param_vector[3]) 521 diffusion_tensor.fold_angles() 522 523 # Parameter index. 524 param_index = param_index + 4 525 526 # Ellipsoidal diffusion. 527 elif cdp.diff_tensor.type == 'ellipsoid': 528 # Values. 529 cdp.diff_tensor.set(param='tm', value=param_vector[0]) 530 cdp.diff_tensor.set(param='Da', value=param_vector[1]) 531 cdp.diff_tensor.set(param='Dr', value=param_vector[2]) 532 cdp.diff_tensor.set(param='alpha', value=param_vector[3]) 533 cdp.diff_tensor.set(param='beta', value=param_vector[4]) 534 cdp.diff_tensor.set(param='gamma', value=param_vector[5]) 535 diffusion_tensor.fold_angles() 536 537 # Parameter index. 538 param_index = param_index + 6 539 540 # Model-free parameters. 541 if model_type != 'diff': 542 # The loop. 543 if spin: 544 loop = [spin] 545 else: 546 loop = spin_loop(spin_id) 547 548 # Loop over the spins. 549 for spin in loop: 550 # Skip deselected spins. 551 if not spin.select: 552 continue 553 554 # Loop over the model-free parameters. 555 for j in range(len(spin.params)): 556 # Local tm. 557 if spin.params[j] == 'local_tm': 558 if sim_index == None: 559 spin.local_tm = param_vector[param_index] 560 else: 561 spin.local_tm_sim[sim_index] = param_vector[param_index] 562 563 # S2. 564 elif spin.params[j] == 's2': 565 if sim_index == None: 566 spin.s2 = param_vector[param_index] 567 else: 568 spin.s2_sim[sim_index] = param_vector[param_index] 569 570 # S2f. 571 elif spin.params[j] == 's2f': 572 if sim_index == None: 573 spin.s2f = param_vector[param_index] 574 else: 575 spin.s2f_sim[sim_index] = param_vector[param_index] 576 577 # S2s. 578 elif spin.params[j] == 's2s': 579 if sim_index == None: 580 spin.s2s = param_vector[param_index] 581 else: 582 spin.s2s_sim[sim_index] = param_vector[param_index] 583 584 # te. 585 elif spin.params[j] == 'te': 586 if sim_index == None: 587 spin.te = param_vector[param_index] 588 else: 589 spin.te_sim[sim_index] = param_vector[param_index] 590 591 # tf. 592 elif spin.params[j] == 'tf': 593 if sim_index == None: 594 spin.tf = param_vector[param_index] 595 else: 596 spin.tf_sim[sim_index] = param_vector[param_index] 597 598 # ts. 599 elif spin.params[j] == 'ts': 600 if sim_index == None: 601 spin.ts = param_vector[param_index] 602 else: 603 spin.ts_sim[sim_index] = param_vector[param_index] 604 605 # Rex. 606 elif spin.params[j] == 'rex': 607 if sim_index == None: 608 spin.rex = param_vector[param_index] 609 else: 610 spin.rex_sim[sim_index] = param_vector[param_index] 611 612 # r. 613 elif spin.params[j] == 'r': 614 if sim_index == None: 615 spin.r = param_vector[param_index] 616 else: 617 spin.r_sim[sim_index] = param_vector[param_index] 618 619 # CSA. 620 elif spin.params[j] == 'csa': 621 if sim_index == None: 622 spin.csa = param_vector[param_index] 623 else: 624 spin.csa_sim[sim_index] = param_vector[param_index] 625 626 # Unknown parameter. 627 else: 628 raise RelaxError("Unknown parameter.") 629 630 # Increment the parameter index. 631 param_index = param_index + 1 632 633 # Calculate all order parameters after unpacking the vector. 634 if model_type != 'diff': 635 # The loop. 636 if spin: 637 loop = [spin] 638 else: 639 loop = spin_loop(spin_id) 640 641 # Loop over the spins. 642 for spin in loop: 643 # Skip deselected residues. 644 if not spin.select: 645 continue 646 647 # Normal values. 648 if sim_index == None: 649 # S2. 650 if 's2' not in spin.params and 's2f' in spin.params and 's2s' in spin.params: 651 spin.s2 = spin.s2f * spin.s2s 652 653 # S2f. 654 if 's2f' not in spin.params and 's2' in spin.params and 's2s' in spin.params: 655 if spin.s2s == 0.0: 656 spin.s2f = 1e99 657 else: 658 spin.s2f = spin.s2 / spin.s2s 659 660 # S2s. 661 if 's2s' not in spin.params and 's2' in spin.params and 's2f' in spin.params: 662 if spin.s2f == 0.0: 663 spin.s2s = 1e99 664 else: 665 spin.s2s = spin.s2 / spin.s2f 666 667 # Simulation values. 668 else: 669 # S2. 670 if 's2' not in spin.params and 's2f' in spin.params and 's2s' in spin.params: 671 spin.s2_sim[sim_index] = spin.s2f_sim[sim_index] * spin.s2s_sim[sim_index] 672 673 # S2f. 674 if 's2f' not in spin.params and 's2' in spin.params and 's2s' in spin.params: 675 if spin.s2s_sim[sim_index] == 0.0: 676 spin.s2f_sim[sim_index] = 1e99 677 else: 678 spin.s2f_sim[sim_index] = spin.s2_sim[sim_index] / spin.s2s_sim[sim_index] 679 680 # S2s. 681 if 's2s' not in spin.params and 's2' in spin.params and 's2f' in spin.params: 682 if spin.s2f_sim[sim_index] == 0.0: 683 spin.s2s_sim[sim_index] = 1e99 684 else: 685 spin.s2s_sim[sim_index] = spin.s2_sim[sim_index] / spin.s2f_sim[sim_index]
686 687
688 -def linear_constraints(num_params, model_type=None, spin=None, spin_id=None, scaling_matrix=None):
689 """Set up the model-free linear constraint matrices A and b. 690 691 Standard notation 692 ================= 693 694 The order parameter constraints are:: 695 696 0 <= S2 <= 1 697 0 <= S2f <= 1 698 0 <= S2s <= 1 699 700 By substituting the formula S2 = S2f.S2s into the above inequalities, the additional two 701 inequalities can be derived:: 702 703 S2 <= S2f 704 S2 <= S2s 705 706 Correlation time constraints are:: 707 708 te >= 0 709 tf >= 0 710 ts >= 0 711 712 tf <= ts 713 714 te, tf, ts <= 2 * tm 715 716 Additional constraints used include:: 717 718 Rex >= 0 719 0.9e-10 <= r <= 2e-10 720 -300e-6 <= CSA <= 0 721 722 723 Rearranged notation 724 =================== 725 726 The above inequality constraints can be rearranged into:: 727 728 S2 >= 0 729 -S2 >= -1 730 S2f >= 0 731 -S2f >= -1 732 S2s >= 0 733 -S2s >= -1 734 S2f - S2 >= 0 735 S2s - S2 >= 0 736 te >= 0 737 tf >= 0 738 ts >= 0 739 ts - tf >= 0 740 Rex >= 0 741 r >= 0.9e-10 742 -r >= -2e-10 743 CSA >= -300e-6 744 -CSA >= 0 745 746 747 Matrix notation 748 =============== 749 750 In the notation A.x >= b, where A is an matrix of coefficients, x is an array of parameter 751 values, and b is a vector of scalars, these inequality constraints are:: 752 753 | 1 0 0 0 0 0 0 0 0 | | 0 | 754 | | | | 755 |-1 0 0 0 0 0 0 0 0 | | -1 | 756 | | | | 757 | 0 1 0 0 0 0 0 0 0 | | 0 | 758 | | | | 759 | 0 -1 0 0 0 0 0 0 0 | | -1 | 760 | | | | 761 | 0 0 1 0 0 0 0 0 0 | | S2 | | 0 | 762 | | | | | | 763 | 0 0 -1 0 0 0 0 0 0 | | S2f | | -1 | 764 | | | | | | 765 |-1 1 0 0 0 0 0 0 0 | | S2s | | 0 | 766 | | | | | | 767 |-1 0 1 0 0 0 0 0 0 | | te | | 0 | 768 | | | | | | 769 | 0 0 0 1 0 0 0 0 0 | . | tf | >= | 0 | 770 | | | | | | 771 | 0 0 0 0 1 0 0 0 0 | | ts | | 0 | 772 | | | | | | 773 | 0 0 0 0 0 1 0 0 0 | | Rex | | 0 | 774 | | | | | | 775 | 0 0 0 0 -1 1 0 0 0 | | r | | 0 | 776 | | | | | | 777 | 0 0 0 0 0 0 1 0 0 | | CSA | | 0 | 778 | | | | 779 | 0 0 0 0 0 0 0 1 0 | | 0.9e-10 | 780 | | | | 781 | 0 0 0 0 0 0 0 -1 0 | | -2e-10 | 782 | | | | 783 | 0 0 0 0 0 0 0 0 1 | | -300e-6 | 784 | | | | 785 | 0 0 0 0 0 0 0 0 -1 | | 0 | 786 787 788 @param num_params: The number of parameters in the model. 789 @type num_params: int 790 @keyword model_type: The model type, one of 'all', 'diff', 'mf', or 'local_tm'. 791 @type model_type: str 792 @keyword spin: The spin data container. If this argument is supplied, then the 793 spin_id argument will be ignored. 794 @type spin: SpinContainer instance 795 @keyword spin_id: The spin identification string. 796 @type spin_id: str 797 @keyword scaling_matrix: The diagonal, square scaling matrix. 798 @type scaling_matrix: numpy diagonal matrix 799 """ 800 801 # Upper limit flag for correlation times. 802 upper_time_limit = 1 803 804 # Initialisation (0..j..m). 805 A = [] 806 b = [] 807 zero_array = zeros(num_params, float64) 808 i = 0 809 j = 0 810 811 # Diffusion tensor parameters. 812 if model_type != 'mf' and diffusion_tensor.diff_data_exists(): 813 # Spherical diffusion. 814 if cdp.diff_tensor.type == 'sphere': 815 # 0 <= tm <= 200 ns. 816 A.append(zero_array * 0.0) 817 A.append(zero_array * 0.0) 818 A[j][i] = 1.0 819 A[j+1][i] = -1.0 820 b.append(0.0 / scaling_matrix[i, i]) 821 b.append(-200.0 * 1e-9 / scaling_matrix[i, i]) 822 i = i + 1 823 j = j + 2 824 825 # Spheroidal diffusion. 826 elif cdp.diff_tensor.type == 'spheroid': 827 # 0 <= tm <= 200 ns. 828 A.append(zero_array * 0.0) 829 A.append(zero_array * 0.0) 830 A[j][i] = 1.0 831 A[j+1][i] = -1.0 832 b.append(0.0 / scaling_matrix[i, i]) 833 b.append(-200.0 * 1e-9 / scaling_matrix[i, i]) 834 i = i + 1 835 j = j + 2 836 837 # Prolate diffusion, Da >= 0. 838 if cdp.diff_tensor.spheroid_type == 'prolate': 839 A.append(zero_array * 0.0) 840 A[j][i] = 1.0 841 b.append(0.0 / scaling_matrix[i, i]) 842 i = i + 1 843 j = j + 1 844 845 # Add two to i for the theta and phi parameters. 846 i = i + 2 847 848 # Oblate diffusion, Da <= 0. 849 elif cdp.diff_tensor.spheroid_type == 'oblate': 850 A.append(zero_array * 0.0) 851 A[j][i] = -1.0 852 b.append(0.0 / scaling_matrix[i, i]) 853 i = i + 1 854 j = j + 1 855 856 # Add two to i for the theta and phi parameters. 857 i = i + 2 858 859 else: 860 # Add three to i for the Da, theta and phi parameters. 861 i = i + 3 862 863 # Ellipsoidal diffusion. 864 elif cdp.diff_tensor.type == 'ellipsoid': 865 # 0 <= tm <= 200 ns. 866 A.append(zero_array * 0.0) 867 A.append(zero_array * 0.0) 868 A[j][i] = 1.0 869 A[j+1][i] = -1.0 870 b.append(0.0 / scaling_matrix[i, i]) 871 b.append(-200.0 * 1e-9 / scaling_matrix[i, i]) 872 i = i + 1 873 j = j + 2 874 875 # Da >= 0. 876 A.append(zero_array * 0.0) 877 A[j][i] = 1.0 878 b.append(0.0 / scaling_matrix[i, i]) 879 i = i + 1 880 j = j + 1 881 882 # 0 <= Dr <= 1. 883 A.append(zero_array * 0.0) 884 A.append(zero_array * 0.0) 885 A[j][i] = 1.0 886 A[j+1][i] = -1.0 887 b.append(0.0 / scaling_matrix[i, i]) 888 b.append(-1.0 / scaling_matrix[i, i]) 889 i = i + 1 890 j = j + 2 891 892 # Add three to i for the alpha, beta, and gamma parameters. 893 i = i + 3 894 895 # Model-free parameters. 896 if model_type != 'diff': 897 # The loop. 898 if spin: 899 loop = [spin] 900 else: 901 loop = spin_loop(spin_id) 902 903 # Loop over the spins. 904 for spin in loop: 905 # Skip deselected spins. 906 if not spin.select: 907 continue 908 909 # Save current value of i. 910 old_i = i 911 912 # Loop over the model-free parameters. 913 for l in range(len(spin.params)): 914 # Local tm. 915 if spin.params[l] == 'local_tm': 916 if upper_time_limit: 917 # 0 <= tm <= 200 ns. 918 A.append(zero_array * 0.0) 919 A.append(zero_array * 0.0) 920 A[j][i] = 1.0 921 A[j+1][i] = -1.0 922 b.append(0.0 / scaling_matrix[i, i]) 923 b.append(-200.0 * 1e-9 / scaling_matrix[i, i]) 924 j = j + 2 925 else: 926 # 0 <= tm. 927 A.append(zero_array * 0.0) 928 A[j][i] = 1.0 929 b.append(0.0 / scaling_matrix[i, i]) 930 j = j + 1 931 932 # Order parameters {S2, S2f, S2s}. 933 elif match('s2', spin.params[l]): 934 # 0 <= S2 <= 1. 935 A.append(zero_array * 0.0) 936 A.append(zero_array * 0.0) 937 A[j][i] = 1.0 938 A[j+1][i] = -1.0 939 b.append(0.0 / scaling_matrix[i, i]) 940 b.append(-1.0 / scaling_matrix[i, i]) 941 j = j + 2 942 943 # S2 <= S2f and S2 <= S2s. 944 if spin.params[l] == 's2': 945 for m in range(len(spin.params)): 946 if spin.params[m] == 's2f' or spin.params[m] == 's2s': 947 A.append(zero_array * 0.0) 948 A[j][i] = -1.0 949 A[j][old_i+m] = 1.0 950 b.append(0.0) 951 j = j + 1 952 953 # Correlation times {te, tf, ts}. 954 elif match('t[efs]', spin.params[l]): 955 # te, tf, ts >= 0. 956 A.append(zero_array * 0.0) 957 A[j][i] = 1.0 958 b.append(0.0 / scaling_matrix[i, i]) 959 j = j + 1 960 961 # tf <= ts. 962 if spin.params[l] == 'ts': 963 for m in range(len(spin.params)): 964 if spin.params[m] == 'tf': 965 A.append(zero_array * 0.0) 966 A[j][i] = 1.0 967 A[j][old_i+m] = -1.0 968 b.append(0.0) 969 j = j + 1 970 971 # te, tf, ts <= 2 * tm. (tf not needed because tf <= ts). 972 if upper_time_limit: 973 if not spin.params[l] == 'tf': 974 if model_type == 'mf': 975 A.append(zero_array * 0.0) 976 A[j][i] = -1.0 977 b.append(-2.0 * cdp.diff_tensor.tm / scaling_matrix[i, i]) 978 else: 979 A.append(zero_array * 0.0) 980 A[j][0] = 2.0 981 A[j][i] = -1.0 982 b.append(0.0) 983 984 j = j + 1 985 986 # Rex. 987 elif spin.params[l] == 'rex': 988 A.append(zero_array * 0.0) 989 A[j][i] = 1.0 990 b.append(0.0 / scaling_matrix[i, i]) 991 j = j + 1 992 993 # Bond length. 994 elif spin.params[l] == 'r': 995 # 0.9e-10 <= r <= 2e-10. 996 A.append(zero_array * 0.0) 997 A.append(zero_array * 0.0) 998 A[j][i] = 1.0 999 A[j+1][i] = -1.0 1000 b.append(0.9e-10 / scaling_matrix[i, i]) 1001 b.append(-2e-10 / scaling_matrix[i, i]) 1002 j = j + 2 1003 1004 # CSA. 1005 elif spin.params[l] == 'csa': 1006 # -300e-6 <= CSA <= 0. 1007 A.append(zero_array * 0.0) 1008 A.append(zero_array * 0.0) 1009 A[j][i] = 1.0 1010 A[j+1][i] = -1.0 1011 b.append(-300e-6 / scaling_matrix[i, i]) 1012 b.append(0.0 / scaling_matrix[i, i]) 1013 j = j + 2 1014 1015 # Increment i. 1016 i = i + 1 1017 1018 # Convert to numpy data structures. 1019 A = array(A, int8) 1020 b = array(b, float64) 1021 1022 return A, b
1023 1024
1025 -def units_rex():
1026 """Return the units for the Rex parameter. 1027 1028 @return: The field strength dependent Rex units. 1029 @rtype: str 1030 """ 1031 1032 # No frequency info. 1033 if not hasattr(cdp, 'frq_labels') or len(cdp.frq_labels) == 0: 1034 return '' 1035 1036 # The units. 1037 return cdp.frq_labels[0] + ' MHz'
1038