"Fossies" - the Fresh Open Source Software Archive

Member "relax-5.0.0/test_suite/system_tests/n_state_model.py" (2 Dec 2019, 58074 Bytes) of package /linux/privat/relax-5.0.0.src.tar.bz2:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) Python source code syntax highlighting (style: standard) with prefixed line numbers. Alternatively you can here view or download the uninterpreted source code file. See also the latest Fossies "Diffs" side-by-side code changes report for "n_state_model.py": 4.1.3_vs_5.0.0.

    1 ###############################################################################
    2 #                                                                             #
    3 # Copyright (C) 2008-2015 Edward d'Auvergne                                   #
    4 #                                                                             #
    5 # This file is part of the program relax (http://www.nmr-relax.com).          #
    6 #                                                                             #
    7 # This program is free software: you can redistribute it and/or modify        #
    8 # it under the terms of the GNU General Public License as published by        #
    9 # the Free Software Foundation, either version 3 of the License, or           #
   10 # (at your option) any later version.                                         #
   11 #                                                                             #
   12 # This program is distributed in the hope that it will be useful,             #
   13 # but WITHOUT ANY WARRANTY; without even the implied warranty of              #
   14 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the               #
   15 # GNU General Public License for more details.                                #
   16 #                                                                             #
   17 # You should have received a copy of the GNU General Public License           #
   18 # along with this program.  If not, see <http://www.gnu.org/licenses/>.       #
   19 #                                                                             #
   20 ###############################################################################
   21 
   22 # Python module imports.
   23 from math import pi
   24 from numpy import array, float64
   25 from numpy.linalg import norm
   26 from os import listdir, sep
   27 from tempfile import mkdtemp
   28 
   29 # relax module imports.
   30 from data_store import Relax_data_store; ds = Relax_data_store()
   31 import dep_check
   32 from pipe_control.align_tensor import calc_chi_tensor
   33 from pipe_control.interatomic import interatomic_loop, return_interatom
   34 from pipe_control.mol_res_spin import return_spin, spin_index_loop, spin_loop
   35 from pipe_control.pipes import get_pipe
   36 from pipe_control.rdc import return_rdc_data
   37 from status import Status; status = Status()
   38 from test_suite.system_tests.base_classes import SystemTestCase
   39 
   40 
   41 class N_state_model(SystemTestCase):
   42     """Class for testing various aspects specific to the N-state model."""
   43 
   44     def __init__(self, methodName='runTest'):
   45         """Skip some tests if scipy is not installed.
   46 
   47         @keyword methodName:    The name of the test.
   48         @type methodName:       str
   49         """
   50 
   51         # Execute the base class method.
   52         super(N_state_model, self).__init__(methodName)
   53 
   54         # Missing module.
   55         if not dep_check.scipy_module and methodName == 'test_frame_order_align_fit':
   56             # Store in the status object. 
   57             status.skipped_tests.append([methodName, 'Scipy', self._skip_type])
   58 
   59 
   60     def check_vectors(self):
   61         """Auxiliary method for checking the correct loading of bond vectors."""
   62 
   63         # The new order.
   64         ds.order_new = array([ds.order_struct[ds.order_model[ds.order_model[0]]],
   65                               ds.order_struct[ds.order_model[ds.order_model[1]]],
   66                               ds.order_struct[ds.order_model[ds.order_model[2]]]])
   67 
   68         # The atom positions.
   69         C_pos = []
   70         C_pos.append(array([6.250,   0.948,   1.968]))
   71         C_pos.append(array([6.438,  -0.139,   1.226]))
   72         C_pos.append(array([6.108,  -0.169,   0.378]))
   73 
   74         H_pos = []
   75         H_pos.append(array([7.243,   0.580,   1.676]))
   76         H_pos.append(array([7.271,  -0.291,   0.525]))
   77         H_pos.append(array([5.735,   0.003,  -0.639]))
   78 
   79         # The real vectors.
   80         vect = []
   81         for i in range(3):
   82             vect.append(H_pos[i] - C_pos[i])
   83 
   84         # Normalise.
   85         for i in range(3):
   86             vect[i] = vect[i] / norm(vect[i])
   87 
   88         # Print out.
   89         print("Structure order: %s" % ds.order_struct)
   90         print("Model order:     %s" % ds.order_model)
   91         print("New order:       %s" % ds.order_new)
   92         for i in range(3):
   93             print("\ni = %i" % i)
   94             print("The real vector:      %s" % vect[i])
   95             print("The reordered vector: %s" % vect[ds.order_new[i]])
   96             print("The loaded vector:    %s" % cdp.interatomic[0].vector[i])
   97 
   98         # Check.
   99         for i in range(3):
  100             self.assertAlmostEqual(norm(C_pos[ds.order_new[i]] - cdp.mol[0].res[0].spin[0].pos[i]), 0.0)
  101             self.assertAlmostEqual(norm(H_pos[ds.order_new[i]] - cdp.mol[0].res[0].spin[1].pos[i]), 0.0)
  102         for i in range(3):
  103             self.assertAlmostEqual(norm(vect[ds.order_new[i]] - cdp.interatomic[0].vector[i]), 0.0)
  104 
  105 
  106     def test_5_state_xz(self):
  107         """A 5-state model in the xz-plane (no pivotting of alpha).
  108 
  109         The 5 states correspond to the Euler angles (z-y-z notation):
  110             - State 1:    {0, pi/4, 0}
  111             - State 2:    {0, pi/8, 0}
  112             - State 3:    {0, 0, 0}
  113             - State 4:    {0, -pi/8, 0}
  114             - State 5:    {0, -pi/4, 0}
  115         """
  116 
  117         # Execute the script.
  118         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'5_state_xz.py')
  119 
  120         # Test the optimised probabilities.
  121         self.assertAlmostEqual(cdp.probs[0], 0.2)
  122         self.assertAlmostEqual(cdp.probs[1], 0.2)
  123         self.assertAlmostEqual(cdp.probs[2], 0.2)
  124         self.assertAlmostEqual(cdp.probs[3], 0.2)
  125         self.assertAlmostEqual(cdp.probs[4], 0.2)
  126 
  127         # Test the optimised alpha Euler angles.
  128         self.assertAlmostEqual(cdp.alpha[0], 0.0)
  129         self.assertAlmostEqual(cdp.alpha[1], 0.0)
  130         self.assertAlmostEqual(cdp.alpha[2], 0.0)
  131         self.assertAlmostEqual(cdp.alpha[3], 0.0)
  132         self.assertAlmostEqual(cdp.alpha[4], 0.0)
  133 
  134         # Test the optimised beta Euler angles.
  135         self.assertAlmostEqual(cdp.beta[0], pi/4)
  136         self.assertAlmostEqual(cdp.beta[1], pi/8)
  137         self.assertAlmostEqual(cdp.beta[2], 0.0)
  138         self.assertAlmostEqual(cdp.beta[3], -pi/8)
  139         self.assertAlmostEqual(cdp.beta[4], -pi/4)
  140 
  141         # Test the optimised gamma Euler angles.
  142         self.assertAlmostEqual(cdp.gamma[0], 0.0)
  143         self.assertAlmostEqual(cdp.gamma[1], 0.0)
  144         self.assertAlmostEqual(cdp.gamma[2], 0.0)
  145         self.assertAlmostEqual(cdp.gamma[3], 0.0)
  146         self.assertAlmostEqual(cdp.gamma[4], 0.0)
  147 
  148         # Test the chi-squared.
  149         self.assertAlmostEqual(cdp.chi2, 3.15009916529e-32)
  150 
  151 
  152     def test_A_to_chi(self):
  153         """Test the conversion of the alignment tensor to the chi tensor."""
  154 
  155         # Execute the script.
  156         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'A_to_chi.py')
  157 
  158         # Test the optimised values.
  159         for i in range(3):
  160             self.assertAlmostEqual(cdp.chi[i, i] * 1e32, cdp.chi_ref[i] * 1e32, 2)
  161 
  162 
  163     def test_CaM_IQ_tensor_fit(self):
  164         """Test the fitting of 5 alignment tensors for the published CaM-IQ RDC and PCS data."""
  165 
  166         # Execute the script.
  167         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'CaM_IQ_tensor_fit.py')
  168 
  169         # Check the optimised alignment tensor.
  170         ids = ['dy', 'tb', 'tm', 'er', 'yb', 'ho']
  171         A = [
  172             [      -5.961228899750e-04,        4.522953035367e-06,        6.168887153854e-04,        7.863257718395e-04,       -2.695142351742e-04],
  173             [      -1.383036217378e-04,       -4.939409871651e-04,        4.130289107370e-04,        7.687580236309e-04,       -3.692016717764e-04],
  174             [      -9.960927192978e-05,        4.477678617640e-04,       -4.062486453226e-04,       -4.332178003608e-04,        3.806525171855e-04],
  175             [      -1.703610649220e-05,        2.102104571469e-04,       -2.836097445400e-04,       -2.989888174606e-04,        1.326155627753e-04],
  176             [       6.087542827320e-05,        1.644473726436e-05,       -1.804561551839e-04,       -1.544765971220e-04,        1.354612889609e-04],
  177             [      -1.902819219985e-04,       -1.267359074456e-04,        2.303672688393e-04,        3.019676781386e-04,       -2.255708108877e-04]
  178         ]
  179         for i in range(len(A)):
  180             print("Checking tensor '%s'." % ids[i])
  181             self.assertAlmostEqual(cdp.align_tensors[i].Axx * 1e4, A[i][0] * 1e4)
  182             self.assertAlmostEqual(cdp.align_tensors[i].Ayy * 1e4, A[i][1] * 1e4)
  183             self.assertAlmostEqual(cdp.align_tensors[i].Axy * 1e4, A[i][2] * 1e4)
  184             self.assertAlmostEqual(cdp.align_tensors[i].Axz * 1e4, A[i][3] * 1e4)
  185             self.assertAlmostEqual(cdp.align_tensors[i].Ayz * 1e4, A[i][4] * 1e4)
  186 
  187         # Check the optimised paramagnetic position.
  188         centre = [  6.328315298868965,   8.951353997015001,  12.311128147383837]
  189         self.assertAlmostEqual(cdp.paramagnetic_centre[0], centre[0])
  190         self.assertAlmostEqual(cdp.paramagnetic_centre[1], centre[1])
  191         self.assertAlmostEqual(cdp.paramagnetic_centre[2], centre[2])
  192 
  193         # Check the minimisation stats.
  194         self.assertAlmostEqual(cdp.chi2, 496.36267335850528)
  195         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.142825408910208)
  196         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.07265356890310298)
  197 
  198 
  199     def test_absolute_rdc(self):
  200         """Test the fitting of signless RDCs."""
  201 
  202         # Execute the script.
  203         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'absolute_rdcs.py')
  204 
  205         # Test the optimised values.
  206         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000)
  207         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000)
  208         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000)
  209         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000)
  210         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000)
  211         self.assertAlmostEqual(cdp.chi2, 0.0)
  212         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  213 
  214         # The signless RDC data.
  215         rdcs = [5.59633342475, 13.31357940769, 7.03826972130, 3.39286328073, 2.09118060289, 11.44314950665, 9.06351706695, 2.33713806872, 5.81432510092, 13.10212128419, 2.52845064335, 4.70528375938, 4.07965480340, 6.28030444828, 4.69179757106, 2.34216201798, 3.89567105101, 5.51427513007, 0.72184322202, 3.81502890358, 10.88354253947, 1.66151988717, 4.29930397984, 4.46950447650, 6.99742077188, 2.27879506276, 3.64303288709, 6.83945430255, 3.19585334782]
  216 
  217         # Back calc.
  218         self.interpreter.rdc.back_calc('abs')
  219 
  220         # Check the spin data.
  221         i = 0
  222         for spin in spin_loop():
  223             # No PCS.
  224             if not hasattr(spin, 'rdc'):
  225                 continue
  226 
  227             # Check the loaded and back-calculated absolute values.
  228             self.assertAlmostEqual(spin.rdc['abs'], abs(rdcs[i]))
  229             self.assertAlmostEqual(spin.rdc_bc['abs'], abs(rdcs[i]))
  230 
  231             # Increment the spin index.
  232             i += 1
  233 
  234 
  235     def test_absolute_rdc_menthol(self):
  236         """Test the fitting of signless RDCs for menthol."""
  237 
  238         # Execute the script.
  239         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'absolute_rdcs_menthol.py')
  240 
  241         # Test the optimised values.
  242         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -2.183595975281820e-04)
  243         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, -7.014379141006286e-05)
  244         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -1.744959310458587e-05)
  245         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 3.646699595026552e-05)
  246         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, 2.592895195459969e-04)
  247         self.assertAlmostEqual(cdp.chi2, 728.32717233107246)
  248         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  249         self.assertAlmostEqual(cdp.q_rdc_norm_squared_sum, 0.7547452273747645)
  250 
  251 
  252     def test_absolute_T(self):
  253         """Test the fitting of signless T values (J+D)."""
  254 
  255         # Execute the script.
  256         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'absolute_T.py')
  257 
  258         # Test the optimised values.
  259         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -1.436586299657e-04)
  260         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, -5.004443735044e-04)
  261         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -5.017832275009e-05)
  262         self.assertAlmostEqual(cdp.align_tensors[0].Axz,  1.366097786433e-04)
  263         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -1.614772175671e-04)
  264         self.assertAlmostEqual(cdp.chi2, 311.70348701353225)
  265         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  266         self.assertAlmostEqual(cdp.q_rdc_norm_squared_sum, 0.086891848854541404)
  267 
  268         # The signless T data.
  269         T = [195.2, 205.9, 109.4, 113.8, 127.4, 156.0, 92.1, 173.2, 183.0, 174.3, 152.2, 97.3, 136.0, 129.2, 145.5, 156.0, 121.6, 128.0, 154.5, 94.2]
  270         T_bc = [
  271             195.009353539915281,
  272             205.456622836526037,
  273             112.285085032859712,
  274             113.628896345578610,
  275             127.440986667041187,
  276             155.505017063790831,
  277             94.332271833299316,
  278             172.408496922639102,
  279             181.972859458051403,
  280             173.655640981746103,
  281             153.402585241137388,
  282             92.115389822570464,
  283             139.743303992644570,
  284             131.399101601878243,
  285             146.219317894376132,
  286             153.945261372587538,
  287             119.541444938794172,
  288             126.620471670822312,
  289             155.940753902549545,
  290             90.813638474619523
  291         ]
  292 
  293         # Back calc.
  294         self.interpreter.rdc.back_calc('Gel')
  295 
  296         # Check the spin data.
  297         i = 0
  298         for interatom in interatomic_loop():
  299             # No PCS.
  300             if not hasattr(interatom, 'rdc'):
  301                 continue
  302 
  303             # Check the loaded and back-calculated absolute values.
  304             self.assertAlmostEqual(interatom.rdc['Gel'], T[i])
  305             self.assertAlmostEqual(interatom.rdc_bc['Gel'], T_bc[i], 5)
  306 
  307             # Increment the spin index.
  308             i += 1
  309 
  310 
  311     def test_align_fit(self):
  312         """Test the use of RDCs and PCSs to find the alignment tensor."""
  313 
  314         # Set the mode to use both RDCs and PCSs.
  315         ds.mode = 'all'
  316         tag = 'synth'
  317 
  318         # Execute the script.
  319         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'align_fit.py')
  320 
  321         # Test the optimised values.
  322         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000)
  323         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000)
  324         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000)
  325         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000)
  326         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000)
  327         self.assertAlmostEqual(cdp.chi2, 0.0)
  328         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  329         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  330 
  331         # The spin data.
  332         pcs = [1.0261275236, 0.75832284646, 0.65377417467, 0.88410306916, 0.83665620282, 1.887881182, 1.6564530832, 1.8489841033, -1.1143070855, -0.52863087918, -0.67600660991, -0.36996952054, -0.50720205688, -0.39889489474, -0.41237130008, -0.71313422816, -0.58642013802, -1.2160818959, -1.3990341569, -1.4084215541, -1.2007391713, -2.1392542193, -2.0165726596, -1.7623442985, -1.6437792517, -1.2415832517, -1.3008765368, -1.5872391105, -1.8060331465, -1.9063640494, -1.9817787999, -0.85264936663, -0.98332177588, -0.13370651687, -0.41762890604, -0.038212181921, -0.37986098085, 0.63582157322, 0.48346482178, 1.7566240094, 1.5694652222, 1.9914499872, 2.5316890107, 1.4559940851, 1.8661428328, 0.65003087965, 0.91690449156, 3.2096229388, 3.5547526651, 3.0579308183, 3.5933428117, 2.9062016872, 3.3750576279, 2.1848555929, 2.4769802024, 1.6466129291, 1.7719619979, 1.1373876736, 1.2182451528]
  333         i = 0
  334         for spin in spin_loop():
  335             # No PCS.
  336             if not hasattr(spin, 'pcs'):
  337                 continue
  338 
  339             # Check the value.
  340             self.assertAlmostEqual(spin.pcs[tag], pcs[i])
  341 
  342             # Check the back-calculated value.
  343             self.assertAlmostEqual(spin.pcs_bc[tag], pcs[i])
  344 
  345             # Check the simulation values.
  346             for sim_index in range(cdp.sim_number):
  347                 self.assert_(abs(spin.pcs_sim[tag][sim_index] - spin.pcs[tag]) < 6.0*spin.pcs_err[tag])
  348 
  349             # Increment the spin index.
  350             i += 1
  351 
  352         # Back calc.
  353         self.interpreter.pcs.back_calc(tag)
  354 
  355         # Check the spin data.
  356         i = 0
  357         for spin in spin_loop():
  358             # No PCS.
  359             if not hasattr(spin, 'pcs'):
  360                 continue
  361 
  362             # Check the back-calculated value.
  363             self.assertAlmostEqual(spin.pcs_bc[tag], pcs[i])
  364 
  365             # Increment the spin index.
  366             i += 1
  367 
  368 
  369     def test_align_fit_rand(self):
  370         """Test the use of randomised RDCs and PCSs to find the alignment tensor."""
  371 
  372         # Set the mode to use both RDCs and PCSs.
  373         ds.mode = 'all'
  374 
  375         # Select the randomised data.
  376         ds.rand = True
  377 
  378         # Execute the script.
  379         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'align_fit.py')
  380 
  381         # The tag.
  382         tag = 'synth'
  383 
  384         # Loop a few times.
  385         for i in range(4):
  386             # Test the optimised values (these values are from relax, so are not 100% reliable as a check).
  387             self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.000189412096996)
  388             self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.000271130087923)
  389             self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.000264898401302)
  390             self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.000284115871058)
  391             self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.000152207413184)
  392             self.assertAlmostEqual(cdp.chi2, 783.530808266)
  393             self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.063345784112045375)
  394             self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.084926009099013003)
  395 
  396             # Get a spin to check.
  397             spin = return_spin(spin_id=':114@N')
  398             spin2 = return_spin(spin_id=':114@H')
  399             interatom = return_interatom(spin_hash1=spin._hash, spin_hash2=spin2._hash)
  400 
  401             # Check the RDC and PCS values.
  402             self.assertAlmostEqual(interatom.rdc[tag], -8.9193269604999994)
  403             self.assertAlmostEqual(interatom.rdc_bc[tag], -9.1030018792821394)
  404             self.assertAlmostEqual(spin.pcs[tag], -0.41430390310999998)
  405             self.assertAlmostEqual(spin.pcs_bc[tag], -0.39723010845807194)
  406 
  407             # MC sims so next round can check if values change.
  408             if i == 0:
  409                 # Set some errors.
  410                 for interatom in interatomic_loop():
  411                     interatom.rdc_err = {tag: 1.0}
  412                 for spin in spin_loop():
  413                     spin.pcs_err = {tag: 0.1}
  414 
  415                 # MC sims.
  416                 self.interpreter.monte_carlo.setup(number=3)
  417                 self.interpreter.monte_carlo.create_data()
  418                 self.interpreter.monte_carlo.initial_values()
  419                 self.interpreter.minimise.execute('simplex', constraints=False, max_iter=5)
  420                 self.interpreter.monte_carlo.error_analysis()
  421 
  422             # Back-calc so next round can check if values change.
  423             if i == 2:
  424                 self.interpreter.rdc.back_calc(tag)
  425                 self.interpreter.pcs.back_calc(tag)
  426 
  427             # Calc Q factors so next round can check if values change.
  428             if i == 1:
  429                 self.interpreter.rdc.calc_q_factors()
  430                 self.interpreter.pcs.calc_q_factors()
  431 
  432 
  433     def test_align_fit_pcs(self):
  434         """Test the use of PCSs to find the alignment tensor."""
  435 
  436         # Set the mode to use PCSs.
  437         ds.mode = 'pcs'
  438 
  439         # Execute the script.
  440         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'align_fit.py')
  441 
  442         # Test the optimised values.
  443         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000)
  444         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000)
  445         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000)
  446         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000)
  447         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000)
  448         self.assertAlmostEqual(cdp.chi2, 0.0)
  449         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  450 
  451 
  452     def test_align_fit_pcs_rand(self):
  453         """Test the use of randomised PCSs to find the alignment tensor."""
  454 
  455         # Set the mode to use PCSs.
  456         ds.mode = 'pcs'
  457 
  458         # Select the randomised data.
  459         ds.rand = True
  460 
  461         # Execute the script.
  462         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'align_fit.py')
  463 
  464         # Test the optimised values (these values are from relax, so are not 100% reliable as a check).
  465         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.000189165581069)
  466         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.000271897288335)
  467         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.000264627388896)
  468         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.000284180080857)
  469         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.00015165641132)
  470         self.assertAlmostEqual(cdp.chi2, 756.268087443)
  471         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.063341567973121266)
  472 
  473 
  474     def test_align_fit_rdc(self):
  475         """Test the use of RDCs to find the alignment tensor."""
  476 
  477         # Set the mode to use RDCs.
  478         ds.mode = 'rdc'
  479 
  480         # Execute the script.
  481         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'align_fit.py')
  482 
  483         # Test the optimised values.
  484         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000)
  485         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000)
  486         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000)
  487         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000)
  488         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000)
  489         self.assertAlmostEqual(cdp.chi2, 0.0)
  490         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  491 
  492 
  493     def test_align_fit_rdc_rand(self):
  494         """Test the use of randomised RDCs to find the alignment tensor."""
  495 
  496         # Set the mode to use RDCs.
  497         ds.mode = 'rdc'
  498 
  499         # Select the randomised data.
  500         ds.rand = True
  501 
  502         # Execute the script.
  503         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'align_fit.py')
  504 
  505         # Test the optimised values (these are about ~10% different from Pales).
  506         # Pales:      S(zz)       S(xx-yy)      S(xy)      S(xz)      S(yz)
  507         # Pales:  -9.8124e-05 -5.2533e-04 -3.4446e-04  3.7369e-04 -1.8949e-04
  508         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.00017045)
  509         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.00024905)
  510         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.00027502)
  511         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.00029833)
  512         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.00015125)
  513         self.assertAlmostEqual(cdp.chi2, 23.5877482365)                 # Pales: 23.709
  514         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.078460000413257444)       # Pales (Q Saupe): 0.079
  515         self.assertAlmostEqual(cdp.q_rdc_norm_squared_sum, 0.14049691097282743)       # Pales (Q RDC_RMS): 0.141
  516 
  517 
  518     def test_data_copying(self):
  519         """The copying of RDC and PCS data from one pipe to another."""
  520 
  521         # Execute the script.
  522         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'data_copying.py')
  523 
  524         # Get the data pipes.
  525         orig = get_pipe('orig')
  526         new = get_pipe('new')
  527 
  528         # Check the data.
  529         self.assertEqual(orig.rdc_ids, new.rdc_ids)
  530         self.assertEqual(orig.pcs_ids, new.pcs_ids)
  531         self.assertEqual(orig.align_ids, new.align_ids)
  532 
  533         # Check the spin data.
  534         for mol_index, res_index, spin_index in spin_index_loop():
  535             # Alias the spin containers.
  536             spin_orig = orig.mol[mol_index].res[res_index].spin[spin_index]
  537             spin_new = new.mol[mol_index].res[res_index].spin[spin_index]
  538 
  539             # Skip deselected spins.
  540             if not spin_orig.select:
  541                 continue
  542 
  543             # Loop over the alignments.
  544             for id in orig.align_ids:
  545                 # RDC checks.
  546                 if hasattr(spin_orig, 'rdc'):
  547                     # Check the keys.
  548                     self.assertEqual(sorted(spin_orig.rdc.keys()), sorted(spin_new.rdc.keys()))
  549 
  550                     # Check the values.
  551                     if id in spin_orig.rdc:
  552                         self.assertEqual(spin_orig.rdc[id], spin_new.rdc[id])
  553 
  554                 # PCS checks.
  555                 if hasattr(spin_orig, 'pcs'):
  556                     # Check the keys.
  557                     self.assertEqual(sorted(spin_orig.pcs.keys()), sorted(spin_new.pcs.keys()))
  558 
  559                     # Check the values.
  560                     if id in spin_orig.pcs:
  561                         self.assertEqual(spin_orig.pcs[id], spin_new.pcs[id])
  562 
  563 
  564     def test_frame_order_align_fit(self):
  565         """Test the use of alignment tensors, RDCs and PCSs from a frame order data pipe for the N-state model."""
  566 
  567         # Execute the script.
  568         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'frame_order_align_fit.py')
  569 
  570         # The actual tensors.
  571         A_5D = []
  572         A_5D.append([1.42219822168827662867e-04, -1.44543001566521341940e-04, -7.07796211648713973798e-04, -6.01619494082773244303e-04, 2.02008007072950861996e-04])
  573         A_5D.append([3.56720663040924505435e-04, -2.68385787902088840916e-04, -1.69361406642305853832e-04, 1.71873715515064501074e-04, -3.05790155096090983822e-04])
  574         A_5D.append([2.32088908680377300801e-07, 2.08076808579168379617e-06, -2.21735465435989729223e-06, -3.74311563209448033818e-06, -2.40784858070560310370e-06])
  575         A_5D.append([-2.62495279588228071048e-04, 7.35617367964106275147e-04, 6.39754192258981332648e-05, 6.27880171180572523460e-05, 2.01197582457700226708e-04])
  576 
  577         # Check the tensors.
  578         for i in range(1):
  579             self.assertAlmostEqual(cdp.align_tensors[i].Axx, A_5D[i][0])
  580             self.assertAlmostEqual(cdp.align_tensors[i].Ayy, A_5D[i][1])
  581             self.assertAlmostEqual(cdp.align_tensors[i].Axy, A_5D[i][2])
  582             self.assertAlmostEqual(cdp.align_tensors[i].Axz, A_5D[i][3])
  583             self.assertAlmostEqual(cdp.align_tensors[i].Ayz, A_5D[i][4])
  584 
  585         # Test the optimised values.
  586         self.assertAlmostEqual(cdp.chi2, 0.0)
  587         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  588         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  589 
  590 
  591     def test_lactose_n_state_fixed(self):
  592         """The 4-state model analysis of lactose using RDCs and PCSs."""
  593 
  594         # The model.
  595         ds.model = 'fixed'
  596 
  597         # Execute the script.
  598         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'lactose_n_state.py')
  599 
  600 
  601     def test_lactose_n_state_population(self):
  602         """The 4-state model analysis of lactose using RDCs and PCSs."""
  603 
  604         # The model.
  605         ds.model = 'population'
  606 
  607         # Execute the script.
  608         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'lactose_n_state.py')
  609 
  610 
  611     def test_mc_sim_failure(self):
  612         """Test the setup of the Monte Carlo simulations
  613         
  614         This failed when this test was added, and is probably due to missing data.
  615         """
  616 
  617         # Reset and load the state.
  618         path = status.install_path + sep+'test_suite'+sep+'shared_data'+sep+'saved_states'+sep+'n_state_model_mc_fail.bz2'
  619         self.interpreter.reset()
  620         self.interpreter.state.load(path)
  621 
  622         # Monte Carlo simulations.
  623         self.interpreter.monte_carlo.setup(number=3)
  624         self.interpreter.monte_carlo.create_data()
  625         self.interpreter.monte_carlo.initial_values()
  626         self.interpreter.minimise.execute('newton', constraints=False)
  627         self.interpreter.monte_carlo.error_analysis()
  628 
  629         # Activate the optimisation of the paramagnetic centre position and try again.
  630         self.interpreter.paramag.centre(fix=False)
  631         self.interpreter.monte_carlo.setup(number=3)
  632         self.interpreter.monte_carlo.create_data()
  633         self.interpreter.monte_carlo.initial_values()
  634         self.interpreter.minimise.execute('bfgs', constraints=False, max_iter=100)
  635         self.interpreter.monte_carlo.error_analysis()
  636 
  637 
  638     def test_metal_pos_opt(self):
  639         """Test a certain algorithm for the optimisation of the lanthanide position using RDCs and PCSs (with missing data)."""
  640 
  641         # Execute the script.
  642         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'metal_pos_opt.py')
  643 
  644         # Check the metal position.
  645         self.assertAlmostEqual(cdp.paramagnetic_centre[0], -14.845)
  646         self.assertAlmostEqual(cdp.paramagnetic_centre[1], 0.969)
  647         self.assertAlmostEqual(cdp.paramagnetic_centre[2], 0.265)
  648 
  649         # The actual tensors.
  650         A_5D = []
  651         A_5D.append([1.42219822168827662867e-04, -1.44543001566521341940e-04, -7.07796211648713973798e-04, -6.01619494082773244303e-04, 2.02008007072950861996e-04])
  652         A_5D.append([3.56720663040924505435e-04, -2.68385787902088840916e-04, -1.69361406642305853832e-04, 1.71873715515064501074e-04, -3.05790155096090983822e-04])
  653         A_5D.append([2.32088908680377300801e-07, 2.08076808579168379617e-06, -2.21735465435989729223e-06, -3.74311563209448033818e-06, -2.40784858070560310370e-06])
  654         A_5D.append([-2.62495279588228071048e-04, 7.35617367964106275147e-04, 6.39754192258981332648e-05, 6.27880171180572523460e-05, 2.01197582457700226708e-04])
  655 
  656         # Check the tensors.
  657         for i in range(len(A_5D)):
  658             self.assertAlmostEqual(cdp.align_tensors[i].Axx, A_5D[i][0])
  659             self.assertAlmostEqual(cdp.align_tensors[i].Ayy, A_5D[i][1])
  660             self.assertAlmostEqual(cdp.align_tensors[i].Axy, A_5D[i][2])
  661             self.assertAlmostEqual(cdp.align_tensors[i].Axz, A_5D[i][3])
  662             self.assertAlmostEqual(cdp.align_tensors[i].Ayz, A_5D[i][4])
  663 
  664         # Test the optimised values.
  665         self.assertAlmostEqual(cdp.chi2, 0.0)
  666         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  667         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  668 
  669         # The spin data.
  670         pcs = {
  671                 'Dy': [ 0.257602207272000,  0.068336164014900,  0.110603133551000,  0.283488541078000,  0.339904760907000,  0.987042929223000,  0.442512226783000,  0.369383947069000,  0.576421001305000,  0.414840770458000,  0.332707187464000, -0.029503579631700,  0.134647155776000,  0.112481223655000,  0.355198596515000, -0.046913037723100, -0.118710220626000,  0.029273074820300,  0.302949891871000,  1.312854817740000,  0.040875625448800,  0.608281919748000,  0.023360076276000,  0.811293913028000,  0.258701382871000,  0.056025797668300,  0.609547727298000,  0.906885191563000],
  672                 'Tb': [0.644225239812000, 0.677797993460000, 0.492912191403000, 0.524119985177000, 0.600493382214000, 0.501470795483000, 0.516761193931000, 0.702967176951000, 0.914836321239000, 1.318114995620000, 1.350055099220000, 0.893477362442000, 0.999911636071000, 1.065232203160000, 0.996395971540000, 0.720982749003000, 1.374318729640000, 2.106902457600000],
  673                 'Tm': [ 0.000125066528687, -0.000564062193363, -0.000607973317902,  0.000090266635200,  0.000174865797403,  0.002488010156480,  0.000830246873289,  0.000762523870219, -0.000096933008248,  0.000742665143642, -0.000215152849719],
  674                 'Er': [-0.481445160767000, -0.361046693640000, -0.370098680342000, -0.413514467402000, -0.410802287329000, -1.081011578870000, -0.963176128222000, -0.745366702244000, -0.674570724880000, -0.751320872646000, -0.684906087274000, -0.461253969271000, -0.443680922437000, -0.344056233315000, -0.328118573270000, -0.395048353548000, -0.356220572284000, -0.324533952261000, -0.411777498713000, -0.511811581196000, -1.018565433020000, -0.959481602761000, -0.734022165690000, -0.660034918889000, -0.709085634512000, -0.878775632044000, -0.553464480425000, -0.765371835675000, -1.548006987460000]
  675         }
  676         spin_deselect_blacklist = ['pcs_bc', 'pcs_sim']
  677 
  678         # Loop over the align IDs.
  679         for tag in ['Dy', 'Tb', 'Tm', 'Er']:
  680             i = 0
  681             for spin, spin_id in spin_loop(return_id=True):
  682                 # Deselected spin.
  683                 if not spin.select:
  684                     print("Checking the deselected spin '%s'" % spin_id)
  685 
  686                     # Check that the container is clean.
  687                     for name in spin_deselect_blacklist:
  688                         print("    Checking the blacklisted object %s." % name)
  689                         self.assert_(not hasattr(spin, name))
  690 
  691                     # Skip the rest of the checks.
  692                     continue
  693 
  694                 # No PCS.
  695                 if not hasattr(spin, 'pcs'):
  696                     print("No PCS for spin '%s'." % spin_id)
  697                     continue
  698                 if not tag in spin.pcs:
  699                     print("No PCS for the '%s' alignment ID for spin '%s'." % (tag, spin_id))
  700                     continue
  701                 print("Checking '%s' PCS data for spin '%s'" % (tag, spin_id))
  702 
  703                 # Check the value.
  704                 self.assertAlmostEqual(spin.pcs[tag], pcs[tag][i])
  705 
  706                 # Check the back-calculated value.
  707                 self.assertAlmostEqual(spin.pcs_bc[tag], pcs[tag][i])
  708 
  709                 # Check the simulation values.
  710                 for sim_index in range(cdp.sim_number):
  711                     self.assert_(abs(spin.pcs_sim[tag][sim_index] - spin.pcs[tag]) < 6.0*spin.pcs_err[tag])
  712 
  713                 # Increment the spin index.
  714                 i += 1
  715 
  716         # The interatomic data.
  717         rdc = {
  718             'Dy': [-30.671893754700001, -31.307246147099999, -29.358121961700000,  -7.235977742280000, -45.562589405399997, -42.816624411200003, -41.019354747700000, -25.361318450599999, -44.1129977796],
  719             'Tm': [-0.057386340972700, -0.045650398398700, -0.074873514450400,  0.099056143214600,  0.021275817005300,  0.037132036464200,  0.047340390362400,  0.128745838536000,  0.010906407989400],
  720             'Er': [ 22.944150028900001,  23.363231565100001,  25.948323446000000,   6.955380304960000,   1.784067087050000,   7.228324193240000,   8.271072502000001,  -7.403618580470000]
  721         }
  722         interatom_deselect_blacklist = ['rdc_bc', 'rdc_sim']
  723 
  724         # Loop over the align IDs.
  725         for tag in ['Dy', 'Tm', 'Er']:
  726             i = 0
  727             for interatom in interatomic_loop():
  728                 # Deselected interatomic container.
  729                 if not interatom.select:
  730                     print("Checking the deselected interatom '%s-%s'" % (interatom.spin_id1, interatom.spin_id2))
  731 
  732                     # Check that the container is clean.
  733                     for name in interatom_deselect_blacklist:
  734                         print("    Checking the blacklisted object %s." % name)
  735                         self.assert_(not hasattr(interatom, name))
  736 
  737                     # Skip the rest of the checks.
  738                     continue
  739 
  740                 # No RDC.
  741                 if not hasattr(interatom, 'rdc'):
  742                     print("No RDC for interatom '%s-%s'." % (interatom.spin_id1, interatom.spin_id2))
  743                     continue
  744                 if not tag in interatom.rdc:
  745                     print("No RDC for the '%s' alignment ID for interatom '%s-%s'." % (tag, interatom.spin_id1, interatom.spin_id2))
  746                     continue
  747                 print("Checking '%s' RDC data for interatom '%s-%s'" % (tag, interatom.spin_id1, interatom.spin_id2))
  748 
  749                 # Check the value.
  750                 self.assertAlmostEqual(interatom.rdc[tag], rdc[tag][i])
  751 
  752                 # Check the back-calculated value.
  753                 self.assertAlmostEqual(interatom.rdc_bc[tag], rdc[tag][i])
  754 
  755                 # Check the simulation values.
  756                 for sim_index in range(cdp.sim_number):
  757                     self.assert_(abs(interatom.rdc_sim[tag][sim_index] - interatom.rdc[tag]) < 6.0*interatom.rdc_err[tag])
  758 
  759                 # Increment the interatom index.
  760                 i += 1
  761 
  762 
  763     def test_missing_data(self):
  764         """Test the use of RDCs and PCSs to find the alignment tensor with missing data."""
  765 
  766         # Execute the script.
  767         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'missing_data_test.py')
  768 
  769         # The actual tensors.
  770         A_5D = []
  771         A_5D.append([1.42219822168827662867e-04, -1.44543001566521341940e-04, -7.07796211648713973798e-04, -6.01619494082773244303e-04, 2.02008007072950861996e-04])
  772         A_5D.append([3.56720663040924505435e-04, -2.68385787902088840916e-04, -1.69361406642305853832e-04, 1.71873715515064501074e-04, -3.05790155096090983822e-04])
  773         A_5D.append([2.32088908680377300801e-07, 2.08076808579168379617e-06, -2.21735465435989729223e-06, -3.74311563209448033818e-06, -2.40784858070560310370e-06])
  774         A_5D.append([-2.62495279588228071048e-04, 7.35617367964106275147e-04, 6.39754192258981332648e-05, 6.27880171180572523460e-05, 2.01197582457700226708e-04])
  775 
  776         # Check the tensors.
  777         for i in range(len(A_5D)):
  778             self.assertAlmostEqual(cdp.align_tensors[i].Axx, A_5D[i][0])
  779             self.assertAlmostEqual(cdp.align_tensors[i].Ayy, A_5D[i][1])
  780             self.assertAlmostEqual(cdp.align_tensors[i].Axy, A_5D[i][2])
  781             self.assertAlmostEqual(cdp.align_tensors[i].Axz, A_5D[i][3])
  782             self.assertAlmostEqual(cdp.align_tensors[i].Ayz, A_5D[i][4])
  783 
  784         # Test the optimised values.
  785         self.assertAlmostEqual(cdp.chi2, 0.0)
  786         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  787         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  788 
  789 
  790     def test_monte_carlo_sims(self):
  791         """Test the Monte Carlo simulation data of fitting RDCs and PCSs."""
  792 
  793         # Execute the script.
  794         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'monte_carlo_testing.py')
  795 
  796         # Test the optimised values.
  797         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000)
  798         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000)
  799         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000)
  800         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000)
  801         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000)
  802         self.assertAlmostEqual(cdp.chi2 / 1e6, 1745860.0485368515 / 1e6, 6)
  803         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0)
  804         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  805 
  806         # The tensor key.
  807         key = 'synth'
  808 
  809         # The spin data.
  810         for spin in spin_loop():
  811             # Print out.
  812             print(spin)
  813 
  814             # Check for simulation data.
  815             self.assert_(hasattr(spin, 'pcs_sim'))
  816             self.assert_(key in spin.pcs_sim)
  817 
  818             # Check the values of the simulated data.
  819             for i in range(cdp.sim_number):
  820                 self.assertAlmostEqual(spin.pcs[key], spin.pcs_sim[key][i])
  821 
  822         # The interatomic data.
  823         for interatom in interatomic_loop():
  824             # Print out.
  825             print(interatom)
  826 
  827             # Check for simulation data.
  828             self.assert_(hasattr(interatom, 'rdc_sim'))
  829             self.assert_(key in interatom.rdc_sim)
  830 
  831             # Check the values of the simulated data.
  832             for i in range(cdp.sim_number):
  833                 self.assertAlmostEqual(interatom.rdc[key], interatom.rdc_sim[key][i], 5)
  834 
  835         # Test the optimised simluation values.
  836         for i in range(cdp.sim_number):
  837             self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000)
  838             self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000)
  839             self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000)
  840             self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000)
  841             self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000)
  842             self.assertAlmostEqual(cdp.chi2 / 1e6, 1745860.0485368515 / 1e6, 5)
  843 
  844         # Test the tensor error values.
  845         self.assertAlmostEqual(cdp.align_tensors[0].Axx_err, 0.0)
  846         self.assertAlmostEqual(cdp.align_tensors[0].Ayy_err, 0.0)
  847         self.assertAlmostEqual(cdp.align_tensors[0].Axy_err, 0.0)
  848         self.assertAlmostEqual(cdp.align_tensors[0].Axz_err, 0.0)
  849         self.assertAlmostEqual(cdp.align_tensors[0].Ayz_err, 0.0)
  850 
  851 
  852     def test_paramag_align_fit(self):
  853         """Test the use of RDCs and PCSs to find the alignment tensor."""
  854 
  855         # Execute the script.
  856         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'paramag_align_fit.py')
  857 
  858         # Test the optimised values.
  859         self.assertAlmostEqual(cdp.align_tensors[1].Axx,  0.001414718232784)
  860         self.assertAlmostEqual(cdp.align_tensors[1].Ayy,  0.001530457843766)
  861         self.assertAlmostEqual(cdp.align_tensors[1].Axy,  0.001689957281873)
  862         self.assertAlmostEqual(cdp.align_tensors[1].Axz,  0.000838692329704)
  863         self.assertAlmostEqual(cdp.align_tensors[1].Ayz, -0.000984302159683)
  864         self.assertAlmostEqual(cdp.q_factors_rdc_norm_tensor_size['Er'], 0.0, 7)
  865         self.assertAlmostEqual(cdp.q_factors_rdc_norm_squared_sum['Er'], 0.0, 7)
  866         self.assertAlmostEqual(cdp.q_factors_pcs_norm_squared_sum['Er'], 0.0, 7)
  867         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0, 7)
  868         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0, 7)
  869 
  870 
  871     def test_paramag_centre_fit(self):
  872         """Test the use of RDCs and PCSs to find the alignment tensor."""
  873 
  874         # Set the mode to use both RDCs and PCSs.
  875         ds.mode = 'all'
  876 
  877         # Execute the script.
  878         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'paramag_centre_fit.py')
  879 
  880         # Check the paramagnetic centre position.
  881         self.assertAlmostEqual(cdp.paramagnetic_centre[0], 32.555, 1)
  882         self.assertAlmostEqual(cdp.paramagnetic_centre[1], -19.130, 1)
  883         self.assertAlmostEqual(cdp.paramagnetic_centre[2], 27.775, 1)
  884 
  885         # Test the optimised values.
  886         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.351261/2000, 5)
  887         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.556994/2000, 5)
  888         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.506392/2000, 5)
  889         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.560544/2000, 5)
  890         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.286367/2000, 5)
  891         self.assertAlmostEqual(cdp.chi2, 0.0, 2)
  892         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0, 2)
  893         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0, 2)
  894 
  895 
  896     def test_pcs_back_calc(self):
  897         """Test the back-calculation of PCSs for ubiquitin."""
  898 
  899         # Execute the script.
  900         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'pcs_back_calc.py')
  901 
  902         # Test the optimised values.
  903         self.assertAlmostEqual(cdp.mol[0].res[0].spin[0].pcs_bc['A'],  0.061941887563792014)
  904         self.assertAlmostEqual(cdp.mol[0].res[1].spin[0].pcs_bc['A'], -0.077886567972081502)
  905         self.assertAlmostEqual(cdp.mol[0].res[2].spin[0].pcs_bc['A'], -0.13928519099517916)
  906 
  907 
  908     def test_pcs_fit_true_pos(self):
  909         """Test the fit of DNA PCSs at the true Ln3+ position."""
  910 
  911         # Set the Ln3+ position.
  912         ds.para_centre = 'true'
  913 
  914         # Execute the script.
  915         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'dna_pcs_fit.py')
  916 
  917         # Test the optimised values.
  918         self.assertAlmostEqual(cdp.align_tensors[0].Axx,  1.42219822168827662867e-04)
  919         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, -1.44543001566521341940e-04)
  920         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -7.07796211648713973798e-04)
  921         self.assertAlmostEqual(cdp.align_tensors[0].Axz, -6.01619494082773244303e-04)
  922         self.assertAlmostEqual(cdp.align_tensors[0].Ayz,  2.02008007072950861996e-04)
  923         self.assertAlmostEqual(cdp.chi2, 0.0)
  924         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0)
  925 
  926 
  927     def test_pcs_fit_zero_pos(self):
  928         """Test the fit of DNA PCSs at a Ln3+ position of [0, 0, 0]."""
  929 
  930         # Set the Ln3+ position.
  931         ds.para_centre = 'zero'
  932 
  933         # Execute the script.
  934         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'dna_pcs_fit.py')
  935 
  936         # Test the optimised values.
  937         self.assertAlmostEqual(cdp.align_tensors[0].Axx,  9.739588118243e-07)
  938         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, -1.077401299806e-05)
  939         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -2.321033328910e-06)
  940         self.assertAlmostEqual(cdp.align_tensors[0].Axz,  5.105903556692e-07)
  941         self.assertAlmostEqual(cdp.align_tensors[0].Ayz,  1.676638764825e-05)
  942         self.assertAlmostEqual(cdp.chi2, 2125.9562247877066)
  943         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.76065986767333704)
  944 
  945         # The chi tensor.
  946         chi_diag = calc_chi_tensor(cdp.align_tensors[0].A_diag, 799.75376122 * 1e6, 298)
  947         chi_diag = chi_diag * 1e33
  948         self.assertAlmostEqual((chi_diag[2, 2] - (chi_diag[0, 0] + chi_diag[1, 1])/2.0), -6.726159808496, 5)
  949         self.assertAlmostEqual((chi_diag[0, 0] - chi_diag[1, 1]), -3.960936794864, 6)
  950 
  951 
  952     def test_pcs_to_rdc(self):
  953         """Test the back-calculation of RDCs from a PCS derived tensor."""
  954 
  955         # Execute the script.
  956         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'pcs_to_rdc.py')
  957 
  958         # Test the values.
  959         self.assertAlmostEqual(cdp.interatomic[0].rdc_bc['A'], 4.1319413321530014)
  960         self.assertAlmostEqual(cdp.interatomic[1].rdc_bc['A'], -9.5802642470087989)
  961         self.assertAlmostEqual(cdp.interatomic[2].rdc_bc['A'], -16.244078605100817)
  962 
  963 
  964     def test_pyrotartaric_anhydride_absT(self):
  965         """Pyrotarctic anhydride alignment tensor optimisation using long range (1J, 2J & 3J) absolute T (J+D) data."""
  966 
  967         # The setup.
  968         ds.abs_data = 'T'
  969 
  970         # Execute the script.
  971         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'pyrotartaric_anhydride.py')
  972 
  973         # Test the optimised values.
  974         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.0001756305, 5)
  975         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.000278497, 5)
  976         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.000253196, 5)
  977         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.000280272, 5)
  978         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.0001431835, 5)
  979         self.assertAlmostEqual(cdp.chi2, 0.0, 2)
  980         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0, 2)
  981 
  982 
  983     def test_pyrotartaric_anhydride_mix(self):
  984         """Pyrotarctic anhydride alignment tensor optimisation using short range RDC and long range (1J, 2J & 3J) absolute T (J+D) data."""
  985 
  986         # The setup.
  987         ds.abs_data = 'mix'
  988 
  989         # Execute the script.
  990         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'pyrotartaric_anhydride.py')
  991 
  992         # Test the optimised values.
  993         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.0001756305, 5)
  994         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.000278497, 5)
  995         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.000253196, 5)
  996         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.000280272, 5)
  997         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.0001431835, 5)
  998         self.assertAlmostEqual(cdp.chi2, 0.0, 2)
  999         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0, 2)
 1000 
 1001 
 1002     def test_pyrotartaric_anhydride_rdcs(self):
 1003         """Pyrotarctic anhydride alignment tensor optimisation using long range (1J, 2J & 3J) RDC data."""
 1004 
 1005         # The setup.
 1006         ds.abs_data = 'D'
 1007 
 1008         # Execute the script.
 1009         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'pyrotartaric_anhydride.py')
 1010 
 1011         # Get the RDC data.
 1012         rdcs, rdc_err, rdc_weight, rdc_vector, rdc_dj, absolute_rdc, T_flags, j_couplings, rdc_pseudo_flags = return_rdc_data()
 1013 
 1014         # The data as it should be.
 1015         ids = ['@9 - @Q9', '@4 - @6', '@5 - @7', '@5 - @8', '@1 - @6', '@3 - @6', '@5 - @6', '@9 - @6', '@1 - @7', '@3 - @7', '@4 - @7', '@9 - @7', '@1 - @8', '@3 - @8', '@4 - @8', '@9 - @8', '@3 - @Q9', '@4 - @Q9', '@5 - @Q9']
 1016         real_rdcs = array([7.051710295953332, 14.64956993990, -9.80224941614,  4.49085022708,  1.16041049951,  0.57071216172,  3.68667449742, -2.26063357144, -4.77232431456, -0.17007443173,  2.37501105989,  1.64523216045, -0.94447557779,  0.06213688971,  1.48958862680,  0.21349779284, -1.2400897128766667, -1.8997427023766667, -0.0129325208733333], float64)
 1017 
 1018         # Check the RDC data.
 1019         for i in range(len(ids)):
 1020             print("Spin pair '%s'." % ids[i])
 1021             self.assertEqual(real_rdcs[i], rdcs[0, i])
 1022 
 1023         # Test the optimised values.
 1024         self.assertAlmostEqual(cdp.align_tensors[0].Axx, -0.0001756305, 5)
 1025         self.assertAlmostEqual(cdp.align_tensors[0].Ayy, 0.000278497, 5)
 1026         self.assertAlmostEqual(cdp.align_tensors[0].Axy, -0.000253196, 5)
 1027         self.assertAlmostEqual(cdp.align_tensors[0].Axz, 0.000280272, 5)
 1028         self.assertAlmostEqual(cdp.align_tensors[0].Ayz, -0.0001431835, 5)
 1029         self.assertAlmostEqual(cdp.chi2, 0.0, 2)
 1030         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0, 2)
 1031 
 1032 
 1033     def test_rdc_tensor(self):
 1034         """Test the calculation of an alignment tensor from RDC data."""
 1035 
 1036         # Execute the script.
 1037         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'rdc_tensor.py')
 1038 
 1039 
 1040     def test_statistics(self):
 1041         """Test the statistics user functions for the N-state model.
 1042         
 1043         This uses the 4-state model analysis of lactose using RDCs and PCSs.
 1044         """
 1045 
 1046         # Execute the script.
 1047         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'lactose_n_state.py')
 1048 
 1049         # Test the optimised chi-squared, then delete the variable.
 1050         self.assertAlmostEqual(cdp.chi2, 1051728.8810718122)
 1051         del cdp.chi2
 1052 
 1053 
 1054         # Statistics at the solution and off the solution.
 1055         dy_tensor = [
 1056                 (-1.809863649202e-05, 1.706576149818e-05, -7.790951217246e-06, -1.314676121261e-05, 4.526806452559e-08),
 1057                 (10.0e-05, 1.706576149818e-05, -7.790951217246e-06, -1.314676121261e-05, 4.526806452559e-08),
 1058         ]
 1059         chi2 = [1051728.8810718122, 5927446.9030144978]
 1060 
 1061         # The number of parameters (4 alignment tensors and 3 probabilities).
 1062         k = 5*4 + 3
 1063 
 1064         # The number of data sets (4 sets of 9 RDCs and 20 PCSs).
 1065         n = 9*4 + 20*4
 1066 
 1067         # Loop over the solution and non-solution.
 1068         for i in range(2):
 1069             # Reset the Dy tensor.
 1070             self.interpreter.align_tensor.init(tensor='Dy', align_id='Dy', params=dy_tensor[i], param_types=3)
 1071 
 1072             # Calculate the statistics.
 1073             self.interpreter.statistics.model()
 1074 
 1075             # Checks.
 1076             self.assertAlmostEqual(cdp.chi2, chi2[i])
 1077             self.assertAlmostEqual(cdp.num_params, k)
 1078             self.assertAlmostEqual(cdp.num_data_points, n)
 1079 
 1080             # Test the AIC statistic.
 1081             self.interpreter.statistics.aic()
 1082             self.assertAlmostEqual(cdp.aic, chi2[i]+2.0*k)
 1083 
 1084 
 1085     def test_stereochem_analysis(self):
 1086         """The full relative stereochemistry analysis."""
 1087 
 1088         # Create a temporary directory for all result files.
 1089         ds.tmpdir = mkdtemp()
 1090 
 1091         # Execute the script.
 1092         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'stereochem_analysis.py')
 1093 
 1094         # Check the base directory files.
 1095         files = listdir(ds.tmpdir)
 1096         for file in files:
 1097             print("Checking file %s." % file)
 1098             self.assert_(file in ['NOE_viol_S_sorted', 'ensembles_superimposed', 'RDC_PAN_dist.agr', 'Q_factors_S', 'NOE_viol_curve.agr', 'NOE_viol_dist.agr', 'RDC_PAN_curve.agr', 'NOE_viol_S', 'Q_factors_R_sorted', 'NOE_results', 'Q_factors_R', 'NOE_viol_R_sorted', 'logs', 'NOE_viol_R', 'Q_factors_S_sorted', 'RDC_PAN_results', 'correlation_plot.agr', 'correlation_plot_scaled.agr'])
 1099 
 1100         # Check the sub-directory files.
 1101         subdirs = ['ensembles_superimposed', 'logs', 'NOE_results', 'RDC_PAN_results']
 1102         files = [['S0.pdb', 'S2.pdb', 'R0.pdb', 'R1.pdb', 'S1.pdb', 'R2.pdb'],
 1103                  ['RDC_PAN_analysis.log', 'NOE_viol.log'],
 1104                  ['S_results_0.bz2', 'S_results_1.bz2', 'R_results_2.bz2', 'R_results_0.bz2', 'S_results_2.bz2', 'R_results_1.bz2'],
 1105                  ['S_results_0.bz2', 'S_results_1.bz2', 'R_results_2.bz2', 'R_results_0.bz2', 'S_results_2.bz2', 'R_results_1.bz2']]
 1106         for i in range(len(subdirs)):
 1107             for file in listdir(ds.tmpdir + sep + subdirs[i]):
 1108                 print("Checking file %s." % file)
 1109                 self.assert_(file in files[i])
 1110 
 1111 
 1112     def test_populations(self):
 1113         """Test the 'population' N-state model optimisation using RDCs and PCSs (with missing data)."""
 1114 
 1115         # Execute the script.
 1116         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'populations.py')
 1117 
 1118         # The actual tensors.
 1119         A_5D = []
 1120         A_5D.append([1.42219822168827662867e-04, -1.44543001566521341940e-04, -7.07796211648713973798e-04, -6.01619494082773244303e-04, 2.02008007072950861996e-04])
 1121         A_5D.append([3.56720663040924505435e-04, -2.68385787902088840916e-04, -1.69361406642305853832e-04, 1.71873715515064501074e-04, -3.05790155096090983822e-04])
 1122         A_5D.append([2.32088908680377300801e-07, 2.08076808579168379617e-06, -2.21735465435989729223e-06, -3.74311563209448033818e-06, -2.40784858070560310370e-06])
 1123         A_5D.append([-2.62495279588228071048e-04, 7.35617367964106275147e-04, 6.39754192258981332648e-05, 6.27880171180572523460e-05, 2.01197582457700226708e-04])
 1124 
 1125         # Check the tensors.
 1126         for i in range(len(A_5D)):
 1127             self.assertAlmostEqual(cdp.align_tensors[i].Axx, A_5D[i][0], 5)
 1128             self.assertAlmostEqual(cdp.align_tensors[i].Ayy, A_5D[i][1], 5)
 1129             self.assertAlmostEqual(cdp.align_tensors[i].Axy, A_5D[i][2], 5)
 1130             self.assertAlmostEqual(cdp.align_tensors[i].Axz, A_5D[i][3], 5)
 1131             self.assertAlmostEqual(cdp.align_tensors[i].Ayz, A_5D[i][4], 5)
 1132 
 1133         # Check the populations.
 1134         self.assertEqual(len(cdp.probs), 3)
 1135         self.assertAlmostEqual(cdp.probs[0], 0.3, 3)
 1136         self.assertAlmostEqual(cdp.probs[1], 0.6, 2)
 1137         self.assertAlmostEqual(cdp.probs[2], 0.1, 2)
 1138 
 1139         # Test the optimised values.
 1140         self.assertAlmostEqual(cdp.chi2, 0.0, 2)
 1141         self.assertAlmostEqual(cdp.q_rdc_norm_tensor_size, 0.0, 2)
 1142         self.assertAlmostEqual(cdp.q_pcs_norm_squared_sum, 0.0, 1)
 1143 
 1144 
 1145     def test_vector_loading1(self):
 1146         """Test the loading of inter-atomic vectors in the 'population' N-state model."""
 1147 
 1148         # Order.
 1149         ds.order_struct = [1, 2, 0]
 1150         ds.order_model  = [0, 1, 2]
 1151 
 1152         # Execute the script.
 1153         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'vector_loading.py')
 1154 
 1155         # Check the vectors.
 1156         self.check_vectors()
 1157 
 1158 
 1159     def test_vector_loading2(self):
 1160         """Test the loading of inter-atomic vectors in the 'population' N-state model."""
 1161 
 1162         # Order.
 1163         ds.order_struct = [0, 1, 2]
 1164         ds.order_model  = [2, 0, 1]
 1165 
 1166         # Execute the script.
 1167         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'vector_loading.py')
 1168 
 1169         # Check the vectors.
 1170         self.check_vectors()
 1171 
 1172 
 1173     def test_vector_loading3(self):
 1174         """Test the loading of inter-atomic vectors in the 'population' N-state model."""
 1175 
 1176         # Order.
 1177         ds.order_struct = [1, 0, 2]
 1178         ds.order_model  = [2, 0, 1]
 1179 
 1180         # Execute the script.
 1181         self.script_exec(status.install_path + sep+'test_suite'+sep+'system_tests'+sep+'scripts'+sep+'n_state_model'+sep+'vector_loading.py')
 1182 
 1183         # Check the vectors.
 1184         self.check_vectors()