Author: bugman
Date: Mon Jun 17 16:16:21 2013
New Revision: 20181
URL: http://svn.gna.org/viewcvs/relax?rev=20181&view=rev
Log:
The N-state model analysis can now handle RDC data of mixed D and T=J+D.
Modified:
trunk/specific_analyses/n_state_model/__init__.py
trunk/target_functions/n_state_model.py
Modified: trunk/specific_analyses/n_state_model/__init__.py
URL:
http://svn.gna.org/viewcvs/relax/trunk/specific_analyses/n_state_model/__init__.py?rev=20181&r1=20180&r2=20181&view=diff
==============================================================================
--- trunk/specific_analyses/n_state_model/__init__.py (original)
+++ trunk/specific_analyses/n_state_model/__init__.py Mon Jun 17 16:16:21 2013
@@ -496,6 +496,7 @@
- vectors, the interatomic vectors.
- rdc_const, the dipolar constants.
- absolute, the absolute value flags (as 1's
and 0's).
+ - T_flags, the flags for T = J+D type data
(as 1's and 0's).
- j_couplings, the J coupling values if the
RDC data type is set to T = J+D.
@rtype: tuple of (numpy rank-2 array, numpy rank-2
array, numpy rank-2 array, numpy rank-3 array, numpy rank-2 array, numpy
rank-2 array)
"""
@@ -507,6 +508,7 @@
unit_vect = []
rdc_const = []
absolute = []
+ T_flags = []
j_couplings = []
# The unit vectors, RDC constants, and J couplings.
@@ -573,6 +575,12 @@
rdc_weight.append([])
absolute.append([])
+ # T-type data.
+ if align_id in cdp.rdc_data_types and
cdp.rdc_data_types[align_id] == 'T':
+ T_flags.append(True)
+ else:
+ T_flags.append(False)
+
# Interatom loop.
for interatom in interatomic_loop():
# Get the spins.
@@ -588,9 +596,6 @@
continue
# Check for J couplings if the RDC data type is T = J+D.
- j_flag = False
- if align_id in cdp.rdc_data_types and
cdp.rdc_data_types[align_id] == 'T':
- j_flag = True
if not hasattr(interatom, 'j_coupling'):
continue
@@ -615,7 +620,7 @@
# The error.
if hasattr(interatom, 'rdc_err') and align_id in
interatom.rdc_err.keys():
# T values.
- if j_flag:
+ if T_flags[-1]:
error = -3.0 *
sqrt(interatom.rdc_err[align_id]**2 + interatom.j_coupling_err**2)
# D values.
@@ -633,7 +638,7 @@
# The error.
if hasattr(interatom, 'rdc_err') and align_id in
interatom.rdc_err.keys():
# T values.
- if j_flag:
+ if T_flags[-1]:
error = sqrt(interatom.rdc_err[align_id]**2 +
interatom.j_coupling_err**2)
# D values.
@@ -665,13 +670,14 @@
unit_vect = array(unit_vect, float64)
rdc_const = array(rdc_const, float64)
absolute = array(absolute, float64)
+ T_flags = array(T_flags, float64)
if opt_uses_j_couplings():
j_couplings = array(j_couplings, float64)
else:
j_couplings = None
# Return the data structures.
- return rdc, rdc_err, rdc_weight, unit_vect, rdc_const, absolute,
j_couplings
+ return rdc, rdc_err, rdc_weight, unit_vect, rdc_const, absolute,
T_flags, j_couplings
def _minimise_setup_tensors(self, sim_index=None):
@@ -832,9 +838,9 @@
pcs, pcs_err, pcs_weight, temp, frq =
self._minimise_setup_pcs(sim_index=sim_index)
# Get the data structures for optimisation using RDCs as base data
sets.
- rdcs, rdc_err, rdc_weight, rdc_vector, rdc_dj, absolute_rdc,
j_couplings = None, None, None, None, None, None, None
+ rdcs, rdc_err, rdc_weight, rdc_vector, rdc_dj, absolute_rdc,
T_flags, j_couplings = None, None, None, None, None, None, None, None
if 'rdc' in data_types:
- rdcs, rdc_err, rdc_weight, rdc_vector, rdc_dj, absolute_rdc,
j_couplings = self._minimise_setup_rdcs(sim_index=sim_index)
+ rdcs, rdc_err, rdc_weight, rdc_vector, rdc_dj, absolute_rdc,
T_flags, j_couplings = self._minimise_setup_rdcs(sim_index=sim_index)
# Get the fixed tensors.
fixed_tensors = None
@@ -863,7 +869,7 @@
centre_fixed = cdp.paramag_centre_fixed
# Set up the class instance containing the target function.
- model = N_state_opt(model=cdp.model, N=cdp.N,
init_params=param_vector, probs=probs, full_tensors=full_tensors,
red_data=red_tensor_elem, red_errors=red_tensor_err,
full_in_ref_frame=full_in_ref_frame, fixed_tensors=fixed_tensors, pcs=pcs,
rdcs=rdcs, pcs_errors=pcs_err, rdc_errors=rdc_err, j_couplings=j_couplings,
pcs_weights=pcs_weight, rdc_weights=rdc_weight, rdc_vect=rdc_vector,
temp=temp, frq=frq, dip_const=rdc_dj, absolute_rdc=absolute_rdc,
atomic_pos=atomic_pos, paramag_centre=paramag_centre,
scaling_matrix=scaling_matrix, centre_fixed=centre_fixed)
+ model = N_state_opt(model=cdp.model, N=cdp.N,
init_params=param_vector, probs=probs, full_tensors=full_tensors,
red_data=red_tensor_elem, red_errors=red_tensor_err,
full_in_ref_frame=full_in_ref_frame, fixed_tensors=fixed_tensors, pcs=pcs,
rdcs=rdcs, pcs_errors=pcs_err, rdc_errors=rdc_err, T_flags=T_flags,
j_couplings=j_couplings, pcs_weights=pcs_weight, rdc_weights=rdc_weight,
rdc_vect=rdc_vector, temp=temp, frq=frq, dip_const=rdc_dj,
absolute_rdc=absolute_rdc, atomic_pos=atomic_pos,
paramag_centre=paramag_centre, scaling_matrix=scaling_matrix,
centre_fixed=centre_fixed)
# Return the data.
return model, param_vector, data_types, scaling_matrix
Modified: trunk/target_functions/n_state_model.py
URL:
http://svn.gna.org/viewcvs/relax/trunk/target_functions/n_state_model.py?rev=20181&r1=20180&r2=20181&view=diff
==============================================================================
--- trunk/target_functions/n_state_model.py (original)
+++ trunk/target_functions/n_state_model.py Mon Jun 17 16:16:21 2013
@@ -38,7 +38,7 @@
class N_state_opt:
"""Class containing the target function of the optimisation of the
N-state model."""
- def __init__(self, model=None, N=None, init_params=None, probs=None,
full_tensors=None, red_data=None, red_errors=None, full_in_ref_frame=None,
fixed_tensors=None, pcs=None, pcs_errors=None, pcs_weights=None, rdcs=None,
rdc_errors=None, rdc_weights=None, rdc_vect=None, j_couplings=None,
temp=None, frq=None, dip_const=None, absolute_rdc=None, atomic_pos=None,
paramag_centre=None, scaling_matrix=None, centre_fixed=True):
+ def __init__(self, model=None, N=None, init_params=None, probs=None,
full_tensors=None, red_data=None, red_errors=None, full_in_ref_frame=None,
fixed_tensors=None, pcs=None, pcs_errors=None, pcs_weights=None, rdcs=None,
rdc_errors=None, rdc_weights=None, rdc_vect=None, T_flags=None,
j_couplings=None, temp=None, frq=None, dip_const=None, absolute_rdc=None,
atomic_pos=None, paramag_centre=None, scaling_matrix=None, centre_fixed=True):
"""Set up the class instance for optimisation.
The N-state models
@@ -97,6 +97,7 @@
- rdcs, the residual dipolar couplings.
- rdc_errors, the optional residual dipolar coupling errors.
- rdc_vect, the interatomic vectors.
+ - T_flags, the array of flags specifying if the data for the
given alignment is of the T = J+D type.
- j_couplings, the J couplings for the case when RDC data is of
the type T = J+D.
- dip_const, the dipolar constants.
- absolute_rdc, the flags specifying whether each RDC is
signless.
@@ -134,6 +135,8 @@
@type rdc_weights: numpy rank-2 array
@keyword rdc_vect: The unit vector lists for the RDC. The
first index must correspond to the spin pair and the second index to each
structure (its size being equal to the number of states).
@type rdc_vect: numpy rank-2 array
+ @keyword T_flags: The array of flags specifying if the
data for the given alignment is of the T = J+D type.
+ @type T_flags: numpy rank-1 array
@keyword j_couplings: The J couplings list, used when the RDC
data is of the type T = J+D. The number of elements must match the second
dimension of the rdcs argument.
@type j_couplings: numpy rank-1 array
@keyword temp: The temperature of each experiment, used
for the PCS.
@@ -163,6 +166,7 @@
self.dip_vect = rdc_vect
self.dip_const = dip_const
self.absolute_rdc = absolute_rdc
+ self.T_flags = T_flags
self.j_couplings = j_couplings
self.temp = temp
self.frq = frq
@@ -182,11 +186,6 @@
self.scaling_flag = True
else:
self.scaling_flag = False
-
- # J coupling flag.
- self.j_flag = True
- if self.j_couplings == None:
- self.j_flag = False
# The 2-domain N-state model.
if model == '2-domain':
@@ -678,7 +677,7 @@
self.rdc_theta[align_index, j] =
ave_rdc_tensor(self.dip_const[j], self.dip_vect[j], self.N,
self.A[align_index], weights=self.probs)
# Add the J coupling to convert into the
back-calculated T = J+D value.
- if self.j_flag:
+ if self.T_flags[align_index]:
self.rdc_theta[align_index, j] +=
self.j_couplings[j]
# Take the absolute value.
r20181 - in /trunk: specific_analyses/n_state_model/__init__.py target_functions/n_state_model.py