GCC Code Coverage Report


source/XpertMassCore/src/
File: source/XpertMassCore/src/IsotopicClusterShaper.cpp
Date: 2025-11-20 01:41:33
Lines:
0/118
0.0%
Functions:
0/19
0.0%
Branches:
0/128
0.0%

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 /////////////////////// StdLib includes
35 #include <cmath>
36 #include <algorithm>
37 #include <limits> // for std::numeric_limits
38
39
40 /////////////////////// Qt includes
41
42
43 /////////////////////// pappsomspp includes
44 #include <pappsomspp/core/massspectrum/massspectrum.h>
45 #include <pappsomspp/core/processing/combiners/massspectrumpluscombiner.h>
46 #include <pappsomspp/core/processing/combiners/tracepluscombiner.h>
47 #include <pappsomspp/core/trace/trace.h>
48 #include <pappsomspp/core/processing/filters/filternormalizeintensities.h>
49 #include <pappsomspp/core/utils.h>
50
51
52 /////////////////////// Local includes
53 #include "MsXpS/libXpertMassCore/globals.hpp"
54 #include "MsXpS/libXpertMassCore/Utils.hpp"
55 #include "MsXpS/libXpertMassCore/MassDataCborMassSpectrumHandler.hpp"
56 #include "MsXpS/libXpertMassCore/IsotopicClusterShaper.hpp"
57 #include "MsXpS/libXpertMassCore/MassPeakShaper.hpp"
58
59 namespace MsXpS
60 {
61
62 namespace libXpertMassCore
63 {
64
65
66 /*!
67 \class MsXpS::libXpertMassCore::IsotopicClusterShaper
68 \inmodule libXpertMassCore
69 \ingroup XpertMassCoreMassCalculations
70 \inheaderfile IsotopicClusterShaper.hpp
71
72 \brief The IsotopicClusterShaper class provides the features needed to shape
73 sets of (peak centroid m/z, intensity) pairs associated to a given charge
74 into a mass spectral pappso;:Trace.
75
76 Each set of (peak centroid m/z, intensity) pairs corresponds to an isotopic
77 cluster that is associated to a charge.
78
79 The configuration of the peak shaping process is held in a specific \l
80 MassPeakShaperConfig class.
81
82 The output of the computation is a pappso::Trace obtained by combining all the
83 different shapes obtained for the different peak centroids of all the sets of
84 (peak centroid m/z, intensity) pairs. If binning was requested, the obtained
85 Trace is the result of a combination accounting for the required bin size,
86 otherwise the obtained Trace is the result of the mere addition of all the
87 points in the different traces.
88
89 \sa MassPeakShaperConfig
90 */
91
92 /*!
93 \variable \
94 MsXpS::libXpertMassCore::IsotopicClusterShaper::m_isotopicClusterChargePairs
95
96 \brief Vector of pappso::Trace instances in pair with charges.
97 */
98
99 /*!
100 \variable
101 MsXpS::libXpertMassCore::IsotopicClusterShaper::mp_massPeakShaperConfig
102
103 \brief The configuration driving the mass peak shaping process.
104 */
105
106 /*!
107 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_mapTrace
108
109 \brief The map relating a m/z value to its intensity.
110
111 This map is a variant of pappso::Trace that is designed to allow for easy mass
112 spectrum combination. It is generally used only for computations and is
113 converted to a pappso::Trace once all the computations have been carried out.
114 */
115
116 /*!
117 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_finalTrace
118
119 \brief The pappso::Trace holding the final results of the computations.
120 */
121
122 /*!
123 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_mzIntegrationParams
124
125 \brief The configuration of the mass spectral combinations (for
126 example, determines the bins, if binning is required).
127 */
128
129 /*!
130 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_mostIntensePeakMz
131
132 \brief The most intense peak encountered during the calculations.
133 */
134
135 /*!
136 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_smallestMz
137
138 \brief The smallest m/z value encountered during the calculations.
139
140 This value is required for the crafting of the bins.
141 */
142
143 /*!
144 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_greatestMz
145
146 \brief The greatest m/z value encountered during the calculations.
147
148 This value is required for the crafting of the bins.
149 */
150
151 /*!
152 \variable MsXpS::libXpertMassCore::IsotopicClusterShaper::m_normalizeIntensity
153
154 \brief The value by which all the peak shapes need to be normalized.
155 */
156
157
158 /*!
159 \brief Constructs a IsotopicClusterShaper instance.
160
161 \list
162 \li \a isotopic_cluster: The peak centroids belonging to an isotopic cluster
163 in the form of a pappso::Trace.
164 \li \a charge: The charge of the analyte below the isotopic cluster peaks.
165 \li \a config: The mass peak shaping process configuration.
166 \endlist
167
168 Each pappso::DataPoint in the \a isotopic_cluster pappso::Trace corresponds to
169 a (m/z, intensity) peak centroid belonging to the isotopic cluster.
170 */
171 IsotopicClusterShaper::IsotopicClusterShaper(
172 const pappso::Trace &isotopic_cluster,
173 int charge,
174 const MassPeakShaperConfig &config)
175 : mp_massPeakShaperConfig(new MassPeakShaperConfig(config))
176 {
177 // No need to reset, we are constructing.
178 setIsotopicCluster(isotopic_cluster, charge, false);
179 }
180
181 /*!
182 \brief Constructs a IsotopicClusterShaper instance.
183
184 \list
185 \li \a isotopic_cluster_charge_pairs: The pairs associating a pappso::Trace
186 instance to its corresponding charge.
187 \li \a config: The mass peak shaping process configuration.
188 \endlist
189
190 In this constructor, a set of sets of (m/z, charge) peak centroids is passed as
191 argument. The pappso::Trace instances in \a isotopic_cluster_charge_pairs might
192 correspond to all the isotopic clusters representing a given analyte at
193 different charges.
194 */
195 IsotopicClusterShaper::IsotopicClusterShaper(
196 const std::vector<IsotopicClusterChargePair> &isotopic_cluster_charge_pairs,
197 const MassPeakShaperConfig &config)
198 : mp_massPeakShaperConfig(new MassPeakShaperConfig(config))
199 {
200 // No need to reset, we are constructing.
201 setIsotopicClusterChargePairs(isotopic_cluster_charge_pairs, false);
202 }
203
204 /*!
205 \brief Destructs this IsotopicClusterShaper instance.
206 */
207 IsotopicClusterShaper::~IsotopicClusterShaper()
208 {
209 }
210
211 /*!
212 \brief Sets the peak shaping process \a{config}uration.
213 */
214 void
215 IsotopicClusterShaper::setConfig(const MassPeakShaperConfig &config)
216 {
217 mp_massPeakShaperConfig->initialize(config);
218 }
219
220 /*!
221 \brief Gets the peak shaping process configuration.
222 */
223 const MassPeakShaperConfig *
224 IsotopicClusterShaper::getConfig() const
225 {
226 return mp_massPeakShaperConfig;
227 }
228
229 /*!
230 \brief Sets the intensity normalization value to \a max_intensity.
231 */
232 void
233 IsotopicClusterShaper::setNormalizeIntensity(int max_intensity)
234 {
235 m_normalizeIntensity = max_intensity;
236 }
237
238 /*!
239 \brief Gets the intensity normalization value.
240 */
241 int
242 IsotopicClusterShaper::getNormalizeIntensity() const
243 {
244 return m_normalizeIntensity;
245 }
246
247 /*!
248 \brief Runs the mass peak shaping process for all the \l
249 IsotopicClusterChargePair objects in \l m_isotopicClusterChargePairs.
250
251 If \a reset is true, the member \l m_mapTrace and \l m_finalTrace are cleared
252 before starting the computations. The \l m_config member is first
253 \l{MassPeakShaperConfig::resolve()}d to check that all the parameters have been
254 properly set and are valid.
255
256 Returns the obtained pappso::Trace corresponding to the combination of all the
257 individual traces obtained for the various isotopic clusters and their
258 corresponding charge.
259 */
260 pappso::Trace &
261 IsotopicClusterShaper::run(bool reset)
262 {
263 if(reset)
264 {
265 // Clear the map trace that will receive the results of the combinations.
266 m_mapTrace.clear();
267 m_finalTrace.clear();
268 }
269
270 ErrorList error_list;
271
272 if(!mp_massPeakShaperConfig->resolve(error_list))
273 {
274 qDebug() << "Failed to resolve the MassPeakShaperConfig with errors:\n"
275 << Utils::joinErrorList(error_list, "\n");
276 return m_finalTrace;
277 }
278
279 // This class works on a vector of pairs containing the following:
280 // 1. a pappso::TraceCstSPtr
281 // 2. a charge
282
283 // We will process each pair in turn. If the integration requires bin, then
284 // each shaped isotopic cluster will be combined into a mass spectrum.
285 // Otherwise a trace combiner will be used.
286
287 // When setting the data (either by construction or using the set<> functions,
288 // we had monitored the smallest and the greatest m/z value over the whole set
289 // of the DataPoint objects in the isotopic clusters (centroid data). This is
290 // because we need, in case binning is required, these values to craft the
291 // bins.
292
293 // We will need to perform combinations, positive combinations.
294 // This mass spectrum combiner is in case we need binning.
295 pappso::MassSpectrumPlusCombiner mass_spectrum_plus_combiner;
296
297 // This trace combiner is in case do *not* need binning.
298 pappso::TracePlusCombiner trace_plus_combiner(-1);
299
300 // Configure the mass spectrum combiner in case we need binning.
301
302 if(mp_massPeakShaperConfig->isWithBins())
303 {
304 // Bins were requested.
305
306 // qDebug() << "Bins are required.";
307
308 // Get the bin size out of the configuration.
309
310 double bin_size = mp_massPeakShaperConfig->getBinSize();
311
312 // qDebug() << "The bin size in the config is:" << bin_size;
313
314 // Because we had monitored the m/z values of all the shapes generated
315 // above, we know the smallest and greatest m/z value that were
316 // encountered in all that peak shaping activity. We thus can create the
317 // bins faithfully.
318
319 Enums::MassPeakWidthLogic logic =
320 mp_massPeakShaperConfig->getMassPeakWidthLogic();
321
322 // The m_smallestMz and m_greatestMz values have been determined by
323 // looking into the unshaped isotopic clusters passed to this object
324 // either by construction or with functions. These two mz values are thus
325 // peak centroids, not data points belonging to a shaped peak since we
326 // have not yet started shaping the peaks. This means that we cannot
327 // create bins start / ending at these values because we would loose the
328 // first half of the first shaped peak centroid and the second half of the
329 // last shaped peak centroid (a shaped peak goes left *and* right of the
330 // peak centroid otherwise there would be no shape).
331 //
332 // This is why we provide a confortable margin fo the bin creation below
333 // by removing 1 Th on the left of the smallest mz and by adding 1 Th on
334 // right of the greatest mz.
335
336 if(logic == Enums::MassPeakWidthLogic::FWHM)
337 {
338 m_mzIntegrationParams.initialize(
339 m_smallestMz - 1,
340 m_greatestMz + 1,
341 pappso::MzIntegrationParams::BinningType::ARBITRARY,
342 pappso::PrecisionFactory::getDaltonInstance(bin_size),
343 /*binSizeDivisor*/ 1,
344 /*decimalPlacesr*/ -1,
345 true,
346 nullptr);
347 }
348 else if(logic == Enums::MassPeakWidthLogic::RESOLUTION)
349 {
350 double resolution = mp_massPeakShaperConfig->getResolution();
351
352 m_mzIntegrationParams.initialize(
353 m_smallestMz - 1,
354 m_greatestMz + 1,
355 pappso::MzIntegrationParams::BinningType::ARBITRARY,
356 pappso::PrecisionFactory::getResInstance(resolution),
357 mp_massPeakShaperConfig->getBinSizeDivisor(),
358 -1,
359 true,
360 nullptr);
361 }
362 else
363 qFatal(
364 "Programming error. By this time, the mass peak width logic should "
365 "have been defined.");
366
367 // qDebug() << "The mz integration params:"
368 //<< m_mzIntegrationParams.toString();
369
370 // Now compute the bins.
371
372 std::vector<double> bins = m_mzIntegrationParams.createBins();
373 // qDebug() << "The bins:" << bins;
374
375 mass_spectrum_plus_combiner.setBins(bins);
376 // qDebug() << "Set bins to the mass spectrum combiner:"
377 //<< mass_spectrum_plus_combiner.getBins().size();
378 }
379
380 std::size_t peak_centroid_count = 0;
381 std::size_t isotopic_cluster_count = 0;
382
383 // Iterate in the isotopic cluster/charge pairs.
384 for(auto isotopic_cluster_charge_pair : m_isotopicClusterChargePairs)
385 {
386 int charge = isotopic_cluster_charge_pair.second;
387
388 // Iterate in the data points of the current centroid data isotopic
389 // cluster.
390 for(auto data_point : *isotopic_cluster_charge_pair.first)
391 {
392 // Note the division by m_charge below!
393
394 pappso::Trace trace = MassPeakShaper::computePeakShape(
395 data_point.x / charge, data_point.y, *mp_massPeakShaperConfig);
396
397 // qDebug() << "The shaped isotopic cluster has size:" <<
398 // trace.size();
399
400 if(trace.size())
401 {
402 if(mp_massPeakShaperConfig->isWithBins())
403 mass_spectrum_plus_combiner.combine(m_mapTrace, trace);
404 else
405 trace_plus_combiner.combine(m_mapTrace, trace);
406
407 // qDebug() << qSetRealNumberPrecision(15)
408 //<< "The map trace for combination has size:"
409 //<< m_mapTrace.size()
410 //<< "and starting m/z:" << m_mapTrace.begin()->first
411 //<< "and ending m/z:"
412 //<< std::prev(m_mapTrace.end())->first;
413
414 ++peak_centroid_count;
415 }
416 }
417 ++isotopic_cluster_count;
418 }
419
420 // qDebug() << QString(
421 //"Successfully processed %1 isotopic clusters for a total of %2 "
422 //"shaped peak centroids")
423 //.arg(isotopic_cluster_count)
424 //.arg(peak_centroid_count);
425
426 // The user might have asked that the most intense m/z peak centroid be used
427 // for normalization. In that case that peak centroid's intensity is brought
428 // to m_normalizeIntensity and the ratio between its current intensity and
429 // m_normalizeIntensity is used to normalize all the other data points in the
430 // trace.
431
432 if(m_normalizeIntensity != 1)
433 {
434
435 // qDebug() << "Now normalizing to intensity = " << m_normalizeIntensity;
436
437 pappso::Trace trace = m_mapTrace.toTrace();
438 m_finalTrace =
439 trace.filter(pappso::FilterNormalizeIntensities(m_normalizeIntensity));
440
441 // double max_int = normalized_trace.maxYDataPoint().y;
442 // qDebug() << "After normalization max int:" << max_int;
443 }
444 else
445 m_finalTrace = m_mapTrace.toTrace();
446
447 // qDebug() << "Returning a trace of size:" << m_finalTrace.size();
448
449 // pappso::Utils::writeToFile(m_finalTrace.toString(), "/tmp/mass/trace.txt");
450
451 return m_finalTrace;
452 }
453
454 /*!
455 \brief Handles the \a isotopic_cluster_sp input isotopic cluster as a
456 pappso::Trace.
457
458 \a isotopic_cluster_sp is associated to a \a charge. If \a reset is true, the
459 member vector of \l IsotopicClusterChargePair instances is cleared before the
460 computations.
461
462 This function is the workhorse for all the functions used to set the initial
463 data for the computations. Its main task is to scrutinize the data in \a
464 isotopic_cluster_sp and update the \l m_smallestMz, \l m_greatestMz and \l
465 m_mostIntensePeakMz values based on the data passed as argument.
466 */
467 void
468 IsotopicClusterShaper::setIsotopicCluster(
469 pappso::TraceCstSPtr isotopic_cluster_sp, int charge, bool reset)
470 {
471 if(reset)
472 m_isotopicClusterChargePairs.clear();
473
474 double min_x = isotopic_cluster_sp->minX();
475 m_smallestMz = std::min(m_smallestMz, min_x);
476
477 double max_x = isotopic_cluster_sp->maxX();
478 m_greatestMz = std::max(m_greatestMz, max_x);
479
480 m_mostIntensePeakMz = isotopic_cluster_sp->maxYDataPoint().x;
481
482 // qDebug() << qSetRealNumberPrecision(15) << "m_smallestMz:" << m_smallestMz
483 //<< "m_greatestMz:" << m_greatestMz
484 //<< "m_mostIntensePeakMz:" << m_mostIntensePeakMz;
485
486 m_isotopicClusterChargePairs.push_back(
487 IsotopicClusterChargePair(isotopic_cluster_sp, charge));
488
489 mp_massPeakShaperConfig->setReferencePeakMz(m_mostIntensePeakMz);
490 }
491
492 /*!
493 \brief Handles the \a isotopic_cluster_sp input isotopic cluster as a
494 pappso::Trace.
495
496 \a isotopic_cluster_sp is associated to a \a charge.
497
498 The member vector of \l IsotopicClusterChargePair instances is cleared before
499 the computations.
500 */
501 void
502 IsotopicClusterShaper::setIsotopicCluster(
503 pappso::TraceCstSPtr isotopic_cluster_sp, int charge)
504 {
505 setIsotopicCluster(isotopic_cluster_sp, charge, true);
506 }
507
508 /*!
509 \brief Handles the \a isotopic_cluster input isotopic cluster as a
510 pappso::Trace.
511
512 \a isotopic_cluster is associated to a \a charge. If \a reset is true, the
513 member vector of \l IsotopicClusterChargePair instances is cleared before the
514 computations.
515 */
516 void
517 IsotopicClusterShaper::setIsotopicCluster(const pappso::Trace &isotopic_cluster,
518 int charge,
519 bool reset)
520 {
521 setIsotopicCluster(
522 std::make_shared<const pappso::Trace>(isotopic_cluster), charge, reset);
523 }
524
525 /*!
526 \brief Handles the \a isotopic_cluster input isotopic cluster as a
527 pappso::Trace.
528
529 \a isotopic_cluster is associated to a \a charge.
530
531 The member vector of \l IsotopicClusterChargePair instances is cleared before
532 the computations.
533 */
534 void
535 IsotopicClusterShaper::setIsotopicCluster(const pappso::Trace &isotopic_cluster,
536 int charge)
537 {
538 setIsotopicCluster(
539 std::make_shared<const pappso::Trace>(isotopic_cluster), charge, true);
540 }
541
542 /*!
543 \brief Handles the \a isotopic_cluster_charge_pair input isotopic cluster as a
544 IsotopicClusterChargePair.
545
546 The member vector of \l IsotopicClusterChargePair instances is cleared before
547 the computations.
548 */
549 void
550 IsotopicClusterShaper::setIsotopicCluster(
551 IsotopicClusterChargePair isotopic_cluster_charge_pair)
552 {
553 setIsotopicCluster(isotopic_cluster_charge_pair.first,
554 isotopic_cluster_charge_pair.second,
555 true);
556 }
557
558 /*!
559 \brief Handles the \a isotopic_cluster_charge_pair input isotopic cluster as a
560 IsotopicClusterChargePair.
561
562 If \a reset is true, the member vector of \l IsotopicClusterChargePair instances
563 is cleared before the computations.
564 */
565 void
566 IsotopicClusterShaper::setIsotopicCluster(
567 IsotopicClusterChargePair isotopic_cluster_charge_pair, bool reset)
568 {
569 setIsotopicCluster(isotopic_cluster_charge_pair.first,
570 isotopic_cluster_charge_pair.second,
571 reset);
572 }
573
574 /*!
575 \brief Handles the \a isotopic_cluster_charge_pairs input isotopic cluster as a
576 vector of IsotopicClusterChargePair instances.
577
578 If \a reset is true, the member vector of \l IsotopicClusterChargePair
579 instances is cleared before the computations.
580 */
581 void
582 IsotopicClusterShaper::setIsotopicClusterChargePairs(
583 const std::vector<IsotopicClusterChargePair> &isotopic_cluster_charge_pairs,
584 bool reset)
585 {
586 for(auto cluster_charge_pair : isotopic_cluster_charge_pairs)
587 setIsotopicCluster(
588 cluster_charge_pair.first, cluster_charge_pair.second, reset);
589
590 // qDebug() << qSetRealNumberPrecision(15) << "m_smallestMz:" << m_smallestMz
591 //<< "m_greatestMz:" << m_greatestMz
592 //<< "m_mostIntensePeakMz:" << m_mostIntensePeakMz;
593 }
594
595 /*!
596 \brief Handles the \a isotopic_cluster_charge_pairs input isotopic cluster as a
597 vector of IsotopicClusterChargePair instances.
598
599 The member vector of \l IsotopicClusterChargePair instances is cleared before
600 the computations.
601 */
602 void
603 IsotopicClusterShaper::setIsotopicClusterChargePairs(
604 const std::vector<IsotopicClusterChargePair> &isotopic_cluster_charge_pairs)
605 {
606 setIsotopicClusterChargePairs(isotopic_cluster_charge_pairs, true);
607
608 // qDebug() << qSetRealNumberPrecision(15) << "m_smallestMz:" << m_smallestMz
609 //<< "m_greatestMz:" << m_greatestMz
610 //<< "m_mostIntensePeakMz:" << m_mostIntensePeakMz;
611 }
612
613 /*!
614 \brief Adds an isotopic cluster in the form of the \a isotopic_cluster
615 pappso::Trace associated to the corresponding \a charge.
616
617 The member vector of \l IsotopicClusterChargePair instances is \e{not} cleared
618 before the computations.
619 */
620 void
621 IsotopicClusterShaper::appendIsotopicCluster(
622 const pappso::Trace &isotopic_cluster, int charge)
623 {
624 // Do not clear the isotopic clusters!
625
626 setIsotopicCluster(isotopic_cluster, charge, false);
627 }
628
629 /*!
630 \brief Adds isotopic clusters in the form of the \a
631 isotopic_cluster_charge_pairs vector of IsotopicClusterChargePair instances.
632
633 The member vector of \l IsotopicClusterChargePair instances is \e{not} cleared
634 before the computations.
635 */
636 void
637 IsotopicClusterShaper::appendIsotopicClusterChargePairs(
638 const std::vector<IsotopicClusterChargePair> &isotopic_cluster_charge_pairs)
639 {
640 setIsotopicClusterChargePairs(isotopic_cluster_charge_pairs, false);
641
642 // qDebug() << "m_smallestMz:" << m_smallestMz << "m_greatestMz:" <<
643 // m_greatestMz
644 //<< "m_mostIntensePeakMz:" << m_mostIntensePeakMz;
645 }
646
647 /*!
648 \brief Returns the final result of all the computations as a string.
649 */
650 QString
651 IsotopicClusterShaper::shapeToString()
652 {
653 return m_finalTrace.toString();
654 }
655
656 } // namespace libXpertMassCore
657
658 } // namespace MsXpS
659