massspectraceplotwidget.cpp

Go to the documentation of this file.

/* This code comes right from the msXpertSuite software project.
 *
 * msXpertSuite - mass spectrometry software suite
 * -----------------------------------------------
 * Copyright(C) 2009,...,2018 Filippo Rusconi
 *
 * http://www.msxpertsuite.org
 *
 * This program 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.
 *
 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>.
 *
 * END software license
 */
 
 
/////////////////////// StdLib includes
#include <vector>
 
 
/////////////////////// Qt includes
#include <QVector>
 
 
/////////////////////// Local includes
 
// For the proton mass
#include "../../core/types.h"
 
#include "massspectraceplotwidget.h"
#include "massspectraceplotcontext.h"
 
 
int massSpecTracePlotContextMetaTypeId =
  qRegisterMetaType<pappso::MassSpecTracePlotContext>("pappso::MassSpecTracePlotContext");
 
int massSpecTracePlotContextPtrMetaTypeId =
  qRegisterMetaType<pappso::MassSpecTracePlotContext *>("pappso::MassSpecTracePlotContext *");
 
 
namespace pappso
{
 
MassSpecTracePlotWidget::MassSpecTracePlotWidget(QWidget *parent) : BaseTracePlotWidget(parent)
{
  BasePlotWidget::m_context.m_dataKind = Enums::DataKind::mz;
  m_context.m_dataKind                 = Enums::DataKind::mz;
 
  // qDebug() << "Data kind:" << static_cast<int>(m_context.m_dataKind);
}
 
MassSpecTracePlotWidget::MassSpecTracePlotWidget(QWidget *parent,
                                                 const QString &x_axis_label,
                                                 const QString &y_axis_label)
  : BaseTracePlotWidget(parent, x_axis_label, y_axis_label)
{
  // Set the base context to be of kind Enums::DataKind::mz;
 
  BasePlotWidget::m_context.m_dataKind = Enums::DataKind::mz;
  m_context.m_dataKind                 = Enums::DataKind::mz;
 
  // qDebug() << "Data kind:" << static_cast<int>(m_context.m_dataKind);
}
 
 
MassSpecTracePlotWidget::~MassSpecTracePlotWidget()
{
}
 
 
//! Set the \c m_pressedKeyCode to the key code in \p event.
void
MassSpecTracePlotWidget::keyPressEvent(QKeyEvent *event)
{
  // qDebug() << "ENTER";
  BasePlotWidget::keyPressEvent(event);
 
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  // qDebug() << "Going to emit keyPressEventSignal(m_context);";
 
  emit keyPressEventSignal(m_context);
}
 
 
//! Handle specific key codes and trigger respective actions.
void
MassSpecTracePlotWidget::keyReleaseEvent(QKeyEvent *event)
{
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  BasePlotWidget::keyReleaseEvent(event);
}
 
 
//! Handle mouse movements, in particular record all the last visited points.
/*!
 
  This function is reponsible for storing at each time the last visited point
  in the graph. Here, point is intended as any x/y coordinate in the plot
  widget viewport, not a graph point.
 
  The stored values are then the basis for a large set of calculations
  throughout all the plot widget.
 
  \param pointer to QMouseEvent from which to retrieve the coordinates of the
  visited viewport points.
  */
void
MassSpecTracePlotWidget::mouseMoveHandler(QMouseEvent *event)
{
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  BasePlotWidget::mouseMoveHandler(event);
}
 
 
void
MassSpecTracePlotWidget::mouseMoveHandlerNotDraggingCursor()
{
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  BasePlotWidget::mouseMoveHandlerNotDraggingCursor();
 
  // qDebug() << "The current m_lastZ value:" << m_context.m_lastZ;
}
 
 
void
MassSpecTracePlotWidget::mouseMoveHandlerDraggingCursor()
{
  BasePlotWidget::mouseMoveHandlerDraggingCursor();
 
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  if(m_context.m_mouseButtonsAtMousePress & Qt::LeftButton)
    {
      if(!m_context.m_isMeasuringDistance)
        return;
 
      // qDebug() << "Mouse Buttons at Mouse press:" << m_context.m_mouseButtonsAtMousePress;
 
      // qDebug() << "The current m_lastZ value:" << m_context.m_lastZ;
 
      deconvolute();
      computeResolvingPower();
    }
}
 
 
//! Record the clicks of the mouse.
void
MassSpecTracePlotWidget::mousePressHandler([[maybe_unused]] QMouseEvent *event)
{
  // qDebug() << "Entering";
 
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  BasePlotWidget::mousePressHandler(event);
 
  // qDebug() << "The current m_lastZ value:" << m_context.m_lastZ;
 
  emit mousePressEventSignal(m_context);
 
  // qDebug() << "Exiting after having emitted mousePressEventSignal with context:"
  // << m_context.toString();
}
 
 
//! React to the release of the mouse buttons.
void
MassSpecTracePlotWidget::mouseReleaseHandler(QMouseEvent *event)
{
  // qDebug() << "Entering";
 
  // Before working on the various data belonging to the base context, we need
  // to get it from the base class and refresh our local context with it.
  refreshBaseContext();
 
  BasePlotWidget::mouseReleaseHandler(event);
 
  // qDebug() << "The current m_lastZ value:" << m_context.m_lastZ;
 
  if(!m_context.m_isMouseDragging)
    {
      // Since we were not dragging, by essence we were not measuring
      // distances and _a fortiori_ we were not deconvoluting. So
      // reset deconvolution data in the context.
 
      m_context.resetDeconvolutionData();
    }
 
  emit mouseReleaseEventSignal(m_context);
 
  // qDebug() << "Exiting after having emitted mouseReleaseEventSignal with context:"
  // << m_context.toString();
}
 
 
const MassSpecTracePlotContext &
MassSpecTracePlotWidget::refreshBaseContext() const
{
  // BasePlotWidget has a member m_context of type BasePlotContext.
  // Here we also have a m_context *distinct* member that is of type
  // MassSpecTracePlotContext.
 
  // While MassSpecTracePlotContext is derived from BasePlotContext, the two
  // m_context members are distinct and there are lots of housekeeping data that
  // are managed by the parent BasePlotWidget class in its m_context member
  // *independently* of what we have in the ::BasePlotContext part of our
  // m_context member that is of type MassSpecTracePlotContext. We thus need to
  // resynchronize the data from BasePlotWidget::m_context to our m_context.
 
  // qDebug().noquote() << "The base plot context:" <<
  // BasePlotWidget::m_context.toString();
 
 
  // qDebug() << "Right before refreshing the context, the current m_lastZ value:"
  //          << m_context.m_lastZ;
 
  m_context.initialize(BasePlotWidget::m_context);
 
 
  // qDebug() << "Right after refreshing the context, the current m_lastZ value:"
  //          << m_context.m_lastZ;
 
  // qDebug().noquote() << "After refreshing the base context, base context:"
  //                    << m_context.toString();
 
  return m_context;
}
 
 
void
MassSpecTracePlotWidget::setChargeMinimalFractionalPart(double charge_fractional_part)
{
  m_chargeMinimalFractionalPart = charge_fractional_part;
}
 
 
double
MassSpecTracePlotWidget::getChargeMinimalFractionalPart() const
{
  return m_chargeMinimalFractionalPart;
}
 
 
void
MassSpecTracePlotWidget::setChargeStateEnvelopePeakSpan(int interval)
{
  m_chargeStateEnvelopePeakSpan = interval;
}
 
 
int
MassSpecTracePlotWidget::getChargeStateEnvelopePeakSpan() const
{
  return m_chargeStateEnvelopePeakSpan;
}
 
 
//! Deconvolute the mass peaks into charge and molecular mass.
bool
MassSpecTracePlotWidget::deconvolute()
{
 
  // qDebug() << "The current m_lastZ value:" << m_context.m_lastZ;
 
  // There are two situations: when the user is deconvoluting on the
  // basis of the distance between two consecutive peaks of a same
  // isotopic cluster or when the user deconvolutes on the basis of two
  // different charged-stated peaks that belong to the same envelope.
 
  // We can tell the difference because in the first case the xDelta
  // should be less than 1. In the other case, of course the difference
  // is much greater than 1.
 
  // In order to do the deconvolutions, we need to know what is the tolerance
  // on the fractional part of the deconvoluted charge value. This value is set
  // in the parent window's double spin box.
 
  if(fabs(m_context.m_xDelta) >= 0 && fabs(m_context.m_xDelta) <= 1.1)
    {
      // qDebug() << "m_xDelta:" << m_context.m_xDelta
      //<< "trying isotope-based deconvolution.";
 
      return deconvoluteIsotopicCluster();
    }
 
  // If not deconvoluting on the basis of the isotopic cluster, then:
 
  return deconvoluteChargedState(m_chargeStateEnvelopePeakSpan);
}
 
 
//! Deconvolute the mass peaks into charge and molecular mass.
/*!
 
  This is one of two methods to deconvolute mass data into a charge value and
  a Mr value. The method implemented in this function is based on the charge
  state envelope offered by the mass spectrum (most often for polymers of a
  reasonable size).
 
  \param span value representing the number of peaks of the charge state
  envelope that are spanned by the user selection. Defaults to 1, that is, the
  span encompasses two \e consecutive mass peaks of a given charge state
  envelope.
 
  Set m_lastMz, m_lastZ and m_lastMass.
 
  \return true if the deconvolution could be performed, false otherwise.
  */
bool
MassSpecTracePlotWidget::deconvoluteChargedState(int span)
{
  // We assume that we are dealing with two successive (if span is 1) mass
  // peaks belonging to a given charge state family.
 
  // We call span the number of intervals in a given charge state envelope
  // that separate the initial peak (lowerMz) from the last peak (upperMz).
  // That parameter defaults to 1, that is the two peaks are immediately
  // consecutive, that is, there is only one interval.
 
  qDebug() << "The span is:" << span;
 
  // The m_context ranges are inherently sorted. This is not something that
  // is favorable to us because we need to know what is the direction of the
  // drag. That should be configurable, but we can establish that by
  // default, the z value will be calculated for the mass peak at which
  // the mouse drag movement started.
 
  // These two values reflect a m/z range:
  // startMz and curMz are always sorted (startMz < curMz)
  double startMz = m_context.m_xRegionRangeStart;
  double curMz   = m_context.m_xRegionRangeEnd;
 
  qDebug() << "startMz:" << startMz << "curMz:" << curMz;
 
  if(startMz == curMz)
    {
      m_context.m_lastZ            = std::numeric_limits<quint16>::max();
      m_context.m_lastMz           = qQNaN();
      m_context.m_lastTicIntensity = qQNaN();
      m_context.m_lastMr           = qQNaN();
 
      qDebug() << "Start and Cur are the same, no deconvolution possible.";
      return false;
    }
 
  // After much discussing, the agreement we reached with users is that the
  // m/z value that is reported has to be that of the clicked peak, that is
  // the start point of the mouse cursor drag operation. The z charge value
  // has to be related to that m/z value.
 
  double chargeTemp = ((startMz * span) - span) / (curMz - startMz);
  // This is the charge of the least charged peak, that is, by necessity
  // the charge of the peak on the right (curMz, because the values
  // are sorted).
 
  // Make a judicious roundup.
 
  double chargeIntPart;
  double chargeFracPart = modf(chargeTemp, &chargeIntPart);
 
  // When calculating the charge of the ion, very rarely does it provide a
  // perfect integer value. Most often (if deconvolution is for bona fide
  // peaks belonging to the same charge state envelope) that value is with
  // either a large fractional part or a very small fractional part. What we
  // test here, it that fractional part. If it is greater than
  // m_chargeMinimalFractionalPart, then we simply round up to the next integer
  // value (that is, chargeIntPart = 27 and chargeFracPart 0.995, then we
  // set charge to 28). If it is lesser or equal to (1 -
  // m_chargeMinimalFractionalPart /* that is >= 0.01 */, then we let
  // chargeIntPart unmodified (that is, chargeIntPart = 29 and
  // chargeFracPart 0.01, then we set charge to 29). If chargeFracPart is in
  // between (1 - m_chargeMinimalFractionalPart) and
  // m_chargeMinimalFractionalPart, then we consider that the peaks do not
  // belong to the same charge state envelope.
 
  // qDebug() << "Charge:" << chargeIntPart << "Charge fractional part: " << chargeFracPart;
 
 
  if(chargeFracPart >= (1 - m_chargeMinimalFractionalPart /* that is >= 0.01 */) &&
     chargeFracPart <= m_chargeMinimalFractionalPart /* that is <= 0.99 */)
    {
      m_context.m_lastZ            = std::numeric_limits<quint16>::max();
      m_context.m_lastMz           = qQNaN();
      m_context.m_lastTicIntensity = qQNaN();
      m_context.m_lastMr           = qQNaN();
 
      // qDebug() << "Not a charge state family peak,"
      //          << "returning from deconvoluteChargeState";
 
      return false;
    }
 
  if(chargeFracPart > m_chargeMinimalFractionalPart)
    m_context.m_lastZ = chargeIntPart + 1;
  else
    m_context.m_lastZ = chargeIntPart;
 
  qDebug() << "Initially computed m_lastZ as " << m_context.m_lastZ;
 
  // We know the that formula above computed the charge of the least charged
  // ion, that corresponds to the greater m/z value of the range: curMz.
 
  // Now, if the user did click on that least-charged peak, nothing is
  // to be modified. But, if the user did click on the most-charged
  // peak, that is on startMz, the we need to fix the values:
 
  if(static_cast<int>(m_context.m_dragDirections) & static_cast<int>(DragDirections::RIGHT_TO_LEFT))
    {
      qDebug() << "Dragging from right to left, that is from higher m/z value.";
      m_context.m_lastMz = curMz;
    }
  else if(static_cast<int>(m_context.m_dragDirections) &
          static_cast<int>(DragDirections::LEFT_TO_RIGHT))
    {
      qDebug() << "Dragging from left to right, we need to report"
               << "the mass for the left peak, that is for z+span and lower m/z value.";
      m_context.m_lastZ += span;
      m_context.m_lastMz = startMz;
    }
 
  qDebug() << "Final m_lastZ is " << m_context.m_lastZ << "m/z: m_context.m_lastMz";
 
  m_context.m_lastMr = (m_context.m_lastMz * m_context.m_lastZ) - (m_context.m_lastZ * MPROTON);
 
  qDebug() << "startMz:" << QString("%1").arg(startMz, 0, 'f', 6)
           << "curMz:" << QString("%1").arg(curMz, 0, 'f', 6)
           << "m_lastMz (make computation for this startMz):"
           << QString("%1").arg(m_context.m_lastMz, 0, 'f', 6)
           << "m_lastZ:" << QString("%1").arg(m_context.m_lastZ)
           << "m_lastMass:" << QString("%1").arg(m_context.m_lastMr, 0, 'f', 6);
 
  // The m_context was refreshed with the base class context in the calling
  // chain.
  emit massDeconvolutionSignal(m_context);
 
  return true;
}
 
 
//! Deconvolute the mass peaks into charge and molecular mass.
/*!
 
  This is one of two methods to deconvolute mass data into a charge value and
  a Mr value. The method implemented in this function is based on the distance
  that separates two immediately consecutive peaks of an isotopic cluster.
  This method can be used as long as the instrument produced data with a
  resolution sufficient to separate reasonably well the different peaks of an
  isotopic cluster.
 
  Set m_lastMz, m_lastZ and m_lastMass.
 
  \return true if the deconvolution could be performed, false otherwise.
  */
bool
MassSpecTracePlotWidget::deconvoluteIsotopicCluster()
{
 
  double startMz = m_context.m_xRegionRangeStart;
  double curMz   = m_context.m_xRegionRangeEnd;
 
  qDebug() << "startMz:" << startMz << "curMz:" << curMz;
 
  if(startMz == curMz)
    {
      m_context.m_lastZ            = std::numeric_limits<quint16>::max();
      m_context.m_lastMz           = qQNaN();
      m_context.m_lastTicIntensity = qQNaN();
      m_context.m_lastMr           = qQNaN();
 
      return false;
    }
 
  double chargeTemp = 1 / fabs(m_context.m_xDelta);
 
  qDebug() << qSetRealNumberPrecision(6) << "curMz - startMz = " << fabs(curMz - startMz)
           << "while m_context.m_xDelta = " << m_context.m_xDelta;
 
  // Sanity check
  if(fabs(m_context.m_xDelta) != fabs(curMz - startMz))
    qFatal() << "Programming error.";
 
  // Make a judicious roundup.
  double chargeIntPart;
  double chargeFracPart = modf(chargeTemp, &chargeIntPart);
 
  // qDebug() << "m_xDelta:" << m_context.m_xDelta
  //          << "chargeTemp:" << QString("%1").arg(chargeTemp, 0, 'f', 6)
  //          << "chargeIntPart:" << chargeIntPart
  //          << "chargeFracPart:" << QString("%1").arg(chargeFracPart, 0, 'f', 6)
  //          << "m_chargeMinimalFractionalPart:" << m_chargeMinimalFractionalPart;
 
  if(chargeFracPart >= (1 - m_chargeMinimalFractionalPart /* that is >= 0.01 */) &&
     chargeFracPart <= m_chargeMinimalFractionalPart /* that is <= 0.99 */)
    {
      m_context.m_lastZ            = std::numeric_limits<quint16>::max();
      m_context.m_lastMz           = qQNaN();
      m_context.m_lastTicIntensity = qQNaN();
      m_context.m_lastMr           = qQNaN();
 
      // qDebug() << "Not in a isotopic cluster peak:"
      //          << "returning from deconvoluteIsotopicCluster";
 
      return false;
    }
 
  if(chargeFracPart > m_chargeMinimalFractionalPart)
    {
      m_context.m_lastZ = chargeIntPart + 1;
 
      // qDebug() << "chargeFracPart > m_chargeMinimalFractionalPart -> m_lastZ
      // = "
      //<< m_context.m_lastZ;
    }
  else
    {
      m_context.m_lastZ = chargeIntPart;
 
      // qDebug()
      //<< "chargeFracPart <=  m_chargeMinimalFractionalPart -> m_lastZ = "
      //<< m_context.m_lastZ;
    }
 
  // qDebug() << "Settled on m_lastZ: " << m_context.m_lastZ;
 
  // Now that we have the charge in the form of an int, we can compute the
  // Mr of the lightest isotopic cluster peak (the one that has the lowest x
  // value). That value is stored in m_xRangeLower.
 
  // We need to sort the xRegionRange before being certain that lower is indeed
  // the left value of the drag span.
 
  m_context.m_lastMz = std::min<double>(m_context.m_xRegionRangeStart, m_context.m_xRegionRangeEnd);
 
  m_context.m_lastMr = (m_context.m_lastMz * m_context.m_lastZ) - (m_context.m_lastZ * MPROTON);
 
  // qDebug() << "m_lastMz :" << QString("%1").arg(m_context.m_lastMz, 0, 'f', 6)
  //          << "m_lastZ:" << QString("%1").arg(m_context.m_lastZ)
  //          << "m_lastMr:" << QString("%1").arg(m_context.m_lastMr, 0, 'f', 6);
 
  // The m_context was refreshed with the base class context in the calling
  // chain.
  emit massDeconvolutionSignal(m_context);
 
  return true;
}
 
 
bool
MassSpecTracePlotWidget::computeResolvingPower()
{
 
  // m_xRangeLower and m_xRangeUpper and m_xDelta (in fabs() form) have been set
  // during mouve movement handling. Note that the range values *are
  // sorted*.
 
  if(!m_context.m_xDelta || m_context.m_xDelta > 1)
    {
      m_context.m_lastResolvingPower = std::numeric_limits<double>::min();
 
      return false;
    }
 
  // Resolving power is m/z / Delta(m/z), for singly-charged species.
 
  // qDebug() << "Calculating the resolving power with these data:"
  //<< "m_context.m_xRegionRangeStart: " << m_context.m_xRegionRangeStart
  //<< "m_context.m_xRegionRangeEnd: " << m_context.m_xRegionRangeEnd
  //<< "m_context.m_xDelta / 2: " << m_context.m_xDelta / 2;
 
  m_context.m_lastResolvingPower =
    (std::min<double>(m_context.m_xRegionRangeStart, m_context.m_xRegionRangeEnd) +
     (m_context.m_xDelta / 2)) /
    m_context.m_xDelta;
 
  // The m_context was refreshed with the base class context in the calling
  // chain.
  emit resolvingPowerComputationSignal(m_context);
 
  return true;
}
 
 
} // namespace pappso