libpappsomspp
Library for mass spectrometry
massspectrum.cpp
Go to the documentation of this file.
1 /**
2  * \file pappsomspp/massspectrum/massspectrum.cpp
3  * \date 15/3/2015
4  * \author Olivier Langella
5  * \brief basic mass spectrum
6  */
7 
8 /*******************************************************************************
9  * Copyright (c) 2015 Olivier Langella <Olivier.Langella@moulon.inra.fr>.
10  *
11  * This file is part of the PAPPSOms++ library.
12  *
13  * PAPPSOms++ is free software: you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License as published by
15  * the Free Software Foundation, either version 3 of the License, or
16  * (at your option) any later version.
17  *
18  * PAPPSOms++ is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with PAPPSOms++. If not, see <http://www.gnu.org/licenses/>.
25  *
26  * Contributors:
27  * Olivier Langella <Olivier.Langella@moulon.inra.fr> - initial API and
28  *implementation
29  ******************************************************************************/
30 
31 #include <cmath>
32 #include <numeric>
33 #include <limits>
34 #include <vector>
35 #include <map>
36 
37 #include <QDebug>
38 
39 #include "../trace/datapoint.h"
40 #include "../trace/trace.h"
41 #include "massspectrum.h"
42 #include "../processing/combiners/massspectrumcombiner.h"
43 #include "../mzrange.h"
44 #include "../pappsoexception.h"
45 
46 #include "../peptide/peptidefragmentionlistbase.h"
47 #include "../exception/exceptionoutofrange.h"
48 #include "../processing/filters/filterresample.h"
49 
50 
51 namespace pappso
52 {
53 
54 
56 {
57 }
58 
59 
61  std::vector<std::pair<pappso_double, pappso_double>> &dataPoints)
62  : Trace::Trace(dataPoints)
63 {
64 }
65 
66 MassSpectrum::MassSpectrum(const Trace &other) : Trace(other)
67 {
68 }
69 
70 MassSpectrum::MassSpectrum(const MapTrace &other) : Trace(other)
71 {
72 }
73 
74 
75 MassSpectrum::MassSpectrum(Trace &&other) : Trace(std::move(other))
76 {
77 }
78 
79 
80 MassSpectrum::MassSpectrum(const MassSpectrum &other) : Trace(other)
81 {
82  // qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << "()";
83 }
84 
85 
86 MassSpectrum::MassSpectrum(MassSpectrum &&other) : Trace(std::move(other))
87 {
88  // Specify std::move so that && reference is passed to the Trace constructor
89  // that takes std::vector<DataPoint> && as rvalue reference.
90 
91  // qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << "()"
92  //<< "Moving MassSpectrum::MassSpectrum(MassSpectrum &&)";
93 }
94 
95 
97 {
98 }
99 
100 
103 {
104  // qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << " ()";
105 
106  assign(other.begin(), other.end());
107  return *this;
108 }
109 
110 
111 MassSpectrum &
113 {
114  vector<DataPoint>::operator=(std::move(other));
115  return *this;
116 }
117 
118 
121 {
122  return std::make_shared<MassSpectrum>(*this);
123 }
124 
125 
128 {
129  return std::make_shared<const MassSpectrum>(*this);
130 }
131 
132 
133 //! Compute the total ion current of this mass spectrum
134 /*!
135  * The sum of all the separate ion currents carried by the ions of different
136  * m/z contributing to a complete mass massSpectrum or in a specified m/z
137  * range of a mass massSpectrum. MS:1000285
138  *
139  * \return <pappso_double> The total ion current.
140  */
143 {
144  return Trace::sumY();
145 }
146 
147 
148 //! Compute the total ion current of this mass spectrum
149 /*!
150  * Convenience function that returns totalIonCurrent();
151  */
153 MassSpectrum::tic() const
154 {
155  return totalIonCurrent();
156 }
157 
158 
160 MassSpectrum::tic(double mzStart, double mzEnd)
161 {
162  return Trace::sumY(mzStart, mzEnd);
163 }
164 
165 
166 //! Find the DataPoint instance having the greatest intensity (y) value.
167 /*!
168  * \return <const DataPoint &> The data point having the maximum intensity (y)
169  * value of the whole mass spectrum.
170  */
171 const DataPoint &
173 {
175 }
176 
177 
178 //! Find the DataPoint instance having the smallest intensity (y) value.
179 /*!
180  * \return <const DataPoint &> The data point having the minimum intensity (y)
181  * value of the whole mass spectrum.
182  */
183 const DataPoint &
185 {
186  return Trace::minYDataPoint();
187 }
188 
189 
190 //! Sort the DataPoint instances of this spectrum.
191 /*!
192  * The DataPoint instances are sorted according to the x value (the m/z
193  * value) and in increasing order.
194  */
195 void
197 {
198  Trace::sortX();
199 }
200 
201 
202 //! Tells if \c this MassSpectrum is equal to \p massSpectrum.
203 /*!
204  * To compare \c this to \p massSpectrum, a tolerance is applied to both the x
205  * and y values, that is defined using \p precision.
206  *
207  * \param massSpectrum Mass spectrum to compare to \c this.
208  *
209  * \param precision Precision to be used to perform the comparison of the x
210  * and y values of the data points in \c this and \massSpectrum mass spectra.
211  */
212 bool
213 MassSpectrum::equals(const MassSpectrum &other, PrecisionPtr precision) const
214 {
215  if(size() != other.size())
216  {
217  qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << "()"
218  << "The other mass spectrum size is not equal to *this size"
219  << "*this size:" << size() << "trace size:" << other.size();
220 
221  return false;
222  }
223 
225 
226  auto trace_it = other.begin();
227 
228  for(auto &&data_point : *this)
229  {
230  qDebug() << "first:" << data_point.x << "," << data_point.y
231  << " second:" << trace_it->x << "," << trace_it->y;
232 
233  if(!MzRange(data_point.x, precision).contains(trace_it->x))
234  {
235  qDebug() << "x:" << data_point.x << " != " << trace_it->x;
236  return false;
237  }
238 
239  if(!MzRange(data_point.y, precint).contains(trace_it->y))
240  {
241  qDebug() << "y:" << data_point.y << " != " << trace_it->y;
242  return false;
243  }
244 
245  trace_it++;
246  }
247 
248  return true;
249 }
250 
251 
252 MassSpectrum
253 MassSpectrum::filterSum(const MzRange &range) const
254 {
255  MassSpectrum massSpectrum;
256 
257  std::vector<DataPoint>::const_iterator it = begin();
258  std::vector<DataPoint>::const_iterator itEnd = end();
259 
260  std::vector<DataPoint>::const_reverse_iterator itRev = rbegin();
261  std::vector<DataPoint>::const_reverse_iterator itRevEnd = rend();
262 
263  pappso_double lower = range.lower();
264  pappso_double upper = range.upper();
265 
266  while((it != itEnd) && (it->x <= itRev->x) && (itRev != itRevEnd))
267  {
268  pappso_double sumX = it->x + itRev->x;
269 
270  if(sumX < lower)
271  {
272  it++;
273  }
274  else if(sumX > upper)
275  {
276  itRev++;
277  }
278  else
279  {
280  massSpectrum.push_back(*it);
281  massSpectrum.push_back(*itRev);
282 
283  std::vector<DataPoint>::const_reverse_iterator itRevIn = itRev;
284  itRevIn++;
285 
286  // FIXME Attention buggy code FR 20180626.
287  sumX = it->x + itRevIn->x;
288  while((sumX > lower) && (it->x <= itRevIn->x) &&
289  (itRevIn != itRevEnd))
290  {
291  sumX = it->x + itRevIn->x;
292  // trace.push_back(*it);
293  massSpectrum.push_back(*itRevIn);
294  itRevIn++;
295  }
296  it++;
297  }
298  }
299 
300  // Sort all the data points in increasing order by x
301  std::sort(massSpectrum.begin(),
302  massSpectrum.end(),
303  [](const DataPoint &a, const DataPoint &b) { return (a.x < b.x); });
304 
305  // Remove all the but the first element of a series of elements that are
306  // considered equal. Sort of deduplication.
307  std::vector<DataPoint>::iterator itEndFix =
308  std::unique(massSpectrum.begin(),
309  massSpectrum.end(),
310  [](const DataPoint &a, const DataPoint &b) {
311  // Return true if both elements should be considered equal.
312  return (a.x == b.x) && (a.y == b.y);
313  });
314 
315  massSpectrum.resize(std::distance(massSpectrum.begin(), itEndFix));
316 
317  return massSpectrum;
318 }
319 
320 
321 void
323 {
324 
325  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__ << size();
326  for(std::size_t i = 0; i < size(); i++)
327  {
328  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__;
329  qDebug() << "mz = " << this->operator[](i).x
330  << ", int = " << this->operator[](i).y;
331  }
332 
333  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__;
334 }
335 
336 
337 QDataStream &
338 operator<<(QDataStream &outstream, const MassSpectrum &massSpectrum)
339 {
340  quint32 vector_size = massSpectrum.size();
341  outstream << vector_size;
342  for(auto &&peak : massSpectrum)
343  {
344  outstream << peak;
345  }
346 
347  return outstream;
348 }
349 
350 
351 QDataStream &
352 operator>>(QDataStream &instream, MassSpectrum &massSpectrum)
353 {
354 
355  quint32 vector_size;
356  DataPoint peak;
357 
358  if(!instream.atEnd())
359  {
360  instream >> vector_size;
361  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
362  << " vector_size=" << vector_size;
363 
364  for(quint32 i = 0; i < vector_size; i++)
365  {
366 
367  if(instream.status() != QDataStream::Ok)
368  {
369  throw PappsoException(
370  QString("error in QDataStream unserialize operator>> of "
371  "massSpectrum :\nread datastream failed status=%1 "
372  "massSpectrum "
373  "i=%2 on size=%3")
374  .arg(instream.status())
375  .arg(i)
376  .arg(vector_size));
377  }
378  instream >> peak;
379  massSpectrum.push_back(peak);
380  }
381  if(instream.status() != QDataStream::Ok)
382  {
383  throw PappsoException(
384  QString(
385  "error in QDataStream unserialize operator>> of massSpectrum "
386  ":\nread datastream failed status=%1")
387  .arg(instream.status()));
388  }
389  }
390 
391  return instream;
392 }
393 
394 
395 MassSpectrum &
396 MassSpectrum::massSpectrumFilter(const MassSpectrumFilterInterface &filter)
397 {
398  return filter.filter(*this);
399 }
400 
401 } // namespace pappso
pappso::Trace::maxYDataPoint
const DataPoint & maxYDataPoint() const
Definition: trace.cpp:655
pappso::MassSpectrum::makeMassSpectrumSPtr
MassSpectrumSPtr makeMassSpectrumSPtr() const
Definition: massspectrum.cpp:141
pappso::pappso_double
double pappso_double
A type definition for doubles.
Definition: types.h:69
pappso::MassSpectrumCstSPtr
std::shared_ptr< const MassSpectrum > MassSpectrumCstSPtr
Definition: massspectrum.h:76
pappso::PrecisionFactory::getPpmInstance
static PrecisionPtr getPpmInstance(pappso_double value)
Definition: precision.cpp:170
massspectrum.h
basic mass spectrum
pappso::MassSpectrum::totalIonCurrent
pappso_double totalIonCurrent() const
Compute the total ion current of this mass spectrum.
Definition: massspectrum.cpp:163
pappso::MassSpectrum::~MassSpectrum
virtual ~MassSpectrum()
Definition: massspectrum.cpp:117
pappso
Definition: aa.cpp:38
pappso::MassSpectrum
Class to represent a mass spectrum.
Definition: massspectrum.h:91
pappso::MassSpectrum::debugPrintValues
void debugPrintValues() const
Definition: massspectrum.cpp:343
pappso::MassSpectrum::operator=
virtual MassSpectrum & operator=(const MassSpectrum &other)
Definition: massspectrum.cpp:123
pappso::MzRange::contains
bool contains(pappso_double) const
Definition: mzrange.cpp:124
pappso::MassSpectrum::sortMz
void sortMz()
Sort the DataPoint instances of this spectrum.
Definition: massspectrum.cpp:217
pappso::DataPoint
Definition: datapoint.h:20
pappso::MassSpectrum::lowestIntensityDataPoint
const DataPoint & lowestIntensityDataPoint() const
Find the DataPoint instance having the smallest intensity (y) value.
Definition: massspectrum.cpp:205
pappso::MassSpectrum::MassSpectrum
MassSpectrum()
Definition: massspectrum.cpp:76
pappso::MassSpectrum::filterSum
MassSpectrum filterSum(const MzRange &mass_range) const
Definition: massspectrum.cpp:274
pappso::Trace::filter
virtual Trace & filter(const FilterInterface &filter) final
apply a filter on this trace
Definition: trace.cpp:803
pappso::Trace::sumY
pappso_double sumY() const
Definition: trace.cpp:688
pappso::MzRange
Definition: mzrange.h:66
pappso::Trace
A simple container of DataPoint instances.
Definition: trace.h:126
pappso::MassSpectrum::massSpectrumFilter
virtual MassSpectrum & massSpectrumFilter(const MassSpectrumFilterInterface &filter) final
apply a filter on this MassSpectrum
Definition: massspectrum.cpp:417
pappso::operator>>
QDataStream & operator>>(QDataStream &instream, MassSpectrum &massSpectrum)
Definition: massspectrum.cpp:373
pappso::operator<<
QDataStream & operator<<(QDataStream &outstream, const MassSpectrum &massSpectrum)
Definition: massspectrum.cpp:359
pappso::Trace::sortX
void sortX()
Definition: trace.cpp:737
pappso::MassSpectrum::tic
pappso_double tic() const
Compute the total ion current of this mass spectrum.
Definition: massspectrum.cpp:174
pappso::PrecisionPtr
const typedef PrecisionBase * PrecisionPtr
Definition: precision.h:143
pappso::MassSpectrum::maxIntensityDataPoint
const DataPoint & maxIntensityDataPoint() const
Find the DataPoint instance having the greatest intensity (y) value.
Definition: massspectrum.cpp:193
pappso::MassSpectrum::equals
bool equals(const MassSpectrum &other, PrecisionPtr precision) const
Tells if this MassSpectrum is equal to massSpectrum.
Definition: massspectrum.cpp:234
pappso::MassSpectrum::makeMassSpectrumCstSPtr
MassSpectrumCstSPtr makeMassSpectrumCstSPtr() const
Definition: massspectrum.cpp:148
pappso::Trace::minYDataPoint
const DataPoint & minYDataPoint() const
Definition: trace.cpp:636
pappso::MassSpectrumSPtr
std::shared_ptr< MassSpectrum > MassSpectrumSPtr
Definition: massspectrum.h:75
pappso::PappsoException
Definition: pappsoexception.h:62