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>> &data_point_vector)
62  : Trace::Trace(data_point_vector)
63 {
64 }
65 
66 MassSpectrum::MassSpectrum(std::vector<DataPoint> &data_point_vector)
67  : Trace::Trace(data_point_vector)
68 {
69 }
70 
71 MassSpectrum::MassSpectrum(const Trace &other) : Trace(other)
72 {
73 }
74 
76 {
77 }
78 
79 
80 MassSpectrum::MassSpectrum(Trace &&other) : Trace(std::move(other))
81 {
82 }
83 
84 
86 {
87  // qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << "()";
88 }
89 
90 
91 MassSpectrum::MassSpectrum(MassSpectrum &&other) : Trace(std::move(other))
92 {
93  // Specify std::move so that && reference is passed to the Trace constructor
94  // that takes std::vector<DataPoint> && as rvalue reference.
95 
96  // qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << "()"
97  //<< "Moving MassSpectrum::MassSpectrum(MassSpectrum &&)";
98 }
99 
100 
102 {
103 }
104 
105 
106 MassSpectrum &
108 {
109  // qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << " ()";
110 
111  assign(other.begin(), other.end());
112  return *this;
113 }
114 
115 
116 MassSpectrum &
118 {
119  vector<DataPoint>::operator=(std::move(other));
120  return *this;
121 }
122 
123 
126 {
127  return std::make_shared<MassSpectrum>(*this);
128 }
129 
130 
133 {
134  return std::make_shared<const MassSpectrum>(*this);
135 }
136 
137 
138 //! Compute the total ion current of this mass spectrum
139 /*!
140  * The sum of all the separate ion currents carried by the ions of different
141  * m/z contributing to a complete mass massSpectrum or in a specified m/z
142  * range of a mass massSpectrum. MS:1000285
143  *
144  * \return <pappso_double> The total ion current.
145  */
148 {
149  return Trace::sumY();
150 }
151 
152 
153 //! Compute the total ion current of this mass spectrum
154 /*!
155  * Convenience function that returns totalIonCurrent();
156  */
159 {
160  return totalIonCurrent();
161 }
162 
163 
165 MassSpectrum::tic(double mzStart, double mzEnd)
166 {
167  return Trace::sumY(mzStart, mzEnd);
168 }
169 
170 
171 //! Find the DataPoint instance having the greatest intensity (y) value.
172 /*!
173  * \return <const DataPoint &> The data point having the maximum intensity (y)
174  * value of the whole mass spectrum.
175  */
176 const DataPoint &
178 {
179  return Trace::maxYDataPoint();
180 }
181 
182 
183 //! Find the DataPoint instance having the smallest intensity (y) value.
184 /*!
185  * \return <const DataPoint &> The data point having the minimum intensity (y)
186  * value of the whole mass spectrum.
187  */
188 const DataPoint &
190 {
191  return Trace::minYDataPoint();
192 }
193 
194 
195 //! Sort the DataPoint instances of this spectrum.
196 /*!
197  * The DataPoint instances are sorted according to the x value (the m/z
198  * value) and in increasing order.
199  */
200 void
202 {
203  Trace::sortX();
204 }
205 
206 
207 //! Tells if \c this MassSpectrum is equal to \p massSpectrum.
208 /*!
209  * To compare \c this to \p massSpectrum, a tolerance is applied to both the x
210  * and y values, that is defined using \p precision.
211  *
212  * \param massSpectrum Mass spectrum to compare to \c this.
213  *
214  * \param precision Precision to be used to perform the comparison of the x
215  * and y values of the data points in \c this and \massSpectrum mass spectra.
216  */
217 bool
218 MassSpectrum::equals(const MassSpectrum &other, PrecisionPtr precision) const
219 {
220  if(size() != other.size())
221  {
222  qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << "()"
223  << "The other mass spectrum size is not equal to *this size"
224  << "*this size:" << size() << "trace size:" << other.size();
225 
226  return false;
227  }
228 
230 
231  auto trace_it = other.begin();
232 
233  for(auto &&data_point : *this)
234  {
235  qDebug() << "first:" << data_point.x << "," << data_point.y
236  << " second:" << trace_it->x << "," << trace_it->y;
237 
238  if(!MzRange(data_point.x, precision).contains(trace_it->x))
239  {
240  qDebug() << "x:" << data_point.x << " != " << trace_it->x;
241  return false;
242  }
243 
244  if(!MzRange(data_point.y, precint).contains(trace_it->y))
245  {
246  qDebug() << "y:" << data_point.y << " != " << trace_it->y;
247  return false;
248  }
249 
250  trace_it++;
251  }
252 
253  return true;
254 }
255 
256 
258 MassSpectrum::filterSum(const MzRange &range) const
259 {
260  MassSpectrum massSpectrum;
261 
262  std::vector<DataPoint>::const_iterator it = begin();
263  std::vector<DataPoint>::const_iterator itEnd = end();
264 
265  std::vector<DataPoint>::const_reverse_iterator itRev = rbegin();
266  std::vector<DataPoint>::const_reverse_iterator itRevEnd = rend();
267 
268  pappso_double lower = range.lower();
269  pappso_double upper = range.upper();
270 
271  while((it != itEnd) && (it->x <= itRev->x) && (itRev != itRevEnd))
272  {
273  pappso_double sumX = it->x + itRev->x;
274 
275  if(sumX < lower)
276  {
277  it++;
278  }
279  else if(sumX > upper)
280  {
281  itRev++;
282  }
283  else
284  {
285  massSpectrum.push_back(*it);
286  massSpectrum.push_back(*itRev);
287 
288  std::vector<DataPoint>::const_reverse_iterator itRevIn = itRev;
289  itRevIn++;
290 
291  // FIXME Attention buggy code FR 20180626.
292  sumX = it->x + itRevIn->x;
293  while((sumX > lower) && (it->x <= itRevIn->x) &&
294  (itRevIn != itRevEnd))
295  {
296  sumX = it->x + itRevIn->x;
297  // trace.push_back(*it);
298  massSpectrum.push_back(*itRevIn);
299  itRevIn++;
300  }
301  it++;
302  }
303  }
304 
305  // Sort all the data points in increasing order by x
306  std::sort(massSpectrum.begin(),
307  massSpectrum.end(),
308  [](const DataPoint &a, const DataPoint &b) { return (a.x < b.x); });
309 
310  // Remove all the but the first element of a series of elements that are
311  // considered equal. Sort of deduplication.
312  std::vector<DataPoint>::iterator itEndFix =
313  std::unique(massSpectrum.begin(),
314  massSpectrum.end(),
315  [](const DataPoint &a, const DataPoint &b) {
316  // Return true if both elements should be considered equal.
317  return (a.x == b.x) && (a.y == b.y);
318  });
319 
320  massSpectrum.resize(std::distance(massSpectrum.begin(), itEndFix));
321 
322  return massSpectrum;
323 }
324 
325 
326 void
328 {
329 
330  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__ << size();
331  for(std::size_t i = 0; i < size(); i++)
332  {
333  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__;
334  qDebug() << "mz = " << this->operator[](i).x
335  << ", int = " << this->operator[](i).y;
336  }
337 
338  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__;
339 }
340 
341 
342 QDataStream &
343 operator<<(QDataStream &outstream, const MassSpectrum &massSpectrum)
344 {
345  quint32 vector_size = massSpectrum.size();
346  outstream << vector_size;
347  for(auto &&peak : massSpectrum)
348  {
349  outstream << peak;
350  }
351 
352  return outstream;
353 }
354 
355 
356 QDataStream &
357 operator>>(QDataStream &instream, MassSpectrum &massSpectrum)
358 {
359 
360  quint32 vector_size;
361  DataPoint peak;
362 
363  if(!instream.atEnd())
364  {
365  instream >> vector_size;
366  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
367  << " vector_size=" << vector_size;
368 
369  for(quint32 i = 0; i < vector_size; i++)
370  {
371 
372  if(instream.status() != QDataStream::Ok)
373  {
374  throw PappsoException(
375  QString("error in QDataStream unserialize operator>> of "
376  "massSpectrum :\nread datastream failed status=%1 "
377  "massSpectrum "
378  "i=%2 on size=%3")
379  .arg(instream.status())
380  .arg(i)
381  .arg(vector_size));
382  }
383  instream >> peak;
384  massSpectrum.push_back(peak);
385  }
386  if(instream.status() != QDataStream::Ok)
387  {
388  throw PappsoException(
389  QString(
390  "error in QDataStream unserialize operator>> of massSpectrum "
391  ":\nread datastream failed status=%1")
392  .arg(instream.status()));
393  }
394  }
395 
396  return instream;
397 }
398 
399 
400 MassSpectrum &
402 {
403  return filter.filter(*this);
404 }
405 
406 } // namespace pappso
generic interface to apply a filter on a MassSpectrum This is the same as FilterInterface,...
Class to represent a mass spectrum.
Definition: massspectrum.h:71
void sortMz()
Sort the DataPoint instances of this spectrum.
void debugPrintValues() const
MassSpectrumCstSPtr makeMassSpectrumCstSPtr() const
bool equals(const MassSpectrum &other, PrecisionPtr precision) const
Tells if this MassSpectrum is equal to massSpectrum.
MassSpectrumSPtr makeMassSpectrumSPtr() const
pappso_double tic() const
Compute the total ion current of this mass spectrum.
const DataPoint & minIntensityDataPoint() const
Find the DataPoint instance having the smallest intensity (y) value.
virtual MassSpectrum & massSpectrumFilter(const MassSpectrumFilterInterface &filter) final
apply a filter on this MassSpectrum
MassSpectrum filterSum(const MzRange &mass_range) const
virtual MassSpectrum & operator=(const MassSpectrum &other)
pappso_double totalIonCurrent() const
Compute the total ion current of this mass spectrum.
const DataPoint & maxIntensityDataPoint() const
Find the DataPoint instance having the greatest intensity (y) value.
pappso_double lower() const
Definition: mzrange.h:71
pappso_double upper() const
Definition: mzrange.h:77
bool contains(pappso_double) const
Definition: mzrange.cpp:120
static PrecisionPtr getPpmInstance(pappso_double value)
get a ppm precision pointer
Definition: precision.cpp:150
A simple container of DataPoint instances.
Definition: trace.h:147
pappso_double sumY() const
Definition: trace.cpp:886
const DataPoint & maxYDataPoint() const
Definition: trace.cpp:853
void sortX()
Definition: trace.cpp:936
virtual Trace & filter(const FilterInterface &filter) final
apply a filter on this trace
Definition: trace.cpp:981
const DataPoint & minYDataPoint() const
Definition: trace.cpp:834
basic mass spectrum
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition: aa.cpp:39
QDataStream & operator>>(QDataStream &instream, MassSpectrum &massSpectrum)
QDataStream & operator<<(QDataStream &outstream, const MassSpectrum &massSpectrum)
double pappso_double
A type definition for doubles.
Definition: types.h:48
std::shared_ptr< const MassSpectrum > MassSpectrumCstSPtr
Definition: massspectrum.h:55
std::shared_ptr< MassSpectrum > MassSpectrumSPtr
Definition: massspectrum.h:54