Author: bugman
Date: Mon Jul 16 18:23:17 2012
New Revision: 17256
URL: http://svn.gna.org/viewcvs/relax?rev=17256&view=rev
Log:
Converted the Frame order specific _minimise_setup_rdcs() to the interatomic
data design.
The _check_rdcs() method was added by directly copying the N-state model
method.
Modified:
branches/frame_order_testing/specific_fns/frame_order.py
Modified: branches/frame_order_testing/specific_fns/frame_order.py
URL:
http://svn.gna.org/viewcvs/relax/branches/frame_order_testing/specific_fns/frame_order.py?rev=17256&r1=17255&r2=17256&view=diff
==============================================================================
--- branches/frame_order_testing/specific_fns/frame_order.py (original)
+++ branches/frame_order_testing/specific_fns/frame_order.py Mon Jul 16
18:23:17 2012
@@ -47,7 +47,7 @@
from maths_fns.coord_transform import spherical_to_cartesian
from maths_fns.rotation_matrix import euler_to_R_zyz, two_vect_to_R
from physical_constants import dipolar_constant, g1H,
return_gyromagnetic_ratio
-from relax_errors import RelaxError, RelaxInfError, RelaxModelError,
RelaxNaNError, RelaxNoModelError, RelaxNoValueError, RelaxProtonTypeError,
RelaxSpinTypeError
+from relax_errors import RelaxError, RelaxInfError, RelaxModelError,
RelaxNaNError, RelaxNoModelError, RelaxNoValueError, RelaxSpinTypeError
from relax_io import open_write_file
from relax_warnings import RelaxWarning, RelaxDeselectWarning
@@ -261,6 +261,56 @@
return list
+ def _check_rdcs(self, interatom, spin1, spin2):
+ """Check if the RDCs for the given interatomic data container should
be used.
+
+ @param interatom: The interatomic data container.
+ @type interatom: InteratomContainer instance
+ @param spin1: The first spin container.
+ @type spin1: SpinContainer instance
+ @param spin2: The second spin container.
+ @type spin2: SpinContainer instance
+ @return: True if the RDCs should be used, False otherwise.
+ """
+
+ # Skip deselected spins.
+ if not spin1.select or not spin2.select:
+ return False
+
+ # Only use interatomic data containers with RDC data.
+ if not hasattr(interatom, 'rdc'):
+ return False
+
+ # RDC data exists but the interatomic vectors are missing?
+ if not hasattr(interatom, 'vector'):
+ # Throw a warning.
+ warn(RelaxWarning("RDC data exists but the interatomic vectors
are missing, skipping the spin pair '%s' and '%s'." % (interatom.spin_id1,
interatom.spin_id2)))
+
+ # Jump to the next spin.
+ return False
+
+ # Skip non-Me pseudo-atoms for the first spin.
+ if hasattr(spin1, 'members') and len(spin1.members) != 3:
+ warn(RelaxWarning("Only methyl group pseudo atoms are supported
due to their fast rotation, skipping the spin pair '%s' and '%s'." %
(interatom.spin_id1, interatom.spin_id2)))
+ return False
+
+ # Skip non-Me pseudo-atoms for the second spin.
+ if hasattr(spin2, 'members') and len(spin2.members) != 3:
+ warn(RelaxWarning("Only methyl group pseudo atoms are supported
due to their fast rotation, skipping the spin pair '%s' and '%s'." %
(interatom.spin_id1, interatom.spin_id2)))
+ return False
+
+ # Checks.
+ if not hasattr(spin1, 'isotope'):
+ raise RelaxSpinTypeError(interatom.spin_id1)
+ if not hasattr(spin2, 'isotope'):
+ raise RelaxSpinTypeError(interatom.spin_id2)
+ if not hasattr(interatom, 'r'):
+ raise RelaxNoValueError("averaged interatomic distance")
+
+ # Everything is ok.
+ return True
+
+
def _cone_pdb(self, size=30.0, file=None, dir=None, inc=40, force=False):
"""Create a PDB file containing a geometric object representing the
Frame Order cone models.
@@ -694,9 +744,10 @@
- rdc, the RDC values.
- rdc_err, the RDC errors.
- rdc_weight, the RDC weights.
- - vectors, the heteronucleus to proton
vectors.
+ - vectors, the interatomic vectors.
- rdc_const, the dipolar constants.
- @rtype: tuple of (numpy rank-2 array, numpy rank-2
array, numpy rank-2 array)
+ - absolute, the absolute value flags (as 1's
and 0's).
+ @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)
"""
# Initialise.
@@ -705,48 +756,52 @@
rdc_weight = []
unit_vect = []
rdc_const = []
+ absolute = []
# The unit vectors and RDC constants.
- for spin, spin_id in spin_loop(return_id=True):
- # Skip deselected spins.
- if not spin.select:
+ for interatom in interatomic_loop():
+ # Get the spins.
+ spin1 = return_spin(interatom.spin_id1)
+ spin2 = return_spin(interatom.spin_id2)
+
+ # RDC checks.
+ if not self._check_rdcs(interatom, spin1, spin2):
continue
- # Only use spins with RDC data.
- if not hasattr(spin, 'rdc'):
+ # Add the vectors.
+ if arg_check.is_float(interatom.vector[0], raise_error=False):
+ unit_vect.append([interatom.vector])
+ else:
+ unit_vect.append(interatom.vector)
+
+ # Gyromagnetic ratios.
+ g1 = return_gyromagnetic_ratio(spin1.isotope)
+ g2 = return_gyromagnetic_ratio(spin2.isotope)
+
+ # Calculate the RDC dipolar constant (in Hertz, and the 3 comes
from the alignment tensor), and append it to the list.
+ rdc_const.append(3.0/(2.0*pi) * dipolar_constant(g1, g2,
interatom.r))
+
+ # Fix the unit vector data structure.
+ num = None
+ for rdc_index in range(len(unit_vect)):
+ # Number of vectors.
+ if num == None:
+ if unit_vect[rdc_index] != None:
+ num = len(unit_vect[rdc_index])
continue
- # RDC data exists but the XH bond vectors are missing?
- if not hasattr(spin, 'xh_vect') and not hasattr(spin,
'bond_vect'):
- warn(RelaxWarning("RDC data exists but the XH bond vectors
are missing, skipping spin %s." % spin_id))
- continue
-
- # Append the RDC and XH vectors to the lists.
- if hasattr(spin, 'xh_vect'):
- vect = getattr(spin, 'xh_vect')
- else:
- vect = getattr(spin, 'bond_vect')
-
- # Add the bond vectors.
- if len(vect) == 1:
- unit_vect.append(vect[0])
- else:
- raise RelaxError("The spin '%s' contains more than one XH
bond vector %s." % (spin_id, vect))
-
- # Checks.
- if not hasattr(spin, 'heteronuc_type'):
- raise RelaxSpinTypeError
- if not hasattr(spin, 'proton_type'):
- raise RelaxProtonTypeError
- if not hasattr(spin, 'r'):
- raise RelaxNoValueError("bond length")
-
- # Gyromagnetic ratios.
- gx = return_gyromagnetic_ratio(spin.heteronuc_type)
- gh = return_gyromagnetic_ratio(spin.proton_type)
-
- # Calculate the RDC dipolar constant (in Hertz, and the 3 comes
from the alignment tensor), and append it to the list.
- rdc_const.append(3.0/(2.0*pi) * dipolar_constant(gx, gh, spin.r))
+ # Check.
+ if unit_vect[rdc_index] != None and len(unit_vect[rdc_index]) !=
num:
+ raise RelaxError, "The number of interatomic vectors for all
no match:\n%s" % unit_vect
+
+ # Missing unit vectors.
+ if num == None:
+ raise RelaxError, "No interatomic vectors could be found."
+
+ # Update None entries.
+ for i in range(len(unit_vect)):
+ if unit_vect[i] == None:
+ unit_vect[i] = [[None, None, None]]*num
# The RDC data.
for align_id in cdp.align_ids:
@@ -758,16 +813,21 @@
rdc.append([])
rdc_err.append([])
rdc_weight.append([])
-
- # Spin loop over the domain.
+ absolute.append([])
+
+ # Interatom loop over the domain.
id = cdp.domain[self._domain_moving()]
- for spin in spin_loop(id):
+ for interatom in interatomic_loop(id):
+ # Get the spins.
+ spin1 = return_spin(interatom.spin_id1)
+ spin2 = return_spin(interatom.spin_id2)
+
# Skip deselected spins.
- if not spin.select:
+ if not spin1.select or not spin2.select:
continue
- # Skip spins without RDC data or XH bond vectors.
- if not hasattr(spin, 'rdc') or (not hasattr(spin, 'members')
and not hasattr(spin, 'xh_vect') and not hasattr(spin, 'bond_vect')):
+ # Only use interatomic data containers with RDC and vector
data.
+ if not hasattr(interatom, 'rdc') or not hasattr(interatom,
'vector'):
continue
# Defaults of None.
@@ -776,13 +836,13 @@
# The RDC.
if sim_index != None:
- value = spin.rdc_sim[align_id][sim_index]
+ value = interatom.rdc_sim[align_id][sim_index]
else:
- value = spin.rdc[align_id]
+ value = interatom.rdc[align_id]
# The error.
- if hasattr(spin, 'rdc_err') and align_id in
spin.rdc_err.keys():
- error = spin.rdc_err[align_id]
+ if hasattr(interatom, 'rdc_err') and align_id in
interatom.rdc_err.keys():
+ error = interatom.rdc_err[align_id]
# Append the RDCs to the list.
rdc[-1].append(value)
@@ -791,10 +851,16 @@
rdc_err[-1].append(error)
# Append the weight.
- if hasattr(spin, 'rdc_weight') and align_id in
spin.rdc_weight.keys():
- rdc_weight[-1].append(spin.rdc_weight[align_id])
+ if hasattr(interatom, 'rdc_weight') and align_id in
interatom.rdc_weight.keys():
+ rdc_weight[-1].append(interatom.rdc_weight[align_id])
else:
rdc_weight[-1].append(1.0)
+
+ # Append the absolute value flag.
+ if hasattr(interatom, 'absolute_rdc') and align_id in
interatom.absolute_rdc.keys():
+ absolute[-1].append(interatom.absolute_rdc[align_id])
+ else:
+ absolute[-1].append(False)
# Convert to numpy objects.
rdc = array(rdc, float64)
@@ -802,9 +868,10 @@
rdc_weight = array(rdc_weight, float64)
unit_vect = array(unit_vect, float64)
rdc_const = array(rdc_const, float64)
+ absolute = array(absolute, float64)
# Return the data structures.
- return rdc, rdc_err, rdc_weight, unit_vect, rdc_const
+ return rdc, rdc_err, rdc_weight, unit_vect, rdc_const, absolute
def _minimise_setup_tensors(self, sim_index=None):
@@ -1118,9 +1185,9 @@
pcs, pcs_err, pcs_weight, pcs_atoms, paramag_centre, 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_vect, rdc_const = None, None, None,
None, None
+ rdcs, rdc_err, rdc_weight, rdc_vect, rdc_const, absolute_rdc = None,
None, None, None, None, None
if 'rdc' in data_types:
- rdcs, rdc_err, rdc_weight, rdc_vect, rdc_const =
self._minimise_setup_rdcs(sim_index=sim_index)
+ rdcs, rdc_err, rdc_weight, rdc_vect, rdc_const, absolute_rdc =
self._minimise_setup_rdcs(sim_index=sim_index)
# The fixed pivot point.
pivot = None
r17256 - /branches/frame_order_testing/specific_fns/frame_order.py