1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23 """The steady-state heteronuclear NOE API object."""
24
25
26 from math import sqrt
27 from warnings import warn
28
29
30 from lib.errors import RelaxError, RelaxNoSequenceError
31 from lib.warnings import RelaxDeselectWarning
32 from pipe_control.mol_res_spin import exists_mol_res_spin_data, spin_loop
33 from pipe_control.pipes import check_pipe
34 from specific_analyses.api_base import API_base
35 from specific_analyses.api_common import API_common
36 from specific_analyses.noe.parameter_object import Noe_params
37
38
39 -class Noe(API_base, API_common):
40 """Specific analysis API class for the steady-state heternuclear NOE analysis."""
41
42
43 instance = None
44
55
56
57 - def calculate(self, spin_id=None, scaling_matrix=None, verbosity=1, sim_index=None):
58 """Calculate the NOE and its error.
59
60 The error for each peak is calculated using the formula::
61 ___________________________________________
62 \/ {sd(sat)*I(unsat)}^2 + {sd(unsat)*I(sat)}^2
63 sd(NOE) = -----------------------------------------------
64 I(unsat)^2
65
66 @keyword spin_id: The spin identification string.
67 @type spin_id: None or str
68 @keyword scaling_matrix: The per-model list of diagonal and square scaling matrices.
69 @type scaling_matrix: list of numpy rank-2, float64 array or list of None
70 @keyword verbosity: The amount of information to print. The higher the value, the greater the verbosity.
71 @type verbosity: int
72 @keyword sim_index: The MC simulation index (unused).
73 @type sim_index: None
74 """
75
76
77 check_pipe()
78
79
80 if not hasattr(cdp, 'spectrum_type'):
81 raise RelaxError("The spectrum types have not been set.")
82
83
84 if not 'ref' in list(cdp.spectrum_type.values()) or not 'sat' in list(cdp.spectrum_type.values()):
85 raise RelaxError("The reference and saturated NOE spectra have not been loaded.")
86
87
88 for spin in spin_loop():
89
90 if not spin.select:
91 continue
92
93
94 sat = 0.0
95 sat_err2 = 0.0
96 sat_count = 0
97 ref = 0.0
98 ref_err2 = 0.0
99 ref_count = 0
100 for id in cdp.spectrum_ids:
101
102 if cdp.spectrum_type[id] == 'sat':
103 sat += spin.peak_intensity[id]
104 sat_err2 += spin.peak_intensity_err[id]**2
105 sat_count += 1
106
107
108 if cdp.spectrum_type[id] == 'ref':
109 ref += spin.peak_intensity[id]
110 ref_err2 += spin.peak_intensity_err[id]**2
111 ref_count += 1
112
113
114 sat = sat / sat_count
115 sat_err2 = sat_err2 / sat_count
116 ref = ref / ref_count
117 ref_err2 = ref_err2 / ref_count
118
119
120 spin.noe = sat / ref
121
122
123 spin.noe_err = sqrt(sat_err2 * ref**2 + ref_err2 * sat**2) / ref**2
124
125
127 """Return a vector of parameter names.
128
129 @keyword model_info: The spin container and the spin ID string from the _model_loop_spin() method.
130 @type model_info: SpinContainer instance, str
131 @return: The vector of parameter names.
132 @rtype: list of str
133 """
134
135
136 return ['noe']
137
138
140 """Deselect spins which have insufficient data to support calculation.
141
142 @keyword data_check: A flag to signal if the presence of base data is to be checked for.
143 @type data_check: bool
144 @keyword verbose: A flag which if True will allow printouts.
145 @type verbose: bool
146 """
147
148
149 if verbose:
150 print("\nOver-fit spin deselection:")
151
152
153 if not exists_mol_res_spin_data():
154 raise RelaxNoSequenceError
155
156
157 deselect_flag = False
158 all_desel = True
159 for spin, spin_id in spin_loop(return_id=True):
160
161 if not spin.select:
162 continue
163
164
165 if not hasattr(spin, 'peak_intensity'):
166 warn(RelaxDeselectWarning(spin_id, 'the absence of intensity data'))
167 spin.select = False
168 deselect_flag = True
169 continue
170
171
172 if not len(spin.peak_intensity) >= 2:
173 warn(RelaxDeselectWarning(spin_id, 'insufficient data (less than two data points)'))
174 spin.select = False
175 deselect_flag = True
176 continue
177
178
179 if not hasattr(spin, 'peak_intensity_err'):
180 warn(RelaxDeselectWarning(spin_id, 'the absence of errors'))
181 spin.select = False
182 deselect_flag = True
183 continue
184
185
186 if not len(spin.peak_intensity_err) >= 2:
187 warn(RelaxDeselectWarning(spin_id, 'missing errors (less than two error points)'))
188 spin.select = False
189 deselect_flag = True
190 continue
191
192
193 all_desel = False
194
195
196 if verbose and not deselect_flag:
197 print("No spins have been deselected.")
198
199
200 if all_desel:
201 raise RelaxError("All spins have been deselected.")
202