timsdatafastmap.cpp

Go to the documentation of this file.

/**
 * \file pappsomspp/vendors/tims/timsdatafastmap.h
 * \date 16/12/2023
 * \author Olivier Langella
 * \brief replacement fot std::map dedicated to tims data for fast computations
 */
 
/*******************************************************************************
 * Copyright (c) 2023 Olivier Langella
 *<Olivier.Langella@universite-paris-saclay.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 "timsdatafastmap.h"
#include <QDebug>
#include "pappsomspp/core/exception/exceptionoutofrange.h"
 
namespace pappso
{
 
 
std::map<QThread *, TimsDataFastMap>
  TimsDataFastMap::m_preAllocatedFastMapPerThread = {};
 
pappso::TimsDataFastMap &
TimsDataFastMap::getTimsDataFastMapInstance()
{
 
  auto it = TimsDataFastMap::m_preAllocatedFastMapPerThread.insert(
    {QThread::currentThread(), TimsDataFastMap()});
 
  if(it.second)
    {
      it.first->second.mapTofIndexIntensity.resize(500000);
    }
  return it.first->second;
}
 
TimsDataFastMap::TimsDataFastMap()
{
  // map.resize(500000);
}
 
std::size_t
TimsDataFastMap::accumulateIntensity(quint32 key, std::size_t intensity)
{
  qDebug();
  try
    {
      TimsDataFastMapElement &map_element = mapTofIndexIntensity.at(key);
 
      if(map_element.first_access)
        {
          map_element.first_access = false;
          map_element.count        = intensity;
          tofIndexList.push_back(key);
        }
      else
        {
          map_element.count += intensity;
        }
      qDebug();
      return map_element.count;
    }
  catch(std::out_of_range &error)
    {
      throw pappso::ExceptionOutOfRange(
        QObject::tr("out of range exception for tof index %1 ").arg(key));
    }
}
 
std::size_t
TimsDataFastMap::readIntensity(quint32 key)
{
  TimsDataFastMapElement &map_element = mapTofIndexIntensity.at(key);
  // FIXME: why set first_access to true below?
  map_element.first_access = true;
  return map_element.count;
}
 
void
TimsDataFastMap::downsizeMzRawMap(std::size_t mzindex_merge_window)
{
  std::vector<std::pair<quint32, std::size_t>> temp_vector;
  for(quint32 tof_index : tofIndexList)
    {
      temp_vector.push_back({tof_index, readIntensity(tof_index)});
    }
 
  tofIndexList.clear();
 
  for(auto &pair_tof_intensity : temp_vector)
    {
 
      quint32 mzkey = (pair_tof_intensity.first / mzindex_merge_window);
      mzkey = (mzkey * mzindex_merge_window) + (mzindex_merge_window / 2);
 
      accumulateIntensity(mzkey, pair_tof_intensity.second);
    }
}
 
void
TimsDataFastMap::builtInCentroid()
{
  qDebug() << "tofIndexList.size()=" << tofIndexList.size();
  if(tofIndexList.size() > 2)
    {
      std::vector<quint32> tof_index_list_tmp = tofIndexList;
      std::sort(tof_index_list_tmp.begin(), tof_index_list_tmp.end());
 
      tofIndexList.clear();
 
      quint32 previous_tof_index     = tof_index_list_tmp[0];
      std::size_t previous_intensity = readIntensity(previous_tof_index);
      for(std::size_t i = 1; i < tof_index_list_tmp.size(); i++)
        {
          quint32 tof_index = tof_index_list_tmp[i];
          if(previous_tof_index == tof_index - 1)
            {
              std::size_t intensity = readIntensity(tof_index);
              if(previous_intensity > intensity)
                {
                  // flush writing current accumulated intensity
                  accumulateIntensity(previous_tof_index,
                                      previous_intensity + intensity);
                  previous_intensity = 0;
                  previous_tof_index = tof_index;
                }
              else
                {
                  // accumulate while intensity increases
                  previous_intensity += intensity;
                  previous_tof_index = tof_index;
                }
            }
          else
            {
              // write accumulated intensity :
              if(previous_intensity > 0)
                {
                  // flush
                  accumulateIntensity(previous_tof_index, previous_intensity);
                }
              previous_tof_index = tof_index;
              previous_intensity = readIntensity(tof_index);
            }
        }
 
      // write remaining intensity :
      if(previous_intensity > 0)
        {
          // flush
          accumulateIntensity(previous_tof_index, previous_intensity);
        }
    }
 
  qDebug() << "tofIndexList.size()=" << tofIndexList.size();
}
 
void
TimsDataFastMap::removeArtefactualSpike()
{
 
  std::vector<quint32> tof_index_list_tmp = tofIndexList;
  tofIndexList.clear();
  for(quint32 tof_index : tof_index_list_tmp)
    {
      if((tof_index > 0) && mapTofIndexIntensity[tof_index - 1].first_access &&
         mapTofIndexIntensity[tof_index + 1].first_access &&
         (mapTofIndexIntensity[tof_index].count == 10))
        {
          // this measure is too small and alone : remove it
          mapTofIndexIntensity[tof_index].first_access = true;
        }
      else
        {
          tofIndexList.push_back(tof_index);
        }
    }
}
 
const std::vector<quint32> &
TimsDataFastMap::getTofIndexList() const
{
  return tofIndexList;
}
 
void
TimsDataFastMap::clear()
{
  tofIndexList.clear();
}
 
 
} // namespace pappso