Source Code for Module lib.dispersion.ns_r1rho_3site

  1  ############################################################################### 
  2  #                                                                             # 
  3  # Copyright (C) 2000-2001 Nikolai Skrynnikov                                  # 
  4  # Copyright (C) 2000-2001 Martin Tollinger                                    # 
  5  # Copyright (C) 2013-2014 Edward d'Auvergne                                   # 
  6  #                                                                             # 
  7  # This file is part of the program relax (http://www.nmr-relax.com).          # 
  8  #                                                                             # 
  9  # This program is free software: you can redistribute it and/or modify        # 
 10  # it under the terms of the GNU General Public License as published by        # 
 11  # the Free Software Foundation, either version 3 of the License, or           # 
 12  # (at your option) any later version.                                         # 
 13  #                                                                             # 
 14  # This program is distributed in the hope that it will be useful,             # 
 15  # but WITHOUT ANY WARRANTY; without even the implied warranty of              # 
 16  # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               # 
 17  # GNU General Public License for more details.                                # 
 18  #                                                                             # 
 19  # You should have received a copy of the GNU General Public License           # 
 20  # along with this program.  If not, see <http://www.gnu.org/licenses/>.       # 
 21  #                                                                             # 
 22  ############################################################################### 
 23   
 24  # Module docstring. 
 25  """The numerical solution for the 3-site Bloch-McConnell equations for R1rho-type data, the U{NS R1rho 3-site linear<http://wiki.nmr-relax.com/NS_R1rho_3-site_linear>} and U{NS R1rho 3-site<http://wiki.nmr-relax.com/NS_R1rho_3-site>} model. 
 26   
 27  Description 
 28  =========== 
 29   
 30  This is the model of the numerical solution for the 3-site Bloch-McConnell equations.  It originates from the funNumrho.m file from the Skrynikov & Tollinger code (the sim_all.tar file U{https://web.archive.org/web/https://gna.org/support/download.php?file_id=18404} attached to U{https://web.archive.org/web/https://gna.org/task/?7712#comment5}). 
 31   
 32   
 33  References 
 34  ========== 
 35   
 36  The solution has been modified to use the from presented in: 
 37   
 38      - Korzhnev, D. M., Orekhov, V. Y., and Kay, L. E. (2005).  Off-resonance R(1rho) NMR studies of exchange dynamics in proteins with low spin-lock fields:  an application to a Fyn SH3 domain.  I{J. Am. Chem. Soc.}, B{127}, 713-721. (U{DOI: 10.1021/ja0446855<http://dx.doi.org/10.1021/ja0446855>}). 
 39   
 40   
 41  Links 
 42  ===== 
 43   
 44  More information on the NS R1rho 3-site linear model can be found in the: 
 45   
 46      - U{relax wiki<http://wiki.nmr-relax.com/NS_R1rho_3-site_linear>}, 
 47      - U{relax manual<http://www.nmr-relax.com/manual/NS_3_site_linear_R1_model.html>}, 
 48      - U{relaxation dispersion page of the relax website<http://www.nmr-relax.com/analyses/relaxation_dispersion.html#NS_R1rho_3-site_linear>}. 
 49   
 50  More information on the NS R1rho 3-site model can be found in the: 
 51   
 52      - U{relax wiki<http://wiki.nmr-relax.com/NS_R1rho_3-site>}, 
 53      - U{relax manual<http://www.nmr-relax.com/manual/NS_3_site_R1_model.html>}, 
 54      - U{relaxation dispersion page of the relax website<http://www.nmr-relax.com/analyses/relaxation_dispersion.html#NS_R1rho_3-site>}. 
 55  """ 
 56   
 57  # Python module imports. 
 58  from math import atan2, cos, log, sin 
 59  from numpy import dot 
 60   
 61  # relax module imports. 
 62  from lib.dispersion.ns_matrices import rr1rho_3d_3site 
 63  from lib.float import isNaN 
 64  from lib.linear_algebra.matrix_exponential import matrix_exponential 
 65   
 66   
 67 -def ns_r1rho_3site(M0=None, matrix=None, r1rho_prime=None, omega=None, offset=None, r1=0.0, pA=None, pB=None, pC=None, dw_AB=None, dw_AC=None, k_AB=None, k_BA=None, k_BC=None, k_CB=None, k_AC=None, k_CA=None, spin_lock_fields=None, relax_time=None, inv_relax_time=None, back_calc=None, num_points=None): 
 68      """The 3-site numerical solution to the Bloch-McConnell equation for R1rho data. 
 69   
 70      This function calculates and stores the R1rho values. 
 71   
 72   
 73      @keyword M0:                This is a vector that contains the initial magnetizations corresponding to the A and B state transverse magnetizations. 
 74      @type M0:                   numpy float64, rank-1, 7D array 
 75      @keyword matrix:            A numpy array to be populated to create the evolution matrix. 
 76      @type matrix:               numpy rank-2, 9D float64 array 
 77      @keyword r1rho_prime:       The R1rho_prime parameter value (R1rho with no exchange). 
 78      @type r1rho_prime:          float 
 79      @keyword omega:             The chemical shift for the spin in rad/s. 
 80      @type omega:                float 
 81      @keyword offset:            The spin-lock offsets for the data. 
 82      @type offset:               numpy rank-1 float array 
 83      @keyword r1:                The R1 relaxation rate. 
 84      @type r1:                   float 
 85      @keyword pA:                The population of state A. 
 86      @type pA:                   float 
 87      @keyword pB:                The population of state B. 
 88      @type pB:                   float 
 89      @keyword pC:                The population of state C. 
 90      @type pC:                   float 
 91      @keyword dw_AB:             The chemical exchange difference between states A and B in rad/s. 
 92      @type dw_AB:                float 
 93      @keyword dw_AC:             The chemical exchange difference between states A and C in rad/s. 
 94      @type dw_AC:                float 
 95      @keyword k_AB:              The rate of exchange from site A to B (rad/s). 
 96      @type k_AB:                 float 
 97      @keyword k_BA:              The rate of exchange from site B to A (rad/s). 
 98      @type k_BA:                 float 
 99      @keyword k_BC:              The rate of exchange from site B to C (rad/s). 
100      @type k_BC:                 float 
101      @keyword k_CB:              The rate of exchange from site C to B (rad/s). 
102      @type k_CB:                 float 
103      @keyword k_AC:              The rate of exchange from site A to C (rad/s). 
104      @type k_AC:                 float 
105      @keyword k_CA:              The rate of exchange from site C to A (rad/s). 
106      @type k_CA:                 float 
107      @keyword spin_lock_fields:  The R1rho spin-lock field strengths (in rad.s^-1). 
108      @type spin_lock_fields:     numpy rank-1 float array 
109      @keyword relax_time:        The total relaxation time period for each spin-lock field strength (in seconds). 
110      @type relax_time:           float 
111      @keyword inv_relax_time:    The inverse of the relaxation time period for each spin-lock field strength (in inverse seconds).  This is used for faster calculations. 
112      @type inv_relax_time:       float 
113      @keyword back_calc:         The array for holding the back calculated R2eff values.  Each element corresponds to one of the CPMG nu1 frequencies. 
114      @type back_calc:            numpy rank-1 float array 
115      @keyword num_points:        The number of points on the dispersion curve, equal to the length of the tcp and back_calc arguments. 
116      @type num_points:           int 
117      """ 
118   
119      # Repetitive calculations (to speed up calculations). 
120      Wa = omega                  # Larmor frequency for state A [s^-1]. 
121      Wb = omega + dw_AB          # Larmor frequency for state B [s^-1]. 
122      Wc = omega + dw_AC          # Larmor frequency for state C [s^-1]. 
123      W = pA*Wa + pB*Wb + pC*Wc   # Population-averaged Larmor frequency [s^-1]. 
124      dA = Wa - offset            # Offset of spin-lock from A. 
125      dB = Wb - offset            # Offset of spin-lock from B. 
126      dC = Wc - offset            # Offset of spin-lock from C. 
127      d = W - offset              # Offset of spin-lock from population-average. 
128   
129      # Loop over the time points, back calculating the R2eff values. 
130      for i in range(num_points): 
131          # The matrix that contains all the contributions to the evolution, i.e. relaxation, exchange and chemical shift evolution. 
132          rr1rho_3d_3site(matrix=matrix, R1=r1, r1rho_prime=r1rho_prime, pA=pA, pB=pB, pC=pC, wA=dA, wB=dB, wC=dC, w1=spin_lock_fields[i], k_AB=k_AB, k_BA=k_BA, k_BC=k_BC, k_CB=k_CB, k_AC=k_AC, k_CA=k_CA) 
133   
134          # The following lines rotate the magnetization previous to spin-lock into the weff frame. 
135          theta = atan2(spin_lock_fields[i], dA) 
136          M0[0] = sin(theta)    # The A state initial X magnetisation. 
137          M0[2] = cos(theta)    # The A state initial Z magnetisation. 
138   
139          # This matrix is a propagator that will evolve the magnetization with the matrix R. 
140          Rexpo = matrix_exponential(matrix*relax_time) 
141   
142          # Magnetization evolution. 
143          MA = dot(M0, dot(Rexpo, M0)) 
144   
145          # The next lines calculate the R1rho using a two-point approximation, i.e. assuming that the decay is mono-exponential. 
146          if MA <= 0.0 or isNaN(MA): 
147              back_calc[i] = 1e99 
148          else: 
149              back_calc[i]= -inv_relax_time * log(MA) 
150