GCC Code Coverage Report

./
File: src/XpertMass/IsotopicDataLibraryHandler.cpp
Date: 2024-08-24 11:26:06
Lines:
74/121
61.2%
Functions:
5/7
71.4%
Branches:
46/150
30.7%
Line Branch Exec Source
1 /* BEGIN software license
2 *
3 * MsXpertSuite - mass spectrometry software suite
4 * -----------------------------------------------
5 * Copyright (C) 2009--2020 Filippo Rusconi
6 *
7 * http://www.msxpertsuite.org
8 *
9 * This file is part of the MsXpertSuite project.
10 *
11 * The MsXpertSuite project is the successor of the massXpert project. This
12 * project now includes various independent modules:
13 *
14 * - massXpert, model polymer chemistries and simulate mass spectrometric data;
15 * - mineXpert, a powerful TIC chromatogram/mass spectrum viewer/miner;
16 *
17 * This program is free software: you can redistribute it and/or modify
18 * it under the terms of the GNU General Public License as published by
19 * the Free Software Foundation, either version 3 of the License, or
20 * (at your option) any later version.
21 *
22 * This program is distributed in the hope that it will be useful,
23 * but WITHOUT ANY WARRANTY; without even the implied warranty of
24 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 * GNU General Public License for more details.
26 *
27 * You should have received a copy of the GNU General Public License
28 * along with this program. If not, see <http://www.gnu.org/licenses/>.
29 *
30 * END software license
31 */
32
33
34 /////////////////////// Qt includes
35 #include <QDebug>
36
37
38 /////////////////////// IsoSpec
39 #include <IsoSpec++/isoSpec++.h>
40 #include <IsoSpec++/element_tables.h>
41
42
43 /////////////////////// Local includes
44 #include "globals.hpp"
45 #include "PeakCentroid.hpp"
46 #include "IsotopicDataLibraryHandler.hpp"
47
48
49 namespace MsXpS
50 {
51
52 namespace libXpertMass
53 {
54
55
56 /*!
57 \class MsXpS::libXpertMass::IsotopicDataLibraryHandler
58 \inmodule libXpertMass
59 \ingroup PolChemDefBuildingdBlocks
60 \inheaderfile IsotopicDataLibraryHandler.hpp
61
62 \brief The IsotopicDataLibraryHandler class handles \l{IsotopicData} from the
63 IsoSpec element data tables directly from the library's data. These are the
64 reference, pristine, \e{unmodified}, isotopic data.
65
66 The IsoSpec element data tables are the following:
67
68 \list
69 \li elem_table_ID
70 \li elem_table_atomicNo
71 \li elem_table_mass
72 \li elem_table_massNo
73 \li elem_table_extraNeutrons
74 \li elem_table_element
75 \li elem_table_symbol
76 \li elem_table_Radioactive
77 \li elem_table_probability
78 \li elem_table_log_probability
79 \endlist
80
81 The data tables are all of the same length and the data in each row of a
82 given table matches the contents of that same row in all the other tables. For
83 example, the first two rows of table elem_table_ID are:
84
85 1
86
87 1
88
89 These two rows match the same rows in elem_table_mass:
90
91 1.00782503227
92
93 2.01410177819
94
95 and the the same rows in elem_table_element:
96
97 "hydrogen"
98
99 "hydrogen"
100
101 By reading, row-by-row, the data from the same row number in each one of the
102 tables, one constructs a fully qualified \l{Isotope}.
103
104 \sa IsotopicDataUserConfigHandler, IsotopicDataManualConfigHandler
105 */
106
107 /*!
108 \typedef IsotopicDataLibraryHandlerSPtr
109 \relates IsotopicDataLibraryHandler
110
111 Synonym for std::shared_ptr<IsotopicDataLibraryHandler>.
112 */
113
114 /*!
115 \typedef IsotopicDataLibraryHandlerCstSPtr
116 \relates IsotopicDataLibraryHandler
117
118 Synonym for std::shared_ptr<const IsotopicDataLibraryHandler>.
119 */
120
121
122 /*!
123 \brief Constructs the \l{IsotopicDataLibraryHandler}.
124
125 The instance will have empty member data.
126 */
127 41 IsotopicDataLibraryHandler::IsotopicDataLibraryHandler()
128
1/2
✓ Branch 2 taken 41 times.
✗ Branch 3 not taken.
 41 : IsotopicDataBaseHandler()
129 {
130 41 }
131
132 /*!
133 \brief Constructs the \l{IsotopicDataLibraryHandler}.
134
135 The instance will have its isotopic data member pointing to \a
136 isotopic_data_sp.
137 */
138 ✗ IsotopicDataLibraryHandler::IsotopicDataLibraryHandler(
139 ✗ IsotopicDataSPtr isotopic_data_sp)
140 ✗ : IsotopicDataBaseHandler(isotopic_data_sp)
141 {
142 ✗ }
143
144
145 // IsotopicDataLibraryHandler::IsotopicDataLibraryHandler(const QString
146 // &file_name) : IsotopicDataBaseHandler(file_name)
147 //{
148 //}
149
150
151 // IsotopicDataLibraryHandler::IsotopicDataLibraryHandler(
152 // IsotopicDataSPtr isotopic_data_sp, const QString &file_name)
153 //: IsotopicDataBaseHandler(isotopic_data_sp, file_name)
154 //{
155 //}
156
157 /*!
158 \brief Destructs the \l{IsotopicDataLibraryHandler}.
159
160 Nothing is explicitely deleted in the destructor.
161 */
162 90 IsotopicDataLibraryHandler::~IsotopicDataLibraryHandler()
163 {
164 // qDebug();
165 90 }
166
167 // NOT documented on purpose.
168 std::size_t
169 ✗ IsotopicDataLibraryHandler::loadData([[maybe_unused]] const QString &filename)
170 {
171 ✗ std::size_t non_isotope_skipped_items = 0;
172 ✗ return loadData(non_isotope_skipped_items);
173 }
174
175 /*!
176 \brief Loads isotopic data directly from IsoSpec library' element data tables.
177
178 The member isotopic data are cleared before setting new data read from the
179 library's element data tables.
180
181 The code iterates, row-by-row, in the all the tables and extracts
182 the data to fill in the Isotope data:
183
184 \code
185 IsotopeSPtr isotope_sp =
186 std::make_shared<Isotope>(IsoSpec::elem_table_ID[iter],
187 QString(IsoSpec::elem_table_element[iter]),
188 QString(IsoSpec::elem_table_symbol[iter]),
189 IsoSpec::elem_table_atomicNo[iter],
190 IsoSpec::elem_table_mass[iter],
191 IsoSpec::elem_table_massNo[iter],
192 IsoSpec::elem_table_extraNeutrons[iter],
193 IsoSpec::elem_table_probability[iter],
194 IsoSpec::elem_table_log_probability[iter],
195 IsoSpec::elem_table_Radioactive[iter]);
196 \endcode
197
198 Note that some rows in the table are populated with non-isotope data, like
199 electron, missing_electron or protonation (twice). When these rows are
200 encountered, the \a non_isotope_skipped_items is incremented.
201
202 Returns the count of \l{Isotope}s that were allocated and stored in the
203 msp_isotopicData member.
204 */
205 std::size_t
206 48 IsotopicDataLibraryHandler::loadData(std::size_t &non_isotope_skipped_items)
207 {
208
209 // We need to allocate one Isotope instance for each element
210 // in the various arrays in the IsoSpec++ source code header file.
211
212 // extern const int elem_table_ID[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
213 // extern const int elem_table_atomicNo[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
214 // extern const double
215 // elem_table_probability[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
216 // extern const double elem_table_mass[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
217 // extern const int elem_table_massNo[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
218 // extern const int
219 // elem_table_extraNeutrons[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
220 // extern const char*
221 // elem_table_element[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES]; extern const
222 // char* elem_table_symbol[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES]; extern const
223 // bool elem_table_Radioactive[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES]; extern
224 // const double
225 // elem_table_log_probability[ISOSPEC_NUMBER_OF_ISOTOPIC_ENTRIES];
226
227 // Big sanity check, all the arrays must be the same length!
228 48 std::size_t array_length = checkConsistency();
229
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(array_length < 1)
230 ✗ return false;
231
232 // qDebug() << "The size of the tables:" << array_length;
233
234 // Clear all the data, since this function might be called multiple times.
235 48 msp_isotopicData->clear();
236
237 48 std::size_t loaded_isotopes = 0;
238
2/2
✓ Branch 0 taken 14016 times.
✓ Branch 1 taken 48 times.
 14064 for(std::size_t iter = 0; iter < array_length; ++iter)
239 {
240
241
1/2
✓ Branch 1 taken 14016 times.
✗ Branch 2 not taken.
 14016 QString elem_element = QString(IsoSpec::elem_table_element[iter]);
242
243 // These are the last items in the various tables. We do not handle them
244 // at the moment.
245
8/8
✓ Branch 1 taken 13968 times.
✓ Branch 2 taken 48 times.
✓ Branch 4 taken 13920 times.
✓ Branch 5 taken 48 times.
✓ Branch 6 taken 96 times.
✓ Branch 7 taken 13824 times.
✓ Branch 8 taken 192 times.
✓ Branch 9 taken 13824 times.
 27936 if(elem_element == "electron" || elem_element == "missing electron" ||
246 13920 elem_element == "protonation" /* occurs twice */)
247 {
248 192 ++non_isotope_skipped_items;
249
250 192 continue;
251 }
252
253 IsotopeSPtr isotope_sp =
254 13824 std::make_shared<Isotope>(IsoSpec::elem_table_ID[iter],
255
1/2
✓ Branch 1 taken 13824 times.
✗ Branch 2 not taken.
 27648 QString(IsoSpec::elem_table_element[iter]),
256 13824 QString(IsoSpec::elem_table_symbol[iter]),
257 13824 IsoSpec::elem_table_atomicNo[iter],
258 13824 IsoSpec::elem_table_mass[iter],
259 13824 IsoSpec::elem_table_massNo[iter],
260 13824 IsoSpec::elem_table_extraNeutrons[iter],
261 13824 IsoSpec::elem_table_probability[iter],
262 13824 IsoSpec::elem_table_log_probability[iter],
263
2/4
✓ Branch 1 taken 13824 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 13824 times.
✗ Branch 5 not taken.
 27648 IsoSpec::elem_table_Radioactive[iter]);
264
265 // We do not want to update the mono/avg maps each time we load an
266 // isotope. We'll call the relevant function later.
267
1/2
✓ Branch 3 taken 13824 times.
✗ Branch 4 not taken.
 13824 msp_isotopicData->appendNewIsotope(isotope_sp, false);
268 13824 ++loaded_isotopes;
269
2/2
✓ Branch 2 taken 13824 times.
✓ Branch 3 taken 192 times.
 14016 }
270
271 // Sanity check
272 48 if((loaded_isotopes + non_isotope_skipped_items) !=
273
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 IsoSpec::isospec_number_of_isotopic_entries)
274 ✗ qFatal("Error loading the isotopic data from IsoSpec++'s tables.");
275
276 // Now ask that the mono/avg mass maps be updated.
277
1/2
✗ Branch 2 not taken.
✓ Branch 3 taken 48 times.
 48 if(!msp_isotopicData->updateMassMaps())
278 ✗ qFatal("Programming error. Failed to update the mass maps.");
279
280 // qDebug() << "Done loading data with :" << msp_isotopicData->size()
281 // << "isotopes in the isotopic data.";
282
283 48 return msp_isotopicData->size();
284 }
285
286 /*!
287 \brief Write all the IsotopicData to \a file_name.
288
289 If \a file_name is empty, m_fileName is tried. If both are empty, the
290 function returns 0. If any one of the file names are correct (file_name takes
291 precedence over m_fileName), then m_fileName is set to that file name.
292
293 The format of the file consists in a single line of data per \l{Isotope} as
294 created using the Isotope::toString() function. Each isotope is output to
295 its own line.
296
297 Returns the count of \l{Isotope}s written to file or 0 if the file does not
298 exist or is not readable.
299
300 \sa Isotope::Isotope(const QString &text)
301 */
302 std::size_t
303 4 IsotopicDataLibraryHandler::writeData(const QString &file_name)
304 {
305 // Although the isotopic data were loaded from the IsoSpec library tables, we
306 // might be willing to store these data to a file.
307
308
2/6
✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 4 times.
 4 if(file_name.isEmpty() && m_fileName.isEmpty())
309 ✗ return 0;
310
311 4 QString temp_file_name;
312
313 // The passed filename takes precedence over the member datum. So copy
314 // that file name to the member datum.
315
316
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 if(!file_name.isEmpty())
317 4 temp_file_name = file_name;
318 else
319 ✗ temp_file_name = m_fileName;
320
321
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 QFile file(temp_file_name);
322
323 // qDebug() << "File name to write to:" << temp_file_name;
324
325
2/4
✓ Branch 2 taken 4 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 4 times.
 4 if(!file.open(QIODevice::WriteOnly | QIODevice::Text))
326 {
327 ✗ qDebug("Failed to open file for writing.");
328 ✗ return 0;
329 }
330
331
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 QTextStream out(&file);
332
333 out
334
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 << "# This file contains isotopic data in a format that can accommodate\n";
335 out
336
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 << "# comments in the form of lines beginning with the '#' character.\n\n";
337
338 4 std::size_t isotope_count = 0;
339
340
2/2
✓ Branch 7 taken 1152 times.
✓ Branch 8 taken 4 times.
 1156 for(auto item : msp_isotopicData->m_isotopes)
341 {
342
2/4
✓ Branch 2 taken 1152 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 1152 times.
✗ Branch 6 not taken.
 1152 out << item->toString();
343 // We need to add it because toString() does not terminate the line with
344 // a new line character.
345
1/2
✓ Branch 1 taken 1152 times.
✗ Branch 2 not taken.
 1152 out << "\n";
346
347 1152 ++isotope_count;
348 1152 }
349
350
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 out.flush();
351
352
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 file.close();
353
354 // Now we know that temp_file_name is fine. Store into m_fileName.
355 4 m_fileName = temp_file_name;
356
357 4 return isotope_count;
358 4 }
359
360 /*!
361 \brief Checks the consistency in all the IsoSpec library's different isotopic
362 data tables.
363
364 This function essentially verifies that each table has the same row count as
365 all the other ones.
366
367 Returns the count of isotopes in the isotopic data.
368 */
369 std::size_t
370 48 IsotopicDataLibraryHandler::checkConsistency()
371 {
372 48 std::size_t array_length = sizeof(IsoSpec::elem_table_atomicNo) /
373 sizeof(IsoSpec::elem_table_atomicNo[0]);
374
375 // qDebug() << "The array length is:" << array_length;
376
377 // All the tables in the header file of the IsoSpec library must
378 // have exactly the same size.
379
380
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(IsoSpec::isospec_number_of_isotopic_entries != array_length)
381 {
382 ✗ qFatal("Found corruption: the size of arrays is not like expected.");
383 }
384
385 // Now test each table one by one.
386 48 std::size_t tested_length = sizeof(IsoSpec::elem_table_probability) /
387 sizeof(IsoSpec::elem_table_probability[0]);
388
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
389 {
390 ✗ qDebug()
391 ✗ << "Found corruption: at least two arrays are not of the same length."
392 ✗ << "tested_length:" << tested_length;
393
394 ✗ return 0;
395 }
396
397 48 tested_length =
398 sizeof(IsoSpec::elem_table_mass) / sizeof(IsoSpec::elem_table_mass[0]);
399
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
400 {
401 ✗ qDebug()
402 ✗ << "Found corruption: at least two arrays are not of the same length."
403 ✗ << "tested_length:" << tested_length;
404
405 ✗ return 0;
406 }
407
408 48 tested_length =
409 sizeof(IsoSpec::elem_table_massNo) / sizeof(IsoSpec::elem_table_massNo[0]);
410
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
411 {
412 ✗ qDebug()
413 ✗ << "Found corruption: at least two arrays are not of the same length."
414 ✗ << "tested_length:" << tested_length;
415
416 ✗ return 0;
417 }
418
419 48 tested_length = sizeof(IsoSpec::elem_table_extraNeutrons) /
420 sizeof(IsoSpec::elem_table_extraNeutrons[0]);
421
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
422 {
423 ✗ qDebug()
424 ✗ << "Found corruption: at least two arrays are not of the same length."
425 ✗ << "tested_length:" << tested_length;
426
427 ✗ return 0;
428 }
429
430 48 tested_length = sizeof(IsoSpec::elem_table_element) /
431 sizeof(IsoSpec::elem_table_element[0]);
432
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
433 {
434 ✗ qDebug()
435 ✗ << "Found corruption: at least two arrays are not of the same length."
436 ✗ << "tested_length:" << tested_length;
437
438 ✗ return 0;
439 }
440
441 48 tested_length =
442 sizeof(IsoSpec::elem_table_symbol) / sizeof(IsoSpec::elem_table_symbol[0]);
443
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
444 {
445 ✗ qDebug()
446 ✗ << "Found corruption: at least two arrays are not of the same length."
447 ✗ << "tested_length:" << tested_length;
448
449 ✗ return 0;
450 }
451
452 48 tested_length = sizeof(IsoSpec::elem_table_Radioactive) /
453 sizeof(IsoSpec::elem_table_Radioactive[0]);
454
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
455 {
456 ✗ qDebug()
457 ✗ << "Found corruption: at least two arrays are not of the same length."
458 ✗ << "tested_length:" << tested_length;
459
460 ✗ return 0;
461 }
462
463 48 tested_length = sizeof(IsoSpec::elem_table_log_probability) /
464 sizeof(IsoSpec::elem_table_log_probability[0]);
465
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 48 times.
 48 if(tested_length != array_length)
466 {
467 ✗ qDebug()
468 ✗ << "Found corruption: at least two arrays are not of the same length."
469 ✗ << "tested_length:" << tested_length;
470
471 ✗ return 0;
472 }
473
474 48 return tested_length;
475 }
476
477
478 } // namespace libXpertMass
479
480 } // namespace MsXpS
481
482
483 #if 0
484
485 Example from IsoSpec.
486
487 const int elementNumber = 2;
488 const int isotopeNumbers[2] = {2,3};
489
490 const int atomCounts[2] = {2,1};
491
492
493 const double hydrogen_masses[2] = {1.00782503207, 2.0141017778};
494 const double oxygen_masses[3] = {15.99491461956, 16.99913170, 17.9991610};
495
496 const double* isotope_masses[2] = {hydrogen_masses, oxygen_masses};
497
498 const double hydrogen_probs[2] = {0.5, 0.5};
499 const double oxygen_probs[3] = {0.5, 0.3, 0.2};
500
501 const double* probs[2] = {hydrogen_probs, oxygen_probs};
502
503 IsoLayeredGenerator iso(Iso(elementNumber, isotopeNumbers, atomCounts,
504 isotope_masses, probs), 0.99);
505
506 #endif
507