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() << " spectrum.size()=" << spectrum.size();
87 
88  auto last_it = findFirstGreaterX(spectrum.begin(), spectrum.end(), m_value);
89  spectrum.erase(spectrum.begin(), last_it);
90 
91  // qDebug() << " spectrum.size()=" << spectrum.size();
92 
93  return spectrum;
94 }
95 
97  double max_x)
98  : m_minX(min_x), m_maxX(max_x)
99 {
100 }
101 
103  const FilterResampleRemoveXRange &other)
104  : FilterResampleRemoveXRange(other.m_minX, other.m_maxX)
105 {
106 }
107 
108 
111 {
112  m_minX = other.m_minX;
113  m_maxX = other.m_maxX;
114 
115  return *this;
116 }
117 
118 
119 Trace &
121 {
122 
123  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
124  // << " m_min_x=" << m_min_x;
125  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
126  // << " m_max_x=" << m_max_x;
127  auto begin_it =
128  findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_minX);
129  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
130  // << " begin_it->x=" << begin_it->x;
131  auto end_it = findFirstGreaterX(begin_it, spectrum.end(), m_maxX);
132  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
133  // << " end_it->x=" << end_it->x;
134  spectrum.erase(begin_it, end_it);
135 
136  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
137  // << " spectrum.size()=" << spectrum.size();
138  return spectrum;
139 }
140 
141 
143  : m_minX(min_x), m_maxX(max_x)
144 {
145 }
146 
148  const FilterResampleKeepXRange &other)
149  : m_minX(other.m_minX), m_maxX(other.m_maxX)
150 {
151 }
152 
153 
156 {
157  if(&other == this)
158  return *this;
159 
160  m_minX = other.m_minX;
161  m_maxX = other.m_maxX;
162 
163  return *this;
164 }
165 
166 
167 Trace &
169 {
170  // qDebug() << "The range to keep:" << m_minX << "-" << m_maxX;
171 
172  auto begin_it =
173  findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_minX);
174 
175  // qDebug() << "Found begin iterator (for m_minX) having:" << begin_it->x
176  //<< "x (m/z) value";
177 
178  auto end_it = findFirstGreaterX(begin_it, spectrum.end(), m_maxX);
179 
180  if(end_it == spectrum.end())
181  {
182  // qDebug() << "The end iterator (for m_maxX) is the end(). The prev "
183  //"iterator has"
184  //<< std::prev(end_it)->x << " x(m / z) value.";
185  }
186  else
187  {
188  // qDebug() << "Found end iterator (for m_maxX) having:" << end_it->x
189  //<< "x (m/z) value";
190  }
191 
192  // qDebug() << "Only keeping range" << begin_it->x << "-"
193  //<< std::prev(end_it)->x;
194 
195  spectrum.erase(end_it, spectrum.end());
196  spectrum.erase(spectrum.begin(), begin_it);
197 
198  return spectrum;
199 }
200 
201 
204  : m_filterRange(mz_range.lower(), mz_range.upper())
205 {
206 }
207 
211  : m_filterRange(other.m_filterRange)
212 {
213 }
214 
215 MassSpectrum &
217 {
218  // qDebug() << m_filterRange.filter(spectrum);
219  m_filterRange.filter(spectrum);
220  return spectrum;
221 }
222 
223 
225  const MzRange &mz_range)
226  : m_filterRange(mz_range.lower(), mz_range.upper())
227 {
228 }
229 
232  : m_filterRange(other.m_filterRange)
233 {
234 }
235 
236 MassSpectrum &
238 {
239  m_filterRange.filter(spectrum);
240  return spectrum;
241 }
242 
243 
245 {
246 }
247 
248 
250  const SelectionPolygon &selection_polygon, DataKind data_kind)
251 {
252  // It is assumed that the selection polygon always has x:MZ and y:DT|RT
253  // depending on the spec data kind.
254 
255  m_selectionPolygonSpecs.push_back(
256  SelectionPolygonSpec(selection_polygon, data_kind));
257 
258  m_lowestMz =
259  m_selectionPolygonSpecs.front().selectionPolygon.getBottomMostPoint().y();
260  m_greatestMz =
261  m_selectionPolygonSpecs.front().selectionPolygon.getTopMostPoint().y();
262 }
263 
264 
266  const SelectionPolygonSpecVector &selection_polygon_specs)
267 {
268  // qDebug();
269 
270  m_selectionPolygonSpecs.assign(selection_polygon_specs.begin(),
271  selection_polygon_specs.end());
272 
273  for(auto &&item : m_selectionPolygonSpecs)
274  {
275  m_lowestMz =
276  std::min(m_lowestMz, item.selectionPolygon.getBottomMostPoint().y());
277 
278  m_greatestMz =
279  std::max(m_greatestMz, item.selectionPolygon.getTopMostPoint().y());
280  }
281 }
282 
283 
286 {
287  // qDebug();
288 
290  other.m_selectionPolygonSpecs.end());
291 
292  for(auto &&item : m_selectionPolygonSpecs)
293  {
294  m_lowestMz =
295  std::min(m_lowestMz, item.selectionPolygon.getBottomMostPoint().y());
296 
297  m_greatestMz =
298  std::max(m_greatestMz, item.selectionPolygon.getTopMostPoint().y());
299  }
300 }
301 
302 
303 void
305  const SelectionPolygonSpec &selection_polygon_spec)
306 {
307  // It is assumed that the selection polygon always has x:MZ and y:DT|RT
308  // depending on the spec data kind.
309 
310  m_selectionPolygonSpecs.push_back(selection_polygon_spec);
311 
312  m_lowestMz = std::min(
313  m_lowestMz,
314  m_selectionPolygonSpecs.back().selectionPolygon.getBottomMostPoint().y());
315 
316  m_greatestMz = std::max(
317  m_greatestMz,
318  m_selectionPolygonSpecs.back().selectionPolygon.getTopMostPoint().y());
319 }
320 
321 
325 {
326  if(this == &other)
327  return *this;
328 
330  other.m_selectionPolygonSpecs.end());
331 
332  m_lowestMz = other.m_lowestMz;
333  m_greatestMz = other.m_greatestMz;
334 
335  return *this;
336 }
337 
338 
339 Trace &
340 FilterResampleKeepPointInPolygon::filter([[maybe_unused]] Trace &trace) const
341 {
342  qFatal("Programming error.");
343  return trace;
344 }
345 
346 
347 Trace &
349  double dt_value,
350  double rt_value) const
351 {
352  // Each time a new selection polygon spec is added, the lowest and greatest
353  // m/z values are computed. We use these values to remove from the spectrum
354  // all the points that are outside of that lowest-gratest range.
355 
356  // Find the iterator to the most front of the DataPoint vector (mass
357  // spectrum).
358 
359  // Note that the m_lowestMz and m_greatestMz are set during construction of
360  // this FilterResampleKeepPointInPolygon filter using the
361  // selection polygon specs.
362 
363  // qDebug() << "The lowest and greatest m/z values:" << m_lowestMz << "and"
364  //<< m_greatestMz;
365 
366  // Start by filtering away all the data points outside of the
367  // [m_lowestMz--m_greatestMz] range.
368 
370 
371  trace = the_filter.filter(trace);
372 
373  // Now iterate in all the data points remaining in the trace and for each
374  // point craft a "virtual" point using the dt|rt value and the m/z value of
375  // the data point (data_point.x).
376 
377  auto begin_it = trace.begin();
378  auto end_it = trace.end();
379 
380  // qDebug() << "Iterating in the m/z range:" << begin_it->x << "-"
381  //<< std::prev(end_it)->x;
382 
383  // Start at the end of the range. The iter is outside of the requested range,
384  // in fact, as shown using iter-- below.
385  auto iter = end_it;
386 
387  while(iter > begin_it)
388  {
389  // Immediately go to the last data point of the desired iteration range of
390  // the trace that we need to filter. Remember that end() is not pointing
391  // to any vector item, but past the last one.
392  iter--;
393 
394  // qDebug() << "Iterating in trace data point with m/z value:" << iter->x;
395 
396  // Now that we have the m/z value, we can check it, in combination with
397  // the value in the selection polygon spec (to make a point) against the
398  // various selection polygon specs in our member vector.
399 
400  double checked_value;
401 
402  for(auto &&spec : m_selectionPolygonSpecs)
403  {
404  // qDebug() << "Iterating in selection polygon spec:" <<
405  // spec.toString();
406 
407  if(spec.dataKind == DataKind::dt)
408  {
409  if(dt_value == -1)
410  qFatal("Programming error.");
411 
412  checked_value = dt_value;
413 
414  // qDebug() << "The data kind: dt.";
415  }
416  else if(spec.dataKind == DataKind::rt)
417  {
418  checked_value = rt_value;
419 
420  // qDebug() << "The data kind: rt.";
421  }
422  else
423  qFatal("Programming error.");
424 
425  // First version doing the whole computation on the basis of the
426  // selection polygon's all faces.
427  //
428  if(!spec.selectionPolygon.contains(QPointF(checked_value, iter->x)))
429  iter = trace.erase(iter);
430 
431 #if 0
432 
433  //This code does not work because depending on the orientation of the
434  //skewed selection polygon (left bottom to right top or left top to
435  //right bottom or or Or depending on the axes of the
436  //bi-dimensional colormap, requiring transposition or not), the
437  //notion of "left line" and of "right line" changes.
438 
439  // Second version checking that point is right of left vertical line
440  // of polygon and left of right vertical line.
441 
442  // double res = sideofline(XX;YY;xA;yA;xB;yB) =
443  // (xB-xA) * (YY-yA) - (yB-yA) * (XX-xA)
444 
445  // If res == 0, the point is on the line
446  // If rest < 0, the point is on the right of the line
447  // If rest > 0, the point is on the left of the line
448 
449  // Left vertical line of the polygon
450 
451  // Bottom point
452  double xA_left =
453  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_LEFT_POINT).x();
454  double yA_left =
455  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_LEFT_POINT).y();
456 
457  // Top point
458  double xB_left =
459  spec.selectionPolygon.getPoint(PointSpecs::TOP_LEFT_POINT).x();
460  double yB_left =
461  spec.selectionPolygon.getPoint(PointSpecs::TOP_LEFT_POINT).y();
462 
463  qDebug() << "The left line goes: (" << xA_left << "," << yA_left
464  << ")->(" << xB_left << "," << yB_left << ")";
465 
466  if((xB_left - xA_left) * (iter->x - yA_left) -
467  (yB_left - yA_left) * (checked_value - xA_left) >
468  0)
469  {
470  // The point is left of the left line. We can remove the point
471  // from the mass spectrum.
472 
473  qDebug() << qSetRealNumberPrecision(10)
474  << "Filtered out point (left of left line):"
475  << checked_value << "-" << iter->x;
476 
477  iter = trace.erase(iter);
478 
479  // No need to go on with the analysis, just go to the next (that
480  // is, previous, since we iterate backwards) data point.
481 
482  continue;
483  }
484  else
485  {
486  qDebug() << qSetRealNumberPrecision(10)
487  << "Kept point (right of left line):" << checked_value
488  << "-" << iter->x;
489  }
490 
491  // Right vertical line of the polygon
492 
493  // Bottom point
494  double xA_right =
495  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_RIGHT_POINT).x();
496  double yA_right =
497  spec.selectionPolygon.getPoint(PointSpecs::BOTTOM_RIGHT_POINT).y();
498 
499  // Top point
500  double xB_right =
501  spec.selectionPolygon.getPoint(PointSpecs::TOP_RIGHT_POINT).x();
502  double yB_right =
503  spec.selectionPolygon.getPoint(PointSpecs::TOP_RIGHT_POINT).y();
504 
505  qDebug() << "The right line goes: (" << xA_right << "," << yA_right
506  << ")->(" << xB_right << "," << yB_right << ")";
507 
508  if((xB_right - xA_right) * (iter->x - yA_right) -
509  (yB_right - yA_right) * (checked_value - xA_right) <
510  0)
511  {
512  qDebug() << qSetRealNumberPrecision(10)
513  << "Filtered out point (right of right line):"
514  << checked_value << "-" << iter->x;
515 
516  // The point is right of the right line. We can remove the point
517  // from the mass spectrum.
518  iter = trace.erase(iter);
519 
520  // Here, continue is implicit.
521  // No need to go on with the analysis, just go to the next (that
522  // is, previous, since we iterate backwards) data point.
523  // continue;
524  }
525  else
526  {
527  if(iter->x >= 449 && iter->x <= 450 && checked_value > 19.5 &&
528  checked_value < 20)
529  qDebug()
530  << qSetRealNumberPrecision(10)
531  << "SHOULD NOT Definitively kept point (left of right line):"
532  << checked_value << "-" << iter->x;
533  else
534  qDebug()
535  << qSetRealNumberPrecision(10)
536  << "MIGHT Definitively kept point (left of right line):"
537  << checked_value << "-" << iter->x;
538  }
539 #endif
540  }
541  // End of
542  // for(auto &&spec : m_selectionPolygonSpecs)
543  }
544  // End of
545  // while(iter > begin_it)
546 
547  return trace;
548 }
549 
550 } // namespace pappso
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:148
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:69
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:97
DataKind
Definition: types.h:172
@ dt
Drift time.
@ rt
Retention time.
std::vector< SelectionPolygonSpec > SelectionPolygonSpecVector