r17254 - /trunk/generic_fns/rdc.py -- July 16, 2012 - 12:02

Author: bugman
Date: Mon Jul 16 12:02:30 2012
New Revision: 17254

URL: http://svn.gna.org/viewcvs/relax?rev=17254&view=rev
Log:
Bug fix for the RDC Q factor calculations for when the dipolar constants are 
not all the same.

Now instead of a RelaxError, a warning is given and the Q factor normalised 
by 2Da^2(4 + 3R)/5 is
skipped.  This allows long range and mixed nuclear pairs to be supported.


Modified:
    trunk/generic_fns/rdc.py

Modified: trunk/generic_fns/rdc.py
URL: 
http://svn.gna.org/viewcvs/relax/trunk/generic_fns/rdc.py?rev=17254&r1=17253&r2=17254&view=diff
==============================================================================
--- trunk/generic_fns/rdc.py (original)
+++ trunk/generic_fns/rdc.py Mon Jul 16 12:02:30 2012
@@ -327,6 +327,7 @@
         interatom_count = 0
         rdc_data = False
         rdc_bc_data = False
+        norm2_flag = True
         for interatom in interatomic_loop():
             # Increment the counter.
             interatom_count += 1
@@ -357,8 +358,9 @@
 
             # Calculate the RDC dipolar constant (in Hertz, and the 3 comes 
from the alignment tensor), and append it to the list.
             dj_new = 3.0/(2.0*pi) * dipolar_constant(g1, g2, interatom.r)
-            if dj and dj_new != dj:
-                raise RelaxError("All the RDCs must come from the same 
nucleus type.")
+            if norm2_flag and dj != None and dj_new != dj:
+                warn(RelaxWarning("The dipolar constant is not the same for 
all RDCs, skipping the Q factor normalised with 2Da^2(4 + 3R)/5."))
+                norm2_flag = False
             else:
                 dj = dj_new
 
@@ -377,20 +379,25 @@
             return
 
         # Normalisation factor of 2Da^2(4 + 3R)/5.
-        D = dj * cdp.align_tensors[cdp.align_ids.index(align_id)].A_diag
-        Da = 1.0/3.0 * (D[2, 2] - (D[0, 0]+D[1, 1])/2.0)
-        Dr = 1.0/3.0 * (D[0, 0] - D[1, 1])
-        if Da == 0:
-            R = nan
+        if norm2_flag:
+            D = dj * cdp.align_tensors[cdp.align_ids.index(align_id)].A_diag
+            Da = 1.0/3.0 * (D[2, 2] - (D[0, 0]+D[1, 1])/2.0)
+            Dr = 1.0/3.0 * (D[0, 0] - D[1, 1])
+            if Da == 0:
+                R = nan
+            else:
+                R = Dr / Da
+            norm = 2.0 * (Da)**2 * (4.0 + 3.0*R**2)/5.0
+            if Da == 0.0:
+                norm = 1e-15
+
+            # The Q-factor for the alignment.
+            cdp.q_factors_rdc[align_id] = sqrt(sse / N / norm)
+
         else:
-            R = Dr / Da
-        norm = 2.0 * (Da)**2 * (4.0 + 3.0*R**2)/5.0
-        if Da == 0.0:
-            norm = 1e-15
-
-        # The Q-factor for the alignment.
-        Q = sqrt(sse / N / norm)
-        cdp.q_factors_rdc[align_id] = Q
+            cdp.q_factors_rdc[align_id] = 0.0
+
+        # The second Q-factor definition.
         cdp.q_factors_rdc_norm2[align_id] = sqrt(sse / D2_sum)
 
     # The total Q-factor.