r20010 - /branches/relax_disp/target_functions/relax_disp.py -- June 10, 2013 - 15:29

Author: bugman
Date: Mon Jun 10 15:29:17 2013
New Revision: 20010

URL: http://svn.gna.org/viewcvs/relax?rev=20010&view=rev
Log:
Created the IT99 model target function.

This is the Ishima and Torchia 1999 2-site model for all timescales with pA 
> > pB.

This commit follows step 4 of the relaxation dispersion model addition 
tutorial
(http://thread.gmane.org/gmane.science.nmr.relax.devel/3907).


Modified:
    branches/relax_disp/target_functions/relax_disp.py

Modified: branches/relax_disp/target_functions/relax_disp.py
URL: 
http://svn.gna.org/viewcvs/relax/branches/relax_disp/target_functions/relax_disp.py?rev=20010&r1=20009&r2=20010&view=diff
==============================================================================
--- branches/relax_disp/target_functions/relax_disp.py (original)
+++ branches/relax_disp/target_functions/relax_disp.py Mon Jun 10 15:29:17 
2013
@@ -29,12 +29,13 @@
 # relax module imports.
 from lib.dispersion.cr72 import r2eff_CR72
 from lib.dispersion.dpl94 import r1rho_DPL94
+from lib.dispersion.it99 import r2eff_IT99
 from lib.dispersion.lm63 import r2eff_LM63
 from lib.dispersion.m61 import r1rho_M61
 from lib.dispersion.m61b import r1rho_M61b
 from lib.errors import RelaxError
 from target_functions.chi2 import chi2
-from specific_analyses.relax_disp.variables import MODEL_CR72, MODEL_DPL94, 
MODEL_LIST_FULL, MODEL_LM63, MODEL_M61, MODEL_M61B, MODEL_NOREX, MODEL_R2EFF
+from specific_analyses.relax_disp.variables import MODEL_CR72, MODEL_DPL94, 
MODEL_IT99, MODEL_LIST_FULL, MODEL_LM63, MODEL_M61, MODEL_M61B, MODEL_NOREX, 
MODEL_R2EFF
 
 
 class Dispersion:
@@ -49,6 +50,7 @@
             - 'No Rex':  The model for no chemical exchange relaxation.
             - 'LM63':  The Luz and Meiboom (1963) 2-site fast exchange model.
             - 'CR72':  The Carver and Richards (1972) 2-site model for all 
time scales.
+            - 'IT99':  The Ishima and Torchia (1999) 2-site model for all 
time scales with skewed populations (pA >> pB).
             - 'M61':  The Meiboom (1961) 2-site fast exchange model for 
R1rho-type experiments.
             - 'DPL94':  The Davis, Perlman and London (1994) 2-site fast 
exchange model for R1rho-type experiments.
             - 'M61 skew':  The Meiboom (1961) on-resonance 2-site model with 
skewed populations (pA >> pB) for R1rho-type experiments.
@@ -123,6 +125,8 @@
             self.func = self.func_LM63
         if model == MODEL_CR72:
             self.func = self.func_CR72
+        if model == MODEL_IT99:
+            self.func = self.func_IT99
         if model == MODEL_M61:
             self.func = self.func_M61
         if model == MODEL_DPL94:
@@ -227,6 +231,54 @@
         return chi2_sum
 
 
+    def func_IT99(self, params):
+        """Target function for the Ishima and Torchia (1999) 2-site model 
for all timescales with pA >> pB.
+
+        @param params:  The vector of parameter values.
+        @type params:   numpy rank-1 float array
+        @return:        The chi-squared value.
+        @rtype:         float
+        """
+
+        # Scaling.
+        if self.scaling_flag:
+            params = dot(params, self.scaling_matrix)
+
+        # Unpack the parameter values.
+        R20 = params[:self.end_index[0]]
+        phi_ex = params[self.end_index[0]:self.end_index[1]]
+        padw2 = params[self.end_index[1]:self.end_index[2]]
+        tex = params[self.end_index[2]]
+
+        # Initialise.
+        chi2_sum = 0.0
+
+        # Loop over the spins.
+        for spin_index in range(self.num_spins):
+            # Loop over the spectrometer frequencies.
+            for frq_index in range(self.num_frq):
+                # The R20 index.
+                r20_index = frq_index + spin_index*self.num_frq
+
+                # Convert phi_ex and pa.dw^2 from ppm^2 to (rad/s)^2.
+                phi_ex_scaled = phi_ex[spin_index] * self.frqs[spin_index, 
frq_index]**2
+                padw2_scaled = padw2[spin_index] * self.frqs[spin_index, 
frq_index]**2
+
+                # Back calculate the R2eff values.
+                r2eff_IT99(r20=R20[r20_index], phi_ex=phi_ex_scaled, 
padw2=padw2_scaled, tex=tex, cpmg_frqs=self.cpmg_frqs, 
back_calc=self.back_calc[spin_index, frq_index], 
num_points=self.num_disp_points)
+
+                # For all missing data points, set the back-calculated value 
to the measured values so that it has no effect on the chi-squared value.
+                for point_index in range(self.num_disp_points):
+                    if self.missing[spin_index, frq_index, point_index]:
+                        self.back_calc[spin_index, frq_index, point_index] = 
self.values[spin_index, frq_index, point_index]
+
+                # Calculate and return the chi-squared value.
+                chi2_sum += chi2(self.values[spin_index, frq_index], 
self.back_calc[spin_index, frq_index], self.errors[spin_index, frq_index])
+
+        # Return the total chi-squared value.
+        return chi2_sum
+
+
     def func_LM63(self, params):
         """Target function for the Luz and Meiboom (1963) fast 2-site 
exchange model.