filterresample.cpp

Go to the documentation of this file.

/**
 * \file pappsomspp/filers/filterresample.cpp
 * \date 28/04/2019
 * \author Olivier Langella
 * \brief collection of filters concerned by X selection
 */
 
/*******************************************************************************
 * Copyright (c) 2019 Olivier Langella <Olivier.Langella@u-psud.fr>.
 *
 * This file is part of the PAPPSOms++ library.
 *
 *     PAPPSOms++ is free software: you can redistribute it and/or modify
 *     it under the terms of the GNU General Public License as published by
 *     the Free Software Foundation, either version 3 of the License, or
 *     (at your option) any later version.
 *
 *     PAPPSOms++ is distributed in the hope that it will be useful,
 *     but WITHOUT ANY WARRANTY; without even the implied warranty of
 *     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *     GNU General Public License for more details.
 *
 *     You should have received a copy of the GNU General Public License
 *     along with PAPPSOms++.  If not, see <http://www.gnu.org/licenses/>.
 *
 ******************************************************************************/
 
#include "filterresample.h"
#include "pappsomspp/core/massspectrum/massspectrum.h"
#include <QDebug>
 
namespace pappso
{
 
 
FilterResampleKeepSmaller::FilterResampleKeepSmaller(double x_value) : m_value(x_value)
{
}
 
FilterResampleKeepSmaller::FilterResampleKeepSmaller(const FilterResampleKeepSmaller &other)
  : FilterResampleKeepSmaller(other.m_value)
{
}
 
 
Trace &
FilterResampleKeepSmaller::filter(Trace &spectrum) const
{
  auto begin_it = findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_value);
  spectrum.erase(begin_it, spectrum.end());
  return spectrum;
}
 
FilterResampleKeepGreater::FilterResampleKeepGreater(double x_value) : m_value(x_value)
{
}
 
FilterResampleKeepGreater::FilterResampleKeepGreater(const FilterResampleKeepGreater &other)
  : FilterResampleKeepGreater(other.m_value)
{
}
 
 
double
FilterResampleKeepGreater::getThresholdX() const
{
  return m_value;
}
 
FilterResampleKeepGreater &
FilterResampleKeepGreater::operator=(const FilterResampleKeepGreater &other)
{
  m_value = other.m_value;
 
  return *this;
}
 
Trace &
FilterResampleKeepGreater::filter(Trace &spectrum) const
{
  // qDebug() << " spectrum.size()=" << spectrum.size();
 
  auto last_it = findFirstGreaterX(spectrum.begin(), spectrum.end(), m_value);
  spectrum.erase(spectrum.begin(), last_it);
 
  // qDebug() <<  " spectrum.size()=" << spectrum.size();
 
  return spectrum;
}
 
FilterResampleRemoveXRange::FilterResampleRemoveXRange(double min_x, double max_x)
  : m_minX(min_x), m_maxX(max_x)
{
}
 
FilterResampleRemoveXRange::FilterResampleRemoveXRange(const FilterResampleRemoveXRange &other)
  : FilterResampleRemoveXRange(other.m_minX, other.m_maxX)
{
}
 
 
FilterResampleRemoveXRange &
FilterResampleRemoveXRange::operator=(const FilterResampleRemoveXRange &other)
{
  m_minX = other.m_minX;
  m_maxX = other.m_maxX;
 
  return *this;
}
 
 
Trace &
FilterResampleRemoveXRange::filter(Trace &spectrum) const
{
 
  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
  //          << " m_min_x=" << m_min_x;
  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
  //          << " m_max_x=" << m_max_x;
  auto begin_it = findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_minX);
  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
  //          << " begin_it->x=" << begin_it->x;
  auto end_it = findFirstGreaterX(begin_it, spectrum.end(), m_maxX);
  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
  //         << " end_it->x=" << end_it->x;
  spectrum.erase(begin_it, end_it);
 
  // qDebug() << __FILE__ << " " << __FUNCTION__ << " " << __LINE__
  //          << " spectrum.size()=" << spectrum.size();
  return spectrum;
}
 
 
FilterResampleKeepXRange::FilterResampleKeepXRange(double min_x, double max_x)
  : m_minX(min_x), m_maxX(max_x)
{
}
 
FilterResampleKeepXRange::FilterResampleKeepXRange(const FilterResampleKeepXRange &other)
  : m_minX(other.m_minX), m_maxX(other.m_maxX)
{
}
 
 
FilterResampleKeepXRange &
FilterResampleKeepXRange::operator=(const FilterResampleKeepXRange &other)
{
  if(&other == this)
    return *this;
 
  m_minX = other.m_minX;
  m_maxX = other.m_maxX;
 
  return *this;
}
 
 
Trace &
FilterResampleKeepXRange::filter(Trace &spectrum) const
{
  // qDebug() << "The range to keep:" << m_minX << "-" << m_maxX;
 
  auto begin_it = findFirstEqualOrGreaterX(spectrum.begin(), spectrum.end(), m_minX);
 
  // qDebug() << "Found begin iterator (for m_minX) having:" << begin_it->x
  //<< "x (m/z) value";
 
  auto end_it = findFirstGreaterX(begin_it, spectrum.end(), m_maxX);
 
  if(end_it == spectrum.end())
    {
      // qDebug() << "The end iterator (for m_maxX) is the end(). The prev "
      //"iterator has"
      //<< std::prev(end_it)->x << " x(m / z) value.";
    }
  else
    {
      // qDebug() << "Found end iterator (for m_maxX) having:" << end_it->x
      //<< "x (m/z) value";
    }
 
  // qDebug() << "Only keeping range" << begin_it->x << "-"
  //<< std::prev(end_it)->x;
 
  spectrum.erase(end_it, spectrum.end());
  spectrum.erase(spectrum.begin(), begin_it);
 
  return spectrum;
}
 
 
MassSpectrumFilterResampleRemoveMzRange::MassSpectrumFilterResampleRemoveMzRange(
  const MzRange &mz_range)
  : m_filterRange(mz_range.lower(), mz_range.upper())
{
}
 
MassSpectrumFilterResampleRemoveMzRange::MassSpectrumFilterResampleRemoveMzRange(
  const MassSpectrumFilterResampleRemoveMzRange &other)
  : m_filterRange(other.m_filterRange)
{
}
 
MassSpectrum &
MassSpectrumFilterResampleRemoveMzRange::filter(MassSpectrum &spectrum) const
{
  // qDebug() << m_filterRange.filter(spectrum);
  m_filterRange.filter(spectrum);
  return spectrum;
}
 
 
MassSpectrumFilterResampleKeepMzRange::MassSpectrumFilterResampleKeepMzRange(
  const MzRange &mz_range)
  : m_filterRange(mz_range.lower(), mz_range.upper())
{
}
 
MassSpectrumFilterResampleKeepMzRange::MassSpectrumFilterResampleKeepMzRange(
  const MassSpectrumFilterResampleKeepMzRange &other)
  : m_filterRange(other.m_filterRange)
{
}
 
MassSpectrum &
MassSpectrumFilterResampleKeepMzRange::filter(MassSpectrum &spectrum) const
{
  m_filterRange.filter(spectrum);
  return spectrum;
}
 
 
FilterResampleKeepPointInPolygon::FilterResampleKeepPointInPolygon()
{
}
 
 
FilterResampleKeepPointInPolygon::FilterResampleKeepPointInPolygon(
  IntegrationScopeBaseSPtr integration_scope_sp, Enums::DataKind data_kind)
{
  m_integrationScopeSpecs.push_back(IntegrationScopeSpec(integration_scope_sp, data_kind));
 
  QPointF point;
  if(!m_integrationScopeSpecs.front().integrationScopeSPtr->getPoint(point))
    qFatal("Failed to get point.");
  m_lowestMz = point.y();
 
  double height;
  if(!m_integrationScopeSpecs.front().integrationScopeSPtr->getHeight(height))
    qFatal("Failed to get height.");
  m_greatestMz = m_lowestMz + height;
}
 
 
FilterResampleKeepPointInPolygon::FilterResampleKeepPointInPolygon(
  const IntegrationScopeSpecVector &integration_scope_specs)
{
  // qDebug();
 
  m_integrationScopeSpecs.assign(integration_scope_specs.begin(), integration_scope_specs.end());
 
  for(auto &&item : m_integrationScopeSpecs)
    {
      if(item.integrationScopeSPtr->isRectangle())
        {
          // qDebug() << "item is Rectangle.";
 
          QPointF point;
          if(!item.integrationScopeSPtr->getPoint(point))
            qFatal("Failed to get point.");
          m_lowestMz = std::min(m_lowestMz, point.y());
 
          double height;
          if(!item.integrationScopeSPtr->getHeight(height))
            qFatal("Failed to get height.");
          m_greatestMz = std::max(m_greatestMz, m_lowestMz + height);
        }
      else
        {
          if(!item.integrationScopeSPtr->isRhomboid())
            qFatal(
              "If integration scope is not Rectangle, then it must be "
              "Rhomboid.");
 
          // qDebug() << "item is Rhomboid.";
 
          QPointF point;
 
          // Because the mz scale is on the coordinate, we ask for
          // bottom-most and top-most points of the rhomboid polygon.
 
          if(!item.integrationScopeSPtr->getBottomMostPoint(point))
            qFatal("Failed to get point.");
          m_lowestMz = std::min(m_lowestMz, point.y());
 
          if(!item.integrationScopeSPtr->getTopMostPoint(point))
            qFatal("Failed to get point.");
          m_greatestMz = std::max(m_greatestMz, point.y());
        }
    }
}
 
 
FilterResampleKeepPointInPolygon::FilterResampleKeepPointInPolygon(
  const FilterResampleKeepPointInPolygon &other)
{
  // qDebug();
 
  m_integrationScopeSpecs.assign(other.m_integrationScopeSpecs.begin(),
                                 other.m_integrationScopeSpecs.end());
 
  for(auto &&item : m_integrationScopeSpecs)
    {
      QPointF point;
      if(!item.integrationScopeSPtr->getPoint(point))
        qFatal("Failed to get point.");
      m_lowestMz = std::min(m_lowestMz, point.y());
 
      double height;
      if(!item.integrationScopeSPtr->getHeight(height))
        qFatal("Failed to get height.");
      m_greatestMz = std::max(m_greatestMz, m_lowestMz + height);
    }
}
 
 
void
FilterResampleKeepPointInPolygon::newSelectionPolygonSpec(
  const IntegrationScopeSpec &integration_scope_spec)
{
 
  m_integrationScopeSpecs.push_back(integration_scope_spec);
 
  QPointF point;
  if(!m_integrationScopeSpecs.back().integrationScopeSPtr->getPoint(point))
    qFatal("Failed to get point.");
  m_lowestMz = std::min(m_lowestMz, point.y());
 
  double height;
  if(!m_integrationScopeSpecs.back().integrationScopeSPtr->getHeight(height))
    qFatal("Failed to get height.");
  m_greatestMz = std::max(m_greatestMz, m_lowestMz + height);
}
 
 
FilterResampleKeepPointInPolygon &
FilterResampleKeepPointInPolygon::operator=(const FilterResampleKeepPointInPolygon &other)
{
  if(this == &other)
    return *this;
 
  m_integrationScopeSpecs.assign(other.m_integrationScopeSpecs.begin(),
                                 other.m_integrationScopeSpecs.end());
 
  m_lowestMz   = other.m_lowestMz;
  m_greatestMz = other.m_greatestMz;
 
  return *this;
}
 
 
Trace &
FilterResampleKeepPointInPolygon::filter([[maybe_unused]] Trace &trace) const
{
  qFatal("Programming error.");
  return trace;
}
 
 
Trace &
FilterResampleKeepPointInPolygon::filter(Trace &trace, double dt_value, double rt_value) const
{
  // Each time a new integration scope spec is added, the lowest and greatest
  // m/z values are computed. We use these values to remove from the spectrum
  // all the points that are outside of that lowest-gratest range.
 
  // Find the iterator to the most front of the DataPoint vector (mass
  // spectrum).
 
  // Note that the m_lowestMz and m_greatestMz are set during construction of
  // this FilterResampleKeepPointInPolygon filter using the
  // integration scope specs.
 
  // qDebug() << "The lowest and greatest m/z values:" << m_lowestMz << "and"
  //<< m_greatestMz;
 
  // Start by filtering away all the data points outside of the
  // [m_lowestMz--m_greatestMz] range.
 
  FilterResampleKeepXRange the_filter(m_lowestMz, m_greatestMz);
 
  trace = the_filter.filter(trace);
 
  // Now iterate in all the data points remaining in the trace and for each
  // point craft a "virtual" point using the dt|rt value and the m/z value of
  // the data point (data_point.x).
 
  auto begin_it = trace.begin();
  auto end_it   = trace.end();
 
  // qDebug() << "Iterating in the m/z range:" << begin_it->x << "-"
  //<< std::prev(end_it)->x;
 
  // Start at the end of the range. The iter is outside of the requested range,
  // in fact, as shown using iter-- below.
  auto iter = end_it;
 
  while(iter > begin_it)
    {
      // Immediately go to the last data point of the desired iteration range of
      // the trace that we need to filter.  Remember that end() is not pointing
      // to any vector item, but past the last one.
      iter--;
 
      // qDebug() << "Iterating in trace data point with m/z value:" << iter->x;
 
      // Now that we have the m/z value, we can check it, in combination with
      // the value in the selection polygon spec (to make a point) against the
      // various selection polygon specs in our member vector.
 
      double checked_value;
 
      for(auto &&spec : m_integrationScopeSpecs)
        {
          // qDebug() << "Iterating in selection polygon spec:" <<
          // spec.toString();
 
          if(spec.dataKind == Enums::DataKind::dt)
            {
              if(dt_value == -1)
                qFatal("Programming error.");
 
              checked_value = dt_value;
 
              // qDebug() << "The data kind: dt.";
            }
          else if(spec.dataKind == Enums::DataKind::rt)
            {
              checked_value = rt_value;
 
              // qDebug() << "The data kind: rt.";
            }
          else
            qFatal("Programming error.");
 
          // First version doing the whole computation on the basis of the
          // selection polygon's all faces.
          //
          if(!spec.integrationScopeSPtr->contains(QPointF(checked_value, iter->x)))
            iter = trace.erase(iter);
        }
      // End of
      // for(auto &&spec : m_selectionPolygonSpecs)
    }
  // End of
  // while(iter > begin_it)
 
  return trace;
}
 
} // namespace pappso