Author: bugman
Date: Fri May 31 11:15:58 2013
New Revision: 19817
URL: http://svn.gna.org/viewcvs/relax?rev=19817&view=rev
Log:
Fixes to the specific_analyses.relax_disp modules to add support for the CR72
dispersion model.
The parameters for the CR72 model are now both correct and correctly handled.
Modified:
branches/relax_disp/specific_analyses/relax_disp/__init__.py
branches/relax_disp/specific_analyses/relax_disp/parameters.py
Modified: branches/relax_disp/specific_analyses/relax_disp/__init__.py
URL:
http://svn.gna.org/viewcvs/relax/branches/relax_disp/specific_analyses/relax_disp/__init__.py?rev=19817&r1=19816&r2=19817&view=diff
==============================================================================
--- branches/relax_disp/specific_analyses/relax_disp/__init__.py (original)
+++ branches/relax_disp/specific_analyses/relax_disp/__init__.py Fri May 31
11:15:58 2013
@@ -93,11 +93,13 @@
self.PARAMS.add('r2eff', scope='spin', default=15.0, desc='The
effective transversal relaxation rate', set='params', py_type=dict,
grace_string='\\qR\\s2,eff\\N\\Q (rad.s\\S-1\\N)', err=True, sim=True)
self.PARAMS.add('i0', scope='spin', default=10000.0, desc='The
initial intensity', py_type=dict, set='params', grace_string='\\qI\\s0\\Q',
err=True, sim=True)
self.PARAMS.add('r2', scope='spin', default=15.0, desc='The
transversal relaxation rate', set='params', py_type=list,
grace_string='\\qR\\s2\\N\\Q (rad.s\\S-1\\N)', err=True, sim=True)
+ self.PARAMS.add('pA', scope='spin', default=0.5, desc='The
population for state A', set='params', py_type=float,
grace_string='\qp\sA\N\Q', err=True, sim=True)
+ self.PARAMS.add('pB', scope='spin', default=0.5, desc='The
population for state B', set='params', py_type=float,
grace_string='\qp\sB\N\Q', err=True, sim=True)
self.PARAMS.add('phi_ex', scope='spin', default=5.0, desc='The
pA.pB.dw**2 value scaled by wH (phi_ex = pA * pB * Delta_omega**2 /
omega_H**2)', set='params', py_type=float, grace_string='\\xF\\B\\sex\\N\\q
(p\\sA\\N.p\\sB\\N.\\xDw\\B\\S2\\N / \\xw\\B\\sH\\N\\S2\\N)', err=True,
sim=True)
+ self.PARAMS.add('dw', scope='spin', default=0.0, desc='The chemical
shift difference between states A and B (in ppm)', set='params',
py_type=float, grace_string='\q\\xDw\f{}\Q (ppm)', err=True, sim=True)
self.PARAMS.add('kex', scope='spin', default=10000.0, desc='The
exchange rate', set='params', py_type=float, grace_string='\\qk\\sex\\N\\Q
(rad.s\\S-1\\N)', err=True, sim=True)
self.PARAMS.add('r2a', scope='spin', default=15.0, desc='The
transversal relaxation rate for state A', set='params', py_type=float,
grace_string='\\qR\\s2,A\\N\\Q (rad.s\\S-1\\N)', err=True, sim=True)
self.PARAMS.add('ka', scope='spin', default=10000.0, desc='The
exchange rate from state A to state B', set='params', py_type=float,
grace_string='\\qk\\sA\\N\\Q (rad.s\\S-1\\N)', err=True, sim=True)
- self.PARAMS.add('dw', scope='spin', default=1000.0, desc='The
chemical shift difference between states A and B', set='params',
py_type=float, grace_string='\\q\\xDw\\f{}\\Q (Hz)', err=True, sim=True)
self.PARAMS.add('params', scope='spin', desc='The model parameters',
py_type=list)
# Add the minimisation data.
@@ -367,35 +369,30 @@
# Only use the parameters of the first spin of the cluster.
spin = spins[0]
for i in range(len(spin.params)):
- # R2 relaxation rate (from 1 to 40 s^-1).
- if spin.params[i] == 'r2':
+ # R2 relaxation rates (from 1 to 40 s^-1).
+ if spin.params[i] in ['r2', 'r2a']:
lower.append(1.0)
upper.append(40.0)
+
+ # The population of state A.
+ elif spin.params[i] == 'pA':
+ lower.append(0.0)
+ upper.append(1.0)
# The pA.pB.dw**2/wH**2 parameter.
elif spin.params[i] == 'phi_ex':
lower.append(1e-20)
upper.append(1e-17)
- # Exchange rate.
- elif spin.params[i] == 'kex':
- lower.append(1.0)
- upper.append(100000.0)
-
- # Transversal relaxation rate for state A.
- elif spin.params[i] == 'r2a':
- lower.append(1.0)
- upper.append(20.0)
-
- # Exchange rate from state A to state B.
- elif spin.params[i] == 'ka':
- lower.append(1.0)
- upper.append(100000.0)
-
# Chemical shift difference between states A and B.
elif spin.params[i] == 'dw':
lower.append(1.0)
upper.append(10000.0)
+
+ # Exchange rates.
+ elif spin.params[i] in ['kex', 'ka']:
+ lower.append(1.0)
+ upper.append(100000.0)
# The full grid size.
grid_size = 1
@@ -828,7 +825,7 @@
params = []
for i in range(cdp.spectrometer_frq_count):
params.append('r2')
- params += ['r2a', 'ka', 'dw']
+ params += ['pa', 'dw', 'kex']
# Invalid model.
else:
@@ -1269,14 +1266,16 @@
return_data_name_doc = Desc_container("Relaxation dispersion curve
fitting data type string matching patterns")
_table = uf_tables.add_table(label="table: dispersion curve-fit data
type patterns", caption="Relaxation dispersion curve fitting data type string
matching patterns.")
_table.add_headings(["Data type", "Object name"])
- _table.add_row(["Transversal relaxation rate", "'r2'"])
- _table.add_row(["The pA.pB.dw**2/wH**2 parameter", "'phi_ex'"])
- _table.add_row(["Exchange rate", "'kex'"])
- _table.add_row(["Transversal relaxation rate for state A", "'r2a'"])
- _table.add_row(["Exchange rate from state A to state B", "'ka'"])
- _table.add_row(["Chemical shift difference between states A and B",
"'dw'"])
+ _table.add_row(["Transversal relaxation rate (rad/s)", "'r2'"])
+ _table.add_row(["Transversal relaxation rate for state A (rad/s)",
"'r2a'"])
+ _table.add_row(["Population of state A", "'pA'"])
+ _table.add_row(["Population of state B", "'pB'"])
+ _table.add_row(["The pA.pB.dw**2 parameter (ppm^2)", "'phi_ex'"])
+ _table.add_row(["Chemical shift difference between states A and B
(ppm)", "'dw'"])
+ _table.add_row(["Exchange rate (rad/s)", "'kex'"])
+ _table.add_row(["Exchange rate from state A to state B (rad/s)", "'ka'"])
_table.add_row(["Peak intensities (series)", "'intensities'"])
- _table.add_row(["CPMG pulse train frequency (series)", "'cpmg_frqs'"])
+ _table.add_row(["CPMG pulse train frequency (series, Hz)",
"'cpmg_frqs'"])
return_data_name_doc.add_table(_table.label)
Modified: branches/relax_disp/specific_analyses/relax_disp/parameters.py
URL:
http://svn.gna.org/viewcvs/relax/branches/relax_disp/specific_analyses/relax_disp/parameters.py?rev=19817&r1=19816&r2=19817&view=diff
==============================================================================
--- branches/relax_disp/specific_analyses/relax_disp/parameters.py (original)
+++ branches/relax_disp/specific_analyses/relax_disp/parameters.py Fri May 31
11:15:58 2013
@@ -116,8 +116,26 @@
else:
param_vector.append(spin.r2[i])
+ # Transversal relaxation rate for state A.
+ elif spin.params[i] == 'r2a':
+ if sim_index != None:
+ param_vector.append(spin.r2a_sim[sim_index])
+ elif spin.r2a == None:
+ param_vector.append(0.0)
+ else:
+ param_vector.append(spin.r2a)
+
+ # The pA parameter.
+ elif spin.params[i] == 'pA':
+ if sim_index != None:
+ param_vector.append(spin.pA_sim[sim_index])
+ elif spin.pA == None:
+ param_vector.append(0.0)
+ else:
+ param_vector.append(spin.pA)
+
# The pA.pB.dw**2/wH**2 parameter.
- if spin.params[i] == 'phi_ex':
+ elif spin.params[i] == 'phi_ex':
if sim_index != None:
param_vector.append(spin.phi_ex_sim[sim_index])
elif spin.phi_ex == None:
@@ -125,6 +143,15 @@
else:
param_vector.append(spin.phi_ex)
+ # Chemical shift difference between states A and B.
+ elif spin.params[i] == 'dw':
+ if sim_index != None:
+ param_vector.append(spin.dw_sim[sim_index])
+ elif spin.dw == None:
+ param_vector.append(0.0)
+ else:
+ param_vector.append(spin.dw)
+
# Exchange rate.
elif spin.params[i] == 'kex':
if sim_index != None:
@@ -134,32 +161,14 @@
else:
param_vector.append(spin.kex)
- # Transversal relaxation rate for state A.
- if spin.params[i] == 'r2a':
- if sim_index != None:
- param_vector.append(spin.r2a_sim[sim_index])
- elif spin.r2a == None:
- param_vector.append(0.0)
- else:
- param_vector.append(spin.r2a)
-
# Exchange rate from state A to state B.
- if spin.params[i] == 'ka':
+ elif spin.params[i] == 'ka':
if sim_index != None:
param_vector.append(spin.ka_sim[sim_index])
elif spin.ka == None:
param_vector.append(0.0)
else:
param_vector.append(spin.ka)
-
- # Chemical shift difference between states A and B.
- if spin.params[i] == 'dw':
- if sim_index != None:
- param_vector.append(spin.dw_sim[sim_index])
- elif spin.dw == None:
- param_vector.append(0.0)
- else:
- param_vector.append(spin.dw)
# Return a numpy array.
return array(param_vector, float64)
@@ -220,28 +229,24 @@
spin = spins[0]
for i in range(len(spin.params)):
# Transversal relaxation rate scaling.
- if spin.params[i] == 'r2':
+ if spin.params[i] in ['r2', 'r2a']:
scaling_matrix[param_index, param_index] = 10
+
+ # The population of state A.
+ elif spin.params[i] == 'pA':
+ scaling_matrix[param_index, param_index] = 1
# The pA.pB.dw**2/wH**2 parameter.
elif spin.params[i] == 'phi_ex':
scaling_matrix[param_index, param_index] = 1e-18
- # Exchange rate scaling.
- elif spin.params[i] == 'kex':
- scaling_matrix[param_index, param_index] = 10000
-
- # Transversal relaxation rate for state A scaling
- elif spin.params[i] == 'r2a':
- scaling_matrix[param_index, param_index] = 10
-
- # Exchange rate from state A to state B scaling.
- elif spin.params[i] == 'ka':
- scaling_matrix[param_index, param_index] = 10000
-
# Chemical shift difference between states A and B scaling.
elif spin.params[i] == 'dw':
- scaling_matrix[param_index, param_index] = 1000
+ scaling_matrix[param_index, param_index] = 1
+
+ # Exchange rate scaling.
+ elif spin.params[i] in ['kex', 'ka']:
+ scaling_matrix[param_index, param_index] = 10000
# Increment the parameter index.
param_index += 1
@@ -346,6 +351,20 @@
else:
spin.r2[i] = param_vector[param_index]
+ # Transversal relaxation rate for state A.
+ if spin.params[i] == 'r2a':
+ if sim_index != None:
+ spin.r2a_sim[sim_index] = param_vector[param_index]
+ else:
+ spin.r2a = param_vector[param_index]
+
+ # The population of state A.
+ if spin.params[i] == 'pA':
+ if sim_index != None:
+ spin.pA_sim[sim_index] = param_vector[param_index]
+ else:
+ spin.pA = param_vector[param_index]
+
# The pA.pB.dw**2/wH**2 parameter.
if spin.params[i] == 'phi_ex':
if sim_index != None:
@@ -353,6 +372,13 @@
else:
spin.phi_ex = param_vector[param_index]
+ # Chemical shift difference between states A and B.
+ if spin.params[i] == 'dw':
+ if sim_index != None:
+ spin.dw_sim[sim_index] = param_vector[param_index]
+ else:
+ spin.dw = param_vector[param_index]
+
# Exchange rate.
elif spin.params[i] == 'kex':
if sim_index != None:
@@ -360,13 +386,6 @@
else:
spin.kex = param_vector[param_index]
- # Transversal relaxation rate for state A.
- if spin.params[i] == 'r2a':
- if sim_index != None:
- spin.r2a_sim[sim_index] = param_vector[param_index]
- else:
- spin.r2a = param_vector[param_index]
-
# Exchange rate from state A to state B.
if spin.params[i] == 'ka':
if sim_index != None:
@@ -374,13 +393,6 @@
else:
spin.ka = param_vector[param_index]
- # Chemical shift difference between states A and B.
- if spin.params[i] == 'dw':
- if sim_index != None:
- spin.dw_sim[sim_index] = param_vector[param_index]
- else:
- spin.dw = param_vector[param_index]
-
# Increment the parameter index.
param_index = param_index + 1
@@ -394,12 +406,14 @@
The different constraints are::
R2 >= 0
- Rex >= 0
+ R2 <= -200
+ R2A >= 0
+ pA >= 0
+ pA >= pB
+ phi_ex >= 0
+ dw >= 0
kex >= 0
-
- R2A >= 0
kA >= 0
- dw >= 0
Matrix notation
@@ -407,19 +421,23 @@
In the notation A.x >= b, where A is a matrix of coefficients, x is an
array of parameter values, and b is a vector of scalars, these inequality
constraints are::
- | 1 0 0 | | R2 | | 0 |
- | | | | | |
- |-1 0 0 | | R2 | | -200 |
- | | | | | |
- | 1 0 0 | | phi | | 0 |
- | | | | | |
- | 1 0 0 | . | kex | >= | 0 |
- | | | | | |
- | 1 0 0 | | R2A | | 0 |
- | | | | | |
- | 1 0 0 | | kA | | 0 |
- | | | | | |
- | 1 0 0 | | dw | | 0 |
+ | 1 0 0 | | R2 | | 0 |
+ | | | | | |
+ |-1 0 0 | | R2 | | -200 |
+ | | | | | |
+ | 1 0 0 | | R2A | | 0 |
+ | | | | | |
+ | 1 0 0 | | pA | | 0 |
+ | | | | | |
+ | 2 0 0 | . | pA | >= | 1 |
+ | | | | | |
+ | 1 0 0 | | phi_ex | | 0 |
+ | | | | | |
+ | 1 0 0 | | dw | | 0 |
+ | | | | | |
+ | 1 0 0 | | kex | | 0 |
+ | | | | | |
+ | 1 0 0 | | kA | | 0 |
@keyword spins: The list of spin data containers for the
block.
@@ -463,8 +481,8 @@
# Only use the parameters of the first spin of the cluster.
spin = spins[0]
for k in range(len(spin.params)):
- # The transversal relaxation rate >= 0.
- if spin.params[k] == 'r2':
+ # The transversal relaxation rates (0 <= r2 <= 200).
+ if spin.params[k] in ['r2', 'r2a']:
A.append(zero_array * 0.0)
A.append(zero_array * 0.0)
A[j][i] = 1.0
@@ -473,25 +491,32 @@
b.append(-200.0 / scaling_matrix[i, i])
j += 2
- # Relaxation rates and phi_ex.
- elif spin.params[k] in ['r2a', 'phi_ex']:
- # phi_ex, R2A >= 0.
+ # The population of state A (pA >= 0 and pA >= pB).
+ elif spin.params[k] == 'pA':
+ A.append(zero_array * 0.0)
+ A.append(zero_array * 0.0)
+ A[j][i] = 1.0
+ A[j+1][i] = 2.0
+ b.append(0.0)
+ b.append(1.0 / scaling_matrix[i, i])
+ j += 2
+
+ # The pA.pB.dw**2/wH**2 parameter (phi_ex >= 0).
+ elif spin.params[k] == 'phi_ex':
A.append(zero_array * 0.0)
A[j][i] = 1.0
b.append(0.0)
j += 1
- # Exchange rates.
- elif search('^k', spin.params[k]):
- # kex, kA >= 0.
+ # Chemical exchange difference (dw >= 0).
+ elif spin.params[k] == 'dw':
A.append(zero_array * 0.0)
A[j][i] = 1.0
b.append(0.0)
j += 1
- # Chemical exchange difference.
- elif spin.params[k] == 'dw':
- # dw >= 0.
+ # Exchange rates (k >= 0).
+ elif spin.params[k] in ['kex', 'ka']:
A.append(zero_array * 0.0)
A[j][i] = 1.0
b.append(0.0)
r19817 - in /branches/relax_disp/specific_analyses/relax_disp: __init__.py parameters.py