Loading...
Searching...
No Matches
massspectrumminuscombiner.cpp
Go to the documentation of this file.
/////////////////////// StdLib includes
#include <numeric>
#include <limits>
#include <vector>
#include <map>
#include <cmath>
#include <iostream>
#include <iomanip>
/////////////////////// Qt includes
#include <QDebug>
#include <QFile>
#include <QThread>
#if 0
// For debugging purposes.
#include <QFile>
#endif
/////////////////////// Local includes
#include "massspectrumminuscombiner.h"
#include "../../types.h"
#include "../../utils.h"
#include "../../pappsoexception.h"
#include "../../exception/exceptionoutofrange.h"
#include "../../exception/exceptionnotpossible.h"
namespace pappso
{
//! Construct an uninitialized instance.
MassSpectrumMinusCombiner::MassSpectrumMinusCombiner()
{
}
MassSpectrumMinusCombiner::MassSpectrumMinusCombiner(int decimal_places)
: MassSpectrumCombiner(decimal_places)
{
}
MassSpectrumMinusCombiner::MassSpectrumMinusCombiner(
const MassSpectrumMinusCombiner &other)
: MassSpectrumCombiner(other)
{
}
MassSpectrumMinusCombiner::MassSpectrumMinusCombiner(
MassSpectrumMinusCombinerCstSPtr other)
: MassSpectrumCombiner(other)
{
}
//! Destruct the instance.
MassSpectrumMinusCombiner::~MassSpectrumMinusCombiner()
{
}
MassSpectrumMinusCombiner &
MassSpectrumMinusCombiner::operator=(const MassSpectrumMinusCombiner &other)
{
if(this == &other)
return *this;
m_decimalPlaces = other.m_decimalPlaces;
m_bins.assign(other.m_bins.begin(), other.m_bins.end());
return *this;
}
MapTrace &
MassSpectrumMinusCombiner::combine(MapTrace &map_trace,
const Trace &trace) const
{
// qDebug();
if(!trace.size())
{
// qDebug() << "Thread:" << QThread::currentThreadId()
//<< "Returning right away because trace is empty.";
return map_trace;
}
// We will need to only use these iterator variables if we do not want to
// loose consistency.
using TraceIter = std::vector<DataPoint>::const_iterator;
TraceIter trace_iter_begin = trace.begin();
TraceIter trace_iter = trace_iter_begin;
TraceIter trace_iter_end = trace.end();
// The destination map trace will be filled-in with the result of the
// combination.
// Sanity check:
if(!m_bins.size())
throw(ExceptionNotPossible("The bin vector cannot be empty."));
using BinIter = std::vector<pappso_double>::const_iterator;
BinIter bin_iter = m_bins.begin();
BinIter bin_end_iter = m_bins.end();
// qDebug() << "initial bins iter at a distance of:"
//<< std::distance(m_bins.begin(), bin_iter)
//<< "bins distance:" << std::distance(m_bins.begin(), m_bins.end())
//<< "bins size:" << m_bins.size() << "first bin:" << m_bins.front()
//<< "last bin:" << m_bins.back();
// Iterate in the vector of bins and for each bin check if there are matching
// data points in the trace.
pappso_double current_bin_mz_value = 0;
pappso_double trace_x = 0;
pappso_double trace_y = 0;
// Lower bound returns an iterator pointing to the first element in the
// range [first, last) that is not less than (i.e. greater or equal to)
// value, or last if no such element is found.
// Get the first bin that would match the very first data point in the trace
// that we need to combine into the map_trace.
auto bin_iter_for_mz = lower_bound(bin_iter, bin_end_iter, trace_iter->x);
if(bin_iter_for_mz != bin_end_iter)
{
// If we found a bin and that the bin is not the first bin of the bin
// vector, then go back one bin to be sure we do not miss data points.
if(bin_iter_for_mz != m_bins.begin())
bin_iter = --bin_iter_for_mz;
}
else
throw(ExceptionNotPossible("The bin vector must match the mz value."));
// Now iterate in the bin vector starting from the std::prev(found bin).
while(bin_iter != bin_end_iter)
{
current_bin_mz_value = *bin_iter;
// qDebug() << "Iterating in new bin:"
//<< QString("%1").arg(current_bin_mz_value, 0, 'f', 15)
//<< "at a distance from the first bin of:"
//<< std::distance(m_bins.begin(), bin_iter);
// For the current bin, we start by instantiating a new DataPoint. By
// essence, each bin will have at most one corresponding DataPoint.
DataPoint bin_data_point;
// Do not set the y value to 0 so that we can actually test if the
// data point is valid later on (try not to push back y=0 data
// points).
bin_data_point.x = current_bin_mz_value;
// qDebug() << "Seed bin_data_point.x with current_bin_mz_value value:"
//<< QString("%1").arg(bin_data_point.x, 0, 'f', 15);
// Now perform a loop over the data points in the mass spectrum.
while(trace_iter != trace_iter_end)
{
bool trace_matched = false;
// If we are not at the end of trace and if the y value of the
// currently iterated trace datapoint is not 0, perform the rounding
// and check if the obtained x value is in the current bin, that is if
// it is less or equal to the current bin.
// qDebug() << "Thread:" << QThread::currentThreadId();
// qDebug() << "Iterating in new trace datapoint:"
//<< trace_iter->toString(15) << "at distance:"
//<< std::distance(trace_iter_begin, trace_iter);
if(trace_iter->y)
{
// qDebug() << "The y value of trace iterated value is not 0.";
// trace_x is the m/z value that we need to combine,
trace_x = trace_iter->x;
trace_y = trace_iter->y;
// Now apply the rounding (if any).
if(m_decimalPlaces != -1)
trace_x = Utils::roundToDecimals(trace_x, m_decimalPlaces);
if(trace_x <= current_bin_mz_value)
{
// qDebug()
//<< "trace_x <= current_bin_mz_value: the match happened.";
// qDebug() << "Going to compute" << bin_data_point.y << "-"
//<< trace_y << "with current_bin_mz_value:"
//<< QString("%1").arg(
// current_bin_mz_value, 0, 'f', 15);
// This is where the bin_data_point actually gets decremented
// by the y value of the currently iterated trace data point.
// So bin_data_point has already the MINUS'ed value in it.
// Relate this code line with the result.first->second +=
// bin_data_point.y; code line elsewhere in this function.
bin_data_point.y -= trace_y;
// qDebug() << "New bin_data_point.y value:" <<
// bin_data_point.y;
// Let's record that we matched.
trace_matched = true;
// Because we matched, we can step-up with the
// iterator.
// qDebug() << "The values matched, increment the trace "
//"itereator but not the bin iterator, as maybe "
//"the next trace datapoint also matches the bin.";
++trace_iter;
}
// else
//{
// We did have a non-0 y value, but that did not
// match. So we do not step-up with the iterator because that
// trace data point will fit into another bin not yet iterated
// into.
//}
}
// End of
// if(trace_iter->y)
else
{
// We iterated into a y=0 data point, so just skip it. Let the
// below code think that we have matched the point and iterate one
// step up.
// qDebug() << "The y value of the currently iterated trace data "
//"point is almost equal to 0, increment the "
//"trace iter but do nothing else.";
trace_matched = true;
++trace_iter;
}
// At this point, check if the trace data point matched the currently
// iterated bin.
if(!trace_matched)
{
// qDebug() << "The currently iterated trace datapoint did not "
//"match the current bin.";
// The trace data point did not match the currently iterated bin.
// All we have to do is go to the next bin. We break and the bin
// vector iterator will be incremented.
// However, if we had a valid new data point, that
// data point needs to be pushed back in the new mass
// spectrum.
if(bin_data_point.isValid())
{
// qDebug() << "But we had a valid bin data point cooking, so
// " "we have to take that into account:"
//<< bin_data_point.toString(15);
// We need to check if that bin value is present already in
// the map_trace object passed as parameter.
std::pair<std::map<pappso_double, pappso_double>::iterator,
bool>
result =
map_trace.insert(std::pair<pappso_double, pappso_double>(
bin_data_point.x, -bin_data_point.y));
if(!result.second)
{
// qDebug()
//<< "The x value (mz) was already present in the "
//"map trace, simply update its y value (intensity).";
// The key already existed! The item was not inserted. We
// need to update the y value.
// qDebug() << "Going to compute" << result.first->second
//<< "+" << bin_data_point.y;
// This might be counterintuitive: one would expect the
// use of the '-=' operator because we are
// MINUS-combining. But no! bin_data_point only get y
// values already minus'ed. So we need to add that
// minus'ed value, otherwise we add the value (- -value is
// the same as +value).
result.first->second += bin_data_point.y;
// qDebug() << "New y value for the item in the map trace
// "
//"(result.first->second value):"
//<< result.first->second;
}
// else
//{
// qDebug()
//<< "The data point has a mz value that was not present "
//"in the map trace."
//<< " Inserted a new data point into the map trace:"
//<< bin_data_point.x << "," << bin_data_point.y;
//}
}
// We need to break this loop! That will increment the
// bin iterator.
break;
}
}
// End of
// while(trace_iter != trace_iter_end)
// Each time we get here, that means that we have consumed all
// the mass spectra data points that matched the current bin.
// So go to the next bin.
if(trace_iter == trace_iter_end)
{
// Make sure a last check is done in case one data point was
// cooking...
if(bin_data_point.isValid())
{
std::pair<std::map<pappso_double, pappso_double>::iterator, bool>
result =
map_trace.insert(std::pair<pappso_double, pappso_double>(
bin_data_point.x, -bin_data_point.y));
if(!result.second)
{
// qDebug()
//<< "The x value (mz) was already present in the "
//"map trace, simply update its y value (intensity).";
// The key already existed! The item was not inserted. We
// need to update the y value.
// qDebug() << "Going to compute" << result.first->second <<
// "+"
//<< bin_data_point.y;
// This might be counterintuitive: one would expect the
// use of the '-=' operator because we are
// MINUS-combining. But no! bin_data_point only get y
// values already minus'ed. So we need to add that
// minus'ed value, otherwise we add the value (- -value is
// the same as +value).
result.first->second += bin_data_point.y;
// qDebug() << "New y value for the item in the map trace "
//"(result.first->second value):"
//<< result.first->second;
}
// else
//{
// qDebug()
//<< "The data point has a mz value that was not present "
//"in the map trace."
//<< " Inserted a new data point into the map trace:"
//<< bin_data_point.x << "," << bin_data_point.y;
//}
}
// Now truly exit the loops...
break;
}
++bin_iter;
}
// End of
// while(bin_iter != bin_end_iter)
// qDebug() << __FILE__ << "@" << __LINE__ << __FUNCTION__ << " ()"
//<< "The combination result mass spectrum being returned has TIC:"
//<< new_trace.sumY();
return map_trace;
}
MapTrace &
MassSpectrumMinusCombiner::combine(MapTrace &map_trace_out,
const MapTrace &map_trace_in) const
{
// qDebug();
if(!map_trace_in.size())
{
// qDebug() << "Thread:" << QThread::currentThreadId()
//<< "Returning right away because map_trace_in is empty.";
return map_trace_out;
}
Trace trace(map_trace_in);
return combine(map_trace_out, trace);
}
} // namespace pappso
pappso::ExceptionNotPossible
Definition exceptionnotpossible.h:29
pappso::MapTrace
Definition maptrace.h:34
pappso::MassDataCombinerInterface::m_decimalPlaces
int m_decimalPlaces
Number of decimals to use for the keys (x values)
Definition massdatacombinerinterface.h:44
pappso::MassSpectrumCombiner
Definition massspectrumcombiner.h:29
pappso::MassSpectrumCombiner::m_bins
std::vector< pappso_double > m_bins
Definition massspectrumcombiner.h:48
pappso::MassSpectrumMinusCombiner
Definition massspectrumminuscombiner.h:26
pappso::MassSpectrumMinusCombiner::combine
virtual MapTrace & combine(MapTrace &map_trace, const Trace &trace) const override
Definition massspectrumminuscombiner.cpp:83
pappso::MassSpectrumMinusCombiner::MassSpectrumMinusCombiner
MassSpectrumMinusCombiner()
Construct an uninitialized instance.
Definition massspectrumminuscombiner.cpp:35
pappso::MassSpectrumMinusCombiner::~MassSpectrumMinusCombiner
virtual ~MassSpectrumMinusCombiner()
Destruct the instance.
Definition massspectrumminuscombiner.cpp:63
pappso::MassSpectrumMinusCombiner::operator=
MassSpectrumMinusCombiner & operator=(const MassSpectrumMinusCombiner &other)
Definition massspectrumminuscombiner.cpp:69
pappso::Trace
A simple container of DataPoint instances.
Definition trace.h:148
pappso::Utils::roundToDecimals
static pappso_double roundToDecimals(pappso_double value, int decimal_places)
Definition utils.cpp:142
massspectrumminuscombiner.h
pappso
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition aa.cpp:39
pappso::MassSpectrumMinusCombinerCstSPtr
std::shared_ptr< const MassSpectrumMinusCombiner > MassSpectrumMinusCombinerCstSPtr
Definition massspectrumminuscombiner.h:20
pappso::pappso_double
double pappso_double
A type definition for doubles.
Definition types.h:50
pappso::DataPoint
Definition datapoint.h:22
pappso::DataPoint::x
pappso_double x
Definition datapoint.h:23
pappso::DataPoint::isValid
bool isValid() const
Definition datapoint.cpp:136
pappso::DataPoint::y
pappso_double y
Definition datapoint.h:24