libpappsomspp
Library for mass spectrometry
filterresample.cpp
Go to the documentation of this file.
1 /**
2  * \file pappsomspp/filers/filterresample.cpp
3  * \date 28/04/2019
4  * \author Olivier Langella
5  * \brief collection of filters concerned by X selection
6  */
7 
8 /*******************************************************************************
9  * Copyright (c) 2019 Olivier Langella <Olivier.Langella@u-psud.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  ******************************************************************************/
27 
28 #include "filterresample.h"
29 #include "../../massspectrum/massspectrum.h"
30 #include <QDebug>
31 
32 namespace pappso
33 {
34 
35 
37  : m_value(x_value)
38 {
39 }
40 
42  const FilterResampleKeepSmaller &other)
43  : FilterResampleKeepSmaller(other.m_value)
44 {
45 }
46 
47 
48 Trace &
50 {
51  auto begin_it =
52  findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_value);
53  spectrum.erase(begin_it, spectrum.end());
54  return spectrum;
55 }
56 
58  : m_value(x_value)
59 {
60 }
61 
63  const FilterResampleKeepGreater &other)
64  : FilterResampleKeepGreater(other.m_value)
65 {
66 }
67 
68 
69 double
71 {
72  return m_value;
73 }
74 
77 {
78  m_value = other.m_value;
79 
80  return *this;
81 }
82 
83 Trace &
85 {
86  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
87  << " spectrum.size()=" << spectrum.size();
88  auto last_it = findFirstGreaterX(spectrum.begin(), spectrum.end(), m_value);
89  spectrum.erase(spectrum.begin(), last_it);
90  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
91  << " spectrum.size()=" << spectrum.size();
92  return spectrum;
93 }
94 
96  double max_x)
97  : m_minX(min_x), m_maxX(max_x)
98 {
99 }
100 
102  const FilterResampleRemoveXRange &other)
103  : FilterResampleRemoveXRange(other.m_minX, other.m_maxX)
104 {
105 }
106 
107 
110 {
111  m_minX = other.m_minX;
112  m_maxX = other.m_maxX;
113 
114  return *this;
115 }
116 
117 
118 Trace &
120 {
121 
122  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
123  // << " m_min_x=" << m_min_x;
124  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
125  // << " m_max_x=" << m_max_x;
126  auto begin_it =
127  findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_minX);
128  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
129  // << " begin_it->x=" << begin_it->x;
130  auto end_it = findFirstGreaterX(begin_it, spectrum.end(), m_maxX);
131  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
132  // << " end_it->x=" << end_it->x;
133  spectrum.erase(begin_it, end_it);
134 
135  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
136  // << " spectrum.size()=" << spectrum.size();
137  return spectrum;
138 }
139 
140 
142  : m_minX(min_x), m_maxX(max_x)
143 {
144 }
145 
147  const FilterResampleKeepXRange &other)
148  : m_minX(other.m_minX), m_maxX(other.m_maxX)
149 {
150 }
151 
152 
155 {
156  if(&other == this)
157  return *this;
158 
159  m_minX = other.m_minX;
160  m_maxX = other.m_maxX;
161 
162  return *this;
163 }
164 
165 
166 Trace &
168 {
169  qDebug() << "The range to keep:" << m_minX << "-" << m_maxX;
170 
171  auto begin_it =
172  findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_minX);
173 
174  qDebug() << "Found begin iterator (for m_minX) having:" << begin_it->x
175  << "x (m/z) value";
176 
177  auto end_it = findFirstGreaterX(begin_it, spectrum.end(), m_maxX);
178 
179  if(end_it == spectrum.end())
180  {
181  qDebug() << "The end iterator (for m_maxX) is the end(). The prev "
182  "iterator has"
183  << std::prev(end_it)->x << " x(m / z) value.";
184  }
185  else
186  qDebug() << "Found end iterator (for m_maxX) having:" << end_it->x
187  << "x (m/z) value";
188 
189  qDebug() << "Only keeping range" << begin_it->x << "-"
190  << std::prev(end_it)->x;
191 
192  spectrum.erase(end_it, spectrum.end());
193  spectrum.erase(spectrum.begin(), begin_it);
194 
195  return spectrum;
196 }
197 
198 
201  : m_filterRange(mz_range.lower(), mz_range.upper())
202 {
203 }
204 
208  : m_filterRange(other.m_filterRange)
209 {
210 }
211 
212 MassSpectrum &
214 {
215  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__;
216  m_filterRange.filter(spectrum);
217  qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__;
218  return spectrum;
219 }
220 
221 
223  const MzRange &mz_range)
224  : m_filterRange(mz_range.lower(), mz_range.upper())
225 {
226 }
227 
230  : m_filterRange(other.m_filterRange)
231 {
232 }
233 
234 MassSpectrum &
236 {
237  m_filterRange.filter(spectrum);
238  return spectrum;
239 }
240 
241 
243 {
244 }
245 
246 
248  const SelectionPolygon &selection_polygon, DataKind data_kind)
249 {
250  // It is assumed that the selection polygon always has x:MZ and y:DT|RT
251  // depending on the spec data kind.
252 
253  m_selectionPolygonSpecs.push_back(
254  SelectionPolygonSpec(selection_polygon, data_kind));
255 
256  m_lowestMz =
257  m_selectionPolygonSpecs.front().selectionPolygon.getBottomMostPoint().y();
258  m_greatestMz =
259  m_selectionPolygonSpecs.front().selectionPolygon.getTopMostPoint().y();
260 }
261 
262 
264  const SelectionPolygonSpecVector &selection_polygon_specs)
265 {
266  qDebug();
267 
268  m_selectionPolygonSpecs.assign(selection_polygon_specs.begin(),
269  selection_polygon_specs.end());
270 
271  for(auto &&item : m_selectionPolygonSpecs)
272  {
273  m_lowestMz =
274  std::min(m_lowestMz, item.selectionPolygon.getBottomMostPoint().y());
275 
276  m_greatestMz =
277  std::max(m_greatestMz, item.selectionPolygon.getTopMostPoint().y());
278  }
279 }
280 
281 
284 {
285  qDebug();
286 
288  other.m_selectionPolygonSpecs.end());
289 
290  for(auto &&item : m_selectionPolygonSpecs)
291  {
292  m_lowestMz =
293  std::min(m_lowestMz, item.selectionPolygon.getBottomMostPoint().y());
294 
295  m_greatestMz =
296  std::max(m_greatestMz, item.selectionPolygon.getTopMostPoint().y());
297  }
298 }
299 
300 
301 void
303  const SelectionPolygonSpec &selection_polygon_spec)
304 {
305  // It is assumed that the selection polygon always has x:MZ and y:DT|RT
306  // depending on the spec data kind.
307 
308  m_selectionPolygonSpecs.push_back(selection_polygon_spec);
309 
310  m_lowestMz = std::min(
311  m_lowestMz,
312  m_selectionPolygonSpecs.back().selectionPolygon.getBottomMostPoint().y());
313 
314  m_greatestMz = std::max(
315  m_greatestMz,
316  m_selectionPolygonSpecs.back().selectionPolygon.getTopMostPoint().y());
317 }
318 
319 
323 {
324  if(this == &other)
325  return *this;
326 
328  other.m_selectionPolygonSpecs.end());
329 
330  m_lowestMz = other.m_lowestMz;
331  m_greatestMz = other.m_greatestMz;
332 
333  return *this;
334 }
335 
336 
337 Trace &
338 FilterResampleKeepPointInPolygon::filter([[maybe_unused]] Trace &trace) const
339 {
340  qFatal("Programming error.");
341  return trace;
342 }
343 
344 
345 Trace &
347  double dt_value,
348  double rt_value) const
349 {
350  // Each time a new selection polygon spec is added, the lowest and greatest
351  // m/z values are computed. We use these values to remove from the spectrum
352  // all the points that are outside of that lowest-gratest range.
353 
354  // Find the iterator to the most front of the DataPoint vector (mass
355  // spectrum).
356 
357  // Note that the m_lowestMz and m_greatestMz are set during construction of
358  // this FilterResampleKeepPointInPolygon filter using the
359  // selection polygon specs.
360 
361  qDebug() << "The lowest and greatest m/z values:" << m_lowestMz << "and"
362  << m_greatestMz;
363 
364  // Start by filtering away all the data points outside of the
365  // [m_lowestMz--m_greatestMz] range.
366 
368 
369  trace = the_filter.filter(trace);
370 
371  // Now iterate in all the data points remaining in the trace and for each
372  // point craft a "virtual" point using the dt|rt value and the m/z value of
373  // the data point (data_point.x).
374 
375  auto begin_it = trace.begin();
376  auto end_it = trace.end();
377 
378  qDebug() << "Iterating in the m/z range:" << begin_it->x << "-"
379  << std::prev(end_it)->x;
380 
381  // Start at the end of the range. The iter is outside of the requested range,
382  // in fact, as shown using iter-- below.
383  auto iter = end_it;
384 
385  while(iter > begin_it)
386  {
387  // Immediately go to the last data point of the desired iteration range of
388  // the trace that we need to filter. Remember that end() is not pointing
389  // to any vector item, but past the last one.
390  iter--;
391 
392  qDebug() << "Iterating in trace data point with m/z value:" << iter->x;
393 
394  // Now that we have the m/z value, we can check it, in combination with
395  // the value in the selection polygon spec (to make a point) against the
396  // various selection polygon specs in our member vector.
397 
398  double checked_value;
399 
400  for(auto &&spec : m_selectionPolygonSpecs)
401  {
402  // qDebug() << "Iterating in selection polygon spec:" <<
403  // spec.toString();
404 
405  if(spec.dataKind == DataKind::dt)
406  {
407  if(dt_value == -1)
408  qFatal("Programming error.");
409 
410  checked_value = dt_value;
411 
412  // qDebug() << "The data kind: dt.";
413  }
414  else if(spec.dataKind == DataKind::rt)
415  {
416  checked_value = rt_value;
417 
418  // qDebug() << "The data kind: rt.";
419  }
420  else
421  qFatal("Programming error.");
422 
423  // First version doing the whole computation on the basis of the
424  // selection polygon's all faces.
425  //
426  if(!spec.selectionPolygon.contains(QPointF(checked_value, iter->x)))
427  iter = trace.erase(iter);
428 
429 #if 0
430 
431  //This code does not work because depending on the orientation of the
432  //skewed selection polygon (left bottom to right top or left top to
433  //right bottom or or Or depending on the axes of the
434  //bi-dimensional colormap, requiring transposition or not), the
435  //notion of "left line" and of "right line" changes.
436 
437  // Second version checking that point is right of left vertical line
438  // of polygon and left of right vertical line.
439 
440  // double res = sideofline(XX;YY;xA;yA;xB;yB) =
441  // (xB-xA) * (YY-yA) - (yB-yA) * (XX-xA)
442 
443  // If res == 0, the point is on the line
444  // If rest < 0, the point is on the right of the line
445  // If rest > 0, the point is on the left of the line
446 
447  // Left vertical line of the polygon
448 
449  // Bottom point
450  double xA_left =
451  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_LEFT_POINT).x();
452  double yA_left =
453  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_LEFT_POINT).y();
454 
455  // Top point
456  double xB_left =
457  spec.selectionPolygon.getPoint(PointSpecs::TOP_LEFT_POINT).x();
458  double yB_left =
459  spec.selectionPolygon.getPoint(PointSpecs::TOP_LEFT_POINT).y();
460 
461  qDebug() << "The left line goes: (" << xA_left << "," << yA_left
462  << ")->(" << xB_left << "," << yB_left << ")";
463 
464  if((xB_left - xA_left) * (iter->x - yA_left) -
465  (yB_left - yA_left) * (checked_value - xA_left) >
466  0)
467  {
468  // The point is left of the left line. We can remove the point
469  // from the mass spectrum.
470 
471  qDebug() << qSetRealNumberPrecision(10)
472  << "Filtered out point (left of left line):"
473  << checked_value << "-" << iter->x;
474 
475  iter = trace.erase(iter);
476 
477  // No need to go on with the analysis, just go to the next (that
478  // is, previous, since we iterate backwards) data point.
479 
480  continue;
481  }
482  else
483  {
484  qDebug() << qSetRealNumberPrecision(10)
485  << "Kept point (right of left line):" << checked_value
486  << "-" << iter->x;
487  }
488 
489  // Right vertical line of the polygon
490 
491  // Bottom point
492  double xA_right =
493  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_RIGHT_POINT).x();
494  double yA_right =
495  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_RIGHT_POINT).y();
496 
497  // Top point
498  double xB_right =
499  spec.selectionPolygon.getPoint(PointSpecs::TOP_RIGHT_POINT).x();
500  double yB_right =
501  spec.selectionPolygon.getPoint(PointSpecs::TOP_RIGHT_POINT).y();
502 
503  qDebug() << "The right line goes: (" << xA_right << "," << yA_right
504  << ")->(" << xB_right << "," << yB_right << ")";
505 
506  if((xB_right - xA_right) * (iter->x - yA_right) -
507  (yB_right - yA_right) * (checked_value - xA_right) <
508  0)
509  {
510  qDebug() << qSetRealNumberPrecision(10)
511  << "Filtered out point (right of right line):"
512  << checked_value << "-" << iter->x;
513 
514  // The point is right of the right line. We can remove the point
515  // from the mass spectrum.
516  iter = trace.erase(iter);
517 
518  // Here, continue is implicit.
519  // No need to go on with the analysis, just go to the next (that
520  // is, previous, since we iterate backwards) data point.
521  // continue;
522  }
523  else
524  {
525  if(iter->x >= 449 && iter->x <= 450 && checked_value > 19.5 &&
526  checked_value < 20)
527  qDebug()
528  << qSetRealNumberPrecision(10)
529  << "SHOULD NOT Definitively kept point (left of right line):"
530  << checked_value << "-" << iter->x;
531  else
532  qDebug()
533  << qSetRealNumberPrecision(10)
534  << "MIGHT Definitively kept point (left of right line):"
535  << checked_value << "-" << iter->x;
536  }
537 #endif
538  }
539  // End of
540  // for(auto &&spec : m_selectionPolygonSpecs)
541  }
542  // End of
543  // while(iter > begin_it)
544 
545  return trace;
546 }
547 
548 } // namespace pappso
549 
FilterResampleKeepGreater & operator=(const FilterResampleKeepGreater &other)
Trace & filter(Trace &trace) const override
void newSelectionPolygonSpec(const SelectionPolygonSpec &selection_polygon_spec)
std::vector< SelectionPolygonSpec > m_selectionPolygonSpecs
virtual Trace & filter(Trace &data_points) const=0
FilterResampleKeepPointInPolygon & operator=(const FilterResampleKeepPointInPolygon &other)
Trace & filter(Trace &trace) const override
FilterResampleKeepXRange & operator=(const FilterResampleKeepXRange &other)
FilterResampleKeepXRange(double min_x=0, double max_x=0)
Trace & filter(Trace &trace) const override
FilterResampleRemoveXRange(double min_x, double max_x)
FilterResampleRemoveXRange & operator=(const FilterResampleRemoveXRange &other)
Trace & filter(Trace &trace) const override
MassSpectrumFilterResampleKeepMzRange(const MzRange &mz_range)
MassSpectrum & filter(MassSpectrum &spectrum) const override
const FilterResampleKeepXRange m_filterRange
const FilterResampleRemoveXRange m_filterRange
MassSpectrumFilterResampleRemoveMzRange(const MzRange &mz_range)
MassSpectrum & filter(MassSpectrum &spectrum) const override
Class to represent a mass spectrum.
Definition: massspectrum.h:71
A simple container of DataPoint instances.
Definition: trace.h:132
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition: aa.cpp:39
std::vector< DataPoint >::iterator findFirstEqualOrGreaterX(std::vector< DataPoint >::iterator begin, std::vector< DataPoint >::iterator end, const double &value)
find the first element in which X is equal or greater than the value searched important : it implies ...
Definition: trace.cpp:32
std::vector< DataPoint >::iterator findFirstGreaterX(std::vector< DataPoint >::iterator begin, std::vector< DataPoint >::iterator end, const double &value)
find the first element in which X is greater than the value searched important : it implies that Trac...
Definition: trace.cpp:60
DataKind
Definition: types.h:171
@ dt
Drift time.
@ rt
Retention time.
std::vector< SelectionPolygonSpec > SelectionPolygonSpecVector