Author: bugman
Date: Mon Nov 18 17:35:22 2013
New Revision: 21504
URL: http://svn.gna.org/viewcvs/relax?rev=21504&view=rev
Log:
Simplified the 'MMQ 2-site' dispersion model target function.
The r2eff_mmq_2site_sq_dq_zq() and r2eff_mmq_2site_mq() functions from
lib.dispersion.mmq_2site are
now aliased by the experiment_type_setup() target function method. Both
functions now have matching
arguments.
Modified:
branches/relax_disp/lib/dispersion/mmq_2site.py
branches/relax_disp/target_functions/relax_disp.py
Modified: branches/relax_disp/lib/dispersion/mmq_2site.py
URL:
http://svn.gna.org/viewcvs/relax/branches/relax_disp/lib/dispersion/mmq_2site.py?rev=21504&r1=21503&r2=21504&view=diff
==============================================================================
--- branches/relax_disp/lib/dispersion/mmq_2site.py (original)
+++ branches/relax_disp/lib/dispersion/mmq_2site.py Mon Nov 18 17:35:22 2013
@@ -47,7 +47,7 @@
@type R20A: float
@keyword R20B: The transverse, spin-spin relaxation rate for
state B.
@type R20B: float
- @keyword dw: The combined chemical exchange difference
parameters between states A and B in rad/s. This can be any combination of
dwX and dwH.
+ @keyword dw: The combined chemical exchange difference
parameters between states A and B in rad/s. This can be any combination of
dw and dwH.
@type dw: float
@keyword k_AB: The rate of exchange from site A to B (rad/s).
@type k_AB: float
@@ -62,7 +62,7 @@
matrix[1, 1] = -k_BA + 1.j*dw - R20B
-def r2eff_mmq_2site_mq(M0=None, m1=None, m2=None, R20A=None, R20B=None,
pA=None, pB=None, dw=None, dwH=None, k_AB=None, k_BA=None, inv_tcpmg=None,
tcp=None, back_calc=None, num_points=None, n=None):
+def r2eff_mmq_2site_mq(M0=None, m1=None, m2=None, R20A=None, R20B=None,
pA=None, pB=None, dw=None, dwH=None, k_AB=None, k_BA=None, inv_tcpmg=None,
tcp=None, back_calc=None, num_points=None, power=None, n=None):
"""The 2-site numerical solution to the Bloch-McConnell equation for MQ
data.
The notation used here comes from:
@@ -108,6 +108,8 @@
@type back_calc: numpy rank-1 float array
@keyword num_points: The number of points on the dispersion curve,
equal to the length of the tcp and back_calc arguments.
@type num_points: int
+ @keyword power: The matrix exponential power array.
+ @type power: numpy int16, rank-1 array
@keyword n: The n value whereby one CPMG block is defined at
2n.
@type n: numpy int16, rank-1 array
"""
@@ -185,7 +187,7 @@
back_calc[i]= -inv_tcpmg * log(Mx / pA)
-def r2eff_mmq_2site_sq_dq_zq(M0=None, F_vector=array([1, 0], float64),
m1=None, m2=None, R20A=None, R20B=None, pA=None, pB=None, dwX=None,
k_AB=None, k_BA=None, inv_tcpmg=None, tcp=None, back_calc=None,
num_points=None, power=None):
+def r2eff_mmq_2site_sq_dq_zq(M0=None, F_vector=array([1, 0], float64),
m1=None, m2=None, R20A=None, R20B=None, pA=None, pB=None, dw=None, dwH=None,
k_AB=None, k_BA=None, inv_tcpmg=None, tcp=None, back_calc=None,
num_points=None, power=None, n=None):
"""The 2-site numerical solution to the Bloch-McConnell equation for SQ,
ZQ, and DQ data.
The notation used here comes from:
@@ -213,6 +215,8 @@
@type pB: float
@keyword dw: The combined chemical exchange difference
between states A and B in rad/s. It should be set to dwH for 1H SQ data, dw
for heteronuclear SQ data, dwH-dw for ZQ data, and dwH+dw for DQ data.
@type dw: float
+ @keyword dwH: Unused - this is simply to match the
r2eff_mmq_2site_mq() function arguments.
+ @type dwH: float
@keyword k_AB: The rate of exchange from site A to B (rad/s).
@type k_AB: float
@keyword k_BA: The rate of exchange from site B to A (rad/s).
@@ -227,11 +231,13 @@
@type num_points: int
@keyword power: The matrix exponential power array.
@type power: numpy int16, rank-1 array
+ @keyword n: The n value whereby one CPMG block is defined at
2n.
+ @type n: numpy int16, rank-1 array
"""
# Populate the m1 and m2 matrices (only once per function call for
speed).
- populate_matrix(matrix=m1, R20A=R20A, R20B=R20B, dw=dwX, k_AB=k_AB,
k_BA=k_BA)
- populate_matrix(matrix=m2, R20A=R20A, R20B=R20B, dw=-dwX, k_AB=k_AB,
k_BA=k_BA)
+ populate_matrix(matrix=m1, R20A=R20A, R20B=R20B, dw=dw, k_AB=k_AB,
k_BA=k_BA)
+ populate_matrix(matrix=m2, R20A=R20A, R20B=R20B, dw=-dw, k_AB=k_AB,
k_BA=k_BA)
# Loop over the time points, back calculating the R2eff values.
for i in range(num_points):
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=21504&r1=21503&r2=21504&view=diff
==============================================================================
--- branches/relax_disp/target_functions/relax_disp.py (original)
+++ branches/relax_disp/target_functions/relax_disp.py Mon Nov 18 17:35:22
2013
@@ -49,7 +49,7 @@
from lib.dispersion.tsmfk01 import r2eff_TSMFK01
from lib.errors import RelaxError
from target_functions.chi2 import chi2
-from specific_analyses.relax_disp.variables import EXP_TYPE_CPMG,
EXP_TYPE_DQ_CPMG, EXP_TYPE_MQ_CPMG, EXP_TYPE_PROTON_MQ_CPMG,
EXP_TYPE_PROTON_SQ_CPMG, EXP_TYPE_R1RHO, EXP_TYPE_ZQ_CPMG, MODEL_CR72,
MODEL_CR72_FULL, MODEL_DPL94, MODEL_IT99, MODEL_LIST_CPMG,
MODEL_LIST_CPMG_NUM, MODEL_LIST_FULL, MODEL_LIST_MQ_CPMG, MODEL_LIST_R1RHO,
MODEL_LM63, MODEL_LM63_3SITE, MODEL_M61, MODEL_M61B, MODEL_MMQ_2SITE,
MODEL_MP05, MODEL_MQ_CR72, MODEL_NOREX, MODEL_NS_CPMG_2SITE_3D,
MODEL_NS_CPMG_2SITE_3D_FULL, MODEL_NS_CPMG_2SITE_EXPANDED,
MODEL_NS_CPMG_2SITE_STAR, MODEL_NS_CPMG_2SITE_STAR_FULL,
MODEL_NS_R1RHO_2SITE, MODEL_R2EFF, MODEL_TAP03, MODEL_TP02, MODEL_TSMFK01
+from specific_analyses.relax_disp.variables import EXP_TYPE_CPMG,
EXP_TYPE_DQ_CPMG, EXP_TYPE_MQ_CPMG, EXP_TYPE_PROTON_MQ_CPMG,
EXP_TYPE_PROTON_SQ_CPMG, EXP_TYPE_R1RHO, EXP_TYPE_ZQ_CPMG, MODEL_CR72,
MODEL_CR72_FULL, MODEL_DPL94, MODEL_IT99, MODEL_LIST_CPMG,
MODEL_LIST_CPMG_NUM, MODEL_LIST_FULL, MODEL_LIST_MMQ, MODEL_LIST_MQ_CPMG,
MODEL_LIST_R1RHO, MODEL_LM63, MODEL_LM63_3SITE, MODEL_M61, MODEL_M61B,
MODEL_MMQ_2SITE, MODEL_MP05, MODEL_MQ_CR72, MODEL_NOREX,
MODEL_NS_CPMG_2SITE_3D, MODEL_NS_CPMG_2SITE_3D_FULL,
MODEL_NS_CPMG_2SITE_EXPANDED, MODEL_NS_CPMG_2SITE_STAR,
MODEL_NS_CPMG_2SITE_STAR_FULL, MODEL_NS_R1RHO_2SITE, MODEL_R2EFF,
MODEL_TAP03, MODEL_TP02, MODEL_TSMFK01
class Dispersion:
@@ -482,8 +482,18 @@
# The number of experiments.
self.num_exp = len(self.exp_types)
- # The non-combined data models.
- if not self.model in [MODEL_MMQ_2SITE]:
+ # The MMQ combined data type models.
+ if self.model in MODEL_LIST_MMQ:
+ # Alias the r2eff functions.
+ self.r2eff_mmq = []
+ for exp_index in range(self.num_exp):
+ if self.exp_types[exp_index] in [EXP_TYPE_CPMG,
EXP_TYPE_PROTON_SQ_CPMG, EXP_TYPE_DQ_CPMG, EXP_TYPE_ZQ_CPMG]:
+ self.r2eff_mmq.append(r2eff_mmq_2site_sq_dq_zq)
+ elif self.exp_types[exp_index] in [EXP_TYPE_MQ_CPMG,
EXP_TYPE_PROTON_MQ_CPMG]:
+ self.r2eff_mmq.append(r2eff_mmq_2site_mq)
+
+ # The single data type models.
+ else:
# Remove the first dimension of the data structures.
self.values = self.values[0]
self.errors = self.errors[0]
@@ -923,6 +933,7 @@
self.M0[1] = pB
# Initialise.
+ aliased_dwH = 0.0
chi2_sum = 0.0
# Loop over the experiment types.
@@ -938,19 +949,24 @@
dw_frq = dw[spin_index] *
self.frqs[exp_index][spin_index][frq_index]
dwH_frq = dwH[spin_index] *
self.frqs[exp_index][spin_index][frq_index]
+ # Alias the dw frequency combinations.
+ if self.exp_types[exp_index] == EXP_TYPE_CPMG:
+ aliased_dw = dw_frq
+ elif self.exp_types[exp_index] ==
EXP_TYPE_PROTON_SQ_CPMG:
+ aliased_dw = dwH_frq
+ elif self.exp_types[exp_index] == EXP_TYPE_DQ_CPMG:
+ aliased_dw = dwH_frq + dw_frq
+ elif self.exp_types[exp_index] == EXP_TYPE_ZQ_CPMG:
+ aliased_dw = dwH_frq - dw_frq
+ elif self.exp_types[exp_index] == EXP_TYPE_MQ_CPMG:
+ aliased_dw = dw_frq
+ aliased_dwH = dwH_frq
+ elif self.exp_types[exp_index] ==
EXP_TYPE_PROTON_MQ_CPMG:
+ aliased_dw = dwH_frq
+ aliased_dwH = dw_frq
+
# Back calculate the R2eff values for each experiment
type.
- if self.exp_types[exp_index] == EXP_TYPE_CPMG:
- r2eff_mmq_2site_sq_dq_zq(M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dwX=dw_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_time,
tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
power=self.power[exp_index][frq_index])
- elif self.exp_types[exp_index] ==
EXP_TYPE_PROTON_SQ_CPMG:
- r2eff_mmq_2site_sq_dq_zq(M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dwX=dwH_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_time,
tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
power=self.power[exp_index][frq_index])
- elif self.exp_types[exp_index] == EXP_TYPE_DQ_CPMG:
- r2eff_mmq_2site_sq_dq_zq(M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dwX=dwH_frq+dw_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_time,
tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
power=self.power[exp_index][frq_index])
- elif self.exp_types[exp_index] == EXP_TYPE_ZQ_CPMG:
- r2eff_mmq_2site_sq_dq_zq(M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dwX=dwH_frq-dw_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_time,
tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
power=self.power[exp_index][frq_index])
- elif self.exp_types[exp_index] == EXP_TYPE_MQ_CPMG:
- r2eff_mmq_2site_mq(M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dw=dw_frq, dwH=dwH_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_time,
tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
n=self.n[exp_index][frq_index])
- elif self.exp_types[exp_index] ==
EXP_TYPE_PROTON_MQ_CPMG:
- r2eff_mmq_2site_mq(M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dw=dwH_frq, dwH=dw_frq, k_AB=k_AB, k_BA=k_BA, inv_tcpmg=self.inv_relax_time,
tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
n=self.n[exp_index][frq_index])
+ self.r2eff_mmq[exp_index](M0=self.M0, m1=self.m1,
m2=self.m2, R20A=R20[r20_index], R20B=R20[r20_index], pA=pA, pB=pB,
dw=aliased_dw, dwH=aliased_dwH, k_AB=k_AB, k_BA=k_BA,
inv_tcpmg=self.inv_relax_time, tcp=self.tau_cpmg[exp_index][frq_index],
back_calc=self.back_calc[exp_index][spin_index][frq_index],
num_points=self.num_disp_points[exp_index][frq_index],
power=self.power[exp_index][frq_index], n=self.n[exp_index][frq_index])
# 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[exp_index][frq_index]):
r21504 - in /branches/relax_disp: lib/dispersion/mmq_2site.py target_functions/relax_disp.py