mzcbormsrunreader.cpp

Go to the documentation of this file.

/**
 * \file pappsomspp/core/msrun/private/mzcbormsrunreader.cpp
 * \date 21/11/2025
 * \author Olivier Langella
 * \brief MSrun file reader for mzcbor
 */
 
/*******************************************************************************
 * Copyright (c) 2025 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 "mzcbormsrunreader.h"
#include "pappsomspp/core/exception/exceptionnotfound.h"
#include "pappsomspp/core/exception/exceptionnotimplemented.h"
#include "pappsomspp/core/processing/cbor/mzcbor/mzcborspectrumcollectionreader.h"
#include "pappsomspp/core/processing/cbor/mzcbor/mzcborbuildindexreader.h"
#include "pappsomspp/core/processing/cbor/mzcbor/mzcborindexreader.h"
#include "pappsomspp/core/processing/cbor/cborstreamreader.h"
#include "pappsomspp/core/processing/uimonitor/uimonitorvoid.h"
#include <qlogging.h>
 
namespace pappso
{
MzcborMsRunReader::MzcborMsRunReader(MsRunIdCstSPtr &msrun_id_csp)
  : pappso::MsRunReader(msrun_id_csp)
{
  initialize();
  qDebug();
}
 
MzcborMsRunReader::~MzcborMsRunReader()
{
  releaseDevice();
}
 
 
const std::vector<qint64> &
MzcborMsRunReader::getSpectrumIndexPositionInFile() const
{
 
  return m_spectrumIndexPositionInFile;
}
 
void
MzcborMsRunReader::initialize()
{
  qDebug();
  m_cborFileInfo.setFile(mcsp_msRunId->getFileName());
 
  if(!m_cborFileInfo.exists())
    throw ExceptionNotFound(m_cborFileInfo.absoluteFilePath());
 
 
  // build the index
  acquireDevice();
 
  QString str_index_file = m_cborFileInfo.absoluteFilePath();
  QFileInfo mzcbor_index_fileinfo(str_index_file.append(".idx"));
  if(mzcbor_index_fileinfo.exists())
    {
      qDebug() << "mzcbor_index_fileinfo.exists()";
      QFile mzcbor_index_file(mzcbor_index_fileinfo.absoluteFilePath());
      mzcbor_index_file.open(QIODevice::ReadOnly);
 
      pappso::cbor::mzcbor::MzcborIndexReader index_reader;
      index_reader.readCbor(&mzcbor_index_file);
      index_reader.close();
      mzcbor_index_file.close();
 
      qDebug();
      std::vector<QString> run_id_list = index_reader.getRunIdList();
      if(run_id_list.size() > 0)
        {
          std::size_t run_position = 0;
          if(run_id_list.size() > 1)
            {
              std::size_t i_run = 0;
              for(const QString &run_id : run_id_list)
                {
                  if(run_id == getMsRunId().get()->getRunId())
                    {
                      run_position = i_run;
                    }
                  i_run++;
                }
            }
 
          m_spectrumIndexPositionInFile =
            index_reader.getRunAndSpectrumOffsetList().at(run_position);
 
 
          if(index_reader.getRunAndSpectrumTotalIonCountList().size() > run_position)
            {
              m_spectrumTotalIonCountList =
                index_reader.getRunAndSpectrumTotalIonCountList().at(run_position);
            }
          if(index_reader.getRunAndSpectrumMsLevelList().size() > run_position)
            {
              m_spectrumMsLevelList = index_reader.getRunAndSpectrumMsLevelList().at(run_position);
            }
 
          if(index_reader.getRunAndSpectrumRtList().size() > run_position)
            {
              m_spectrumRtList = index_reader.getRunAndSpectrumRtList().at(run_position);
            }
          // qFatal() << m_spectrumTotalIonCountList.size();
 
          std::size_t i = 0;
          for(const QString &native_id : index_reader.getRunAndSpectrumIdList().at(run_position))
            {
              m_nativeId2SpectrumIndexMap.insert({native_id, i});
              i++;
            }
        }
    }
 
  if(m_spectrumIndexPositionInFile.size() == 0)
    {
 
      pappso::UiMonitorVoid monitor;
      pappso::cbor::mzcbor::MzcborBuildIndexReader mzcbor_build_index_reader(getMsRunId());
      mzcbor_build_index_reader.readCbor(mpa_mzcborFileDevice, monitor);
 
 
      qDebug();
      mzcbor_build_index_reader.close();
      // releaseDevice();
      qDebug();
 
      m_spectrumIndexPositionInFile = mzcbor_build_index_reader.getSpectrumIndexPositionInFile();
      m_nativeId2SpectrumIndexMap   = mzcbor_build_index_reader.getNativeId2SpectrumIndexMap();
    }
}
 
bool
MzcborMsRunReader::hasScanNumbers() const
{
  if(m_nativeId2SpectrumIndexMap.size() > 0)
    {
      if(m_nativeId2SpectrumIndexMap.begin()->first.contains("scan="))
        return true;
    }
  return false;
}
 
std::size_t
MzcborMsRunReader::scanNumber2SpectrumIndex(std::size_t scan_number)
{
  if(m_scan2SpectrumIndexMap.size() == 0)
    {
      for(auto &index_pair : m_nativeId2SpectrumIndexMap)
        {
          QStringList native_id_list = index_pair.first.split("=");
          if(native_id_list.size() < 2)
            {
            }
          else
            {
              std::size_t scan_number = native_id_list.back().toULong();
              m_scan2SpectrumIndexMap.insert(
                std::pair<std::size_t, std::size_t>(scan_number, index_pair.second));
            }
        }
    }
 
  auto it = m_scan2SpectrumIndexMap.find(scan_number);
 
  if(it == m_scan2SpectrumIndexMap.end())
    {
      throw ExceptionNotFound(QObject::tr("error reading file %1 : scan %2 not found")
                                .arg(mcsp_msRunId.get()->getFileName())
                                .arg(scan_number));
    }
  return it->second;
}
 
bool
MzcborMsRunReader::accept(const QString &file_name [[maybe_unused]]) const
{
 
  throw pappso::ExceptionNotImplemented(
    QObject::tr("%1 %2 %3 not implemented").arg(__FILE__).arg(__FUNCTION__).arg(__LINE__));
}
 
 
const OboPsiModTerm
MzcborMsRunReader::getOboPsiModTermInstrumentModelName() const
{
  throw pappso::ExceptionNotImplemented(
    QObject::tr("%1 %2 %3 not implemented").arg(__FILE__).arg(__FUNCTION__).arg(__LINE__));
}
 
bool
pappso::MzcborMsRunReader::acquireDevice()
{
  if(mpa_mzcborFileDevice == nullptr)
    {
      mpa_mzcborFileDevice = new QFile(m_cborFileInfo.absoluteFilePath());
      mpa_mzcborFileDevice->open(QIODevice::ReadOnly);
    }
  return true;
}
 
 
bool
pappso::MzcborMsRunReader::releaseDevice()
{
  if(mpa_mzcborFileDevice != nullptr)
    {
      mpa_mzcborFileDevice->close();
      delete mpa_mzcborFileDevice;
      mpa_mzcborFileDevice = nullptr;
    }
  return true;
}
 
 
pappso::XicCoordSPtr
pappso::MzcborMsRunReader::newXicCoordSPtrFromQualifiedMassSpectrum(
  const pappso::QualifiedMassSpectrum &mass_spectrum, pappso::PrecisionPtr precision) const
{
 
  XicCoordSPtr xic_coord = std::make_shared<XicCoord>();
 
  xic_coord.get()->rtTarget = mass_spectrum.getRtInSeconds();
  xic_coord.get()->mzRange  = MzRange(mass_spectrum.getPrecursorMz(), precision);
  return xic_coord;
}
 
pappso::XicCoordSPtr
pappso::MzcborMsRunReader::newXicCoordSPtrFromSpectrumIndex(std::size_t spectrum_index,
                                                            pappso::PrecisionPtr precision) const
{
 
  pappso::cbor::mzcbor::Spectrum cbor_spectrum;
 
  fillMzcborSpectrum(spectrum_index, cbor_spectrum, false);
 
  XicCoordSPtr xic_coord = std::make_shared<XicCoord>();
 
  xic_coord.get()->rtTarget = cbor_spectrum.getRtInSeconds();
 
  bool get_mz = false;
  if(cbor_spectrum.precursorList.size() > 0)
    {
 
      for(auto &precursor : cbor_spectrum.precursorList)
        {
          for(auto &ion : precursor.selectedIonList)
            {
              xic_coord.get()->mzRange = MzRange(ion.getMz(), precision);
              get_mz                   = true;
            }
        }
      if(!get_mz)
        {
          throw pappso::ExceptionNotFound(
            QObject::tr("precursor m/z not found for this spectrum index %1").arg(spectrum_index));
        }
    }
  else
    {
      throw pappso::ExceptionNotFound(
        QObject::tr("no precursor found for this spectrum index %1").arg(spectrum_index));
    }
 
  return xic_coord;
}
 
MassSpectrumCstSPtr
pappso::MzcborMsRunReader::massSpectrumCstSPtr(std::size_t spectrum_index)
{
  return massSpectrumSPtr(spectrum_index);
}
 
MassSpectrumSPtr
pappso::MzcborMsRunReader::massSpectrumSPtr(std::size_t spectrum_index)
{
  MassSpectrumSPtr mass_spectrum_sp;
  try
    {
      acquireDevice();
      pappso::cbor::mzcbor::Spectrum cbor_spectrum;
 
      fillMzcborSpectrum(spectrum_index, cbor_spectrum, true);
      if(cbor_spectrum.binaryDataArrayList.size() == 2)
        {
          mass_spectrum_sp = std::make_shared<MassSpectrum>();
          cbor_spectrum.decodeTrace(*(mass_spectrum_sp.get()));
        }
      else
        {
          throw pappso::PappsoException(
            QObject::tr("cbor_spectrum.binaryDataArrayList.size() != 2"));
        }
    }
  catch(const pappso::PappsoException &pappso_error)
    {
      qDebug() << "Going to throw";
 
      throw pappso::PappsoException(QObject::tr("Error reading data (massSpectrumSPtr) using the "
                                                "mzcbor reader: %1")
                                      .arg(pappso_error.what()));
    }
  return mass_spectrum_sp;
}
 
QualifiedMassSpectrum
pappso::MzcborMsRunReader::qualifiedMassSpectrum(std::size_t spectrum_index,
                                                 bool want_binary_data) const
{
  qDebug();
  try
    {
      pappso::cbor::mzcbor::Spectrum cbor_spectrum;
 
      fillMzcborSpectrum(spectrum_index, cbor_spectrum, want_binary_data);
      qDebug() << cbor_spectrum.index;
      QualifiedMassSpectrum qualified_mass_spectrum;
 
      MassSpectrumId spectrum_id(mcsp_msRunId, spectrum_index);
      spectrum_id.setNativeId(cbor_spectrum.id);
 
      spectrum_id.setSpectrumIndex(cbor_spectrum.index);
 
      qualified_mass_spectrum.setMassSpectrumId(spectrum_id);
      qualified_mass_spectrum.setRtInSeconds(cbor_spectrum.getRtInSeconds());
      if(cbor_spectrum.precursorList.size() > 0)
        {
          qualified_mass_spectrum.setPrecursorNativeId(
            cbor_spectrum.precursorList.at(0).spectrumRef);
          qualified_mass_spectrum.setPrecursorSpectrumIndex(
            m_nativeId2SpectrumIndexMap.at(cbor_spectrum.precursorList.at(0).spectrumRef));
          // qWarning() << "pappso::MzcborMsRunReader::qualifiedMassSpectrum";
 
          for(auto &precursor : cbor_spectrum.precursorList)
            {
 
              // *****************************************************
              // get isolation window (very simplified)
              // *****************************************************
 
              try
                {
                  // qWarning() << "pappso::MzcborMsRunReader::qualifiedMassSpectrum";
                  qualified_mass_spectrum.setParameterValue(
                    QualifiedMassSpectrumParameter::IsolationWindowTargetMz,
                    precursor.isolationWindow.getTargetMz());
 
 
                  auto it = precursor.isolationWindow.cvParamMap.find("MS:1000828");
 
                  if(it == precursor.isolationWindow.cvParamMap.end())
                    {
                      throw pappso::ExceptionNotFound(QObject::tr(
                        "isolation window lower offset MS:1000828 not found in cvParam map"));
                    }
                  qualified_mass_spectrum.setParameterValue(
                    QualifiedMassSpectrumParameter::IsolationWindowLowerOffset,
                    it->second.valueDouble);
 
 
                  it = precursor.isolationWindow.cvParamMap.find("MS:1000829");
 
                  if(it == precursor.isolationWindow.cvParamMap.end())
                    {
                      throw pappso::ExceptionNotFound(QObject::tr(
                        "isolation window upper offset MS:1000829 not found in cvParam map"));
                    }
 
                  qualified_mass_spectrum.setParameterValue(
                    QualifiedMassSpectrumParameter::IsolationWindowUpperOffset,
                    it->second.valueDouble);
                }
              catch(const pappso::ExceptionNotFound &not_found)
                {
                  // isolation window not available
                  // throw not_found;
                }
              // *****************************************************
              // end get isolation window (very simplified)
              // *****************************************************
 
              for(auto &ion : precursor.selectedIonList)
                {
                  PrecursorIonData precursor_ion_data;
                  precursor_ion_data.charge = ion.getChargeState();
                  bool is_ok;
                  precursor_ion_data.intensity = ion.getIntensity(&is_ok);
                  precursor_ion_data.mz        = ion.getMz();
                  qualified_mass_spectrum.appendPrecursorIonData(precursor_ion_data);
                }
            }
        }
 
 
      qDebug();
      qualified_mass_spectrum.setMsLevel(cbor_spectrum.getMsLevel());
 
      qDebug();
      qualified_mass_spectrum.setEmptyMassSpectrum(!cbor_spectrum.defaultArrayLength);
      if(cbor_spectrum.binaryDataArrayList.size() == 2)
        {
          MassSpectrumSPtr mass_spectrum_sp = std::make_shared<MassSpectrum>();
          cbor_spectrum.decodeTrace(*(mass_spectrum_sp.get()));
          qualified_mass_spectrum.setMassSpectrumSPtr(mass_spectrum_sp);
        }
 
      qDebug() << "spectrum id=" << cbor_spectrum.id;
 
      return qualified_mass_spectrum;
    }
 
  catch(const pappso::PappsoException &pappso_error)
    {
      qDebug() << "Going to throw";
 
      throw pappso::PappsoException(
        QObject::tr("Error reading data (qualifiedMassSpectrum) using the "
                    "mzcbor reader: %1")
          .arg(pappso_error.what()));
    }
}
 
void
pappso::MzcborMsRunReader::readSpectrumCollection(SpectrumCollectionHandlerInterface &handler)
{
 
  throw pappso::ExceptionNotImplemented(
    QObject::tr("%1 %2 %3 not implemented").arg(__FILE__).arg(__FUNCTION__).arg(__LINE__));
  MsRunReadConfig config;
  std::vector<size_t> ms_levels;
  for(std::size_t i = 1; i < 9; i++)
    {
      if(handler.needMsLevelPeakList(i))
        {
          ms_levels.push_back(i);
        }
    }
  config.setMsLevels(ms_levels);
  config.setNeedPeakList(handler.needPeakList());
  // readSpectrumCollectionByMsLevel(reader_timeline, 1);
 
  try
    {
      readSpectrumCollectionWithMsrunReadConfig(config, handler);
    }
  catch(const pappso::PappsoException &pappso_error)
    {
      qDebug() << "Going to throw";
 
      throw pappso::PappsoException(
        QObject::tr("Error reading data (spectrum collection2) using the "
                    "mzcbor reader: %1")
          .arg(pappso_error.what()));
    }
  catch(std::exception &error)
    {
      qDebug() << "Going to throw";
 
      throw pappso::PappsoException(
        QObject::tr("Error reading data (spectrum collection) using the "
                    "mzcbor reader: %1")
          .arg(error.what()));
    }
}
 
void
pappso::MzcborMsRunReader::readSpectrumCollection2(const MsRunReadConfig &config,
                                                   SpectrumCollectionHandlerInterface &handler)
{
  qDebug();
  try
    {
      readSpectrumCollectionWithMsrunReadConfig(config, handler);
    }
  catch(const pappso::PappsoException &pappso_error)
    {
      qDebug() << "Going to throw";
 
      throw pappso::PappsoException(
        QObject::tr("Error reading data (spectrum collection2) using the "
                    "mzcbor reader: %1")
          .arg(pappso_error.what()));
    }
  catch(std::exception &error)
    {
      qDebug() << "Going to throw";
 
      throw pappso::PappsoException(
        QObject::tr("Error reading data (spectrum collection2) using the "
                    "mzcbor reader: %1")
          .arg(error.what()));
    }
}
 
void
pappso::MzcborMsRunReader::readSpectrumCollectionByMsLevel(
  SpectrumCollectionHandlerInterface &handler, unsigned int ms_level)
{
  MsRunReadConfig config;
  config.setMsLevels({ms_level});
  config.setNeedPeakList(handler.needPeakList());
  // readSpectrumCollectionByMsLevel(reader_timeline, 1);
  readSpectrumCollection2(config, handler);
}
 
std::size_t
pappso::MzcborMsRunReader::spectrumListSize() const
{
  return m_spectrumIndexPositionInFile.size();
}
 
std::size_t
pappso::MzcborMsRunReader::spectrumStringIdentifier2SpectrumIndex(
  const QString &spectrum_identifier)
{
  auto it = m_nativeId2SpectrumIndexMap.find(spectrum_identifier);
  if(it == m_nativeId2SpectrumIndexMap.end())
    {
      throw pappso::ExceptionNotFound(
        QObject::tr("spectrum identifier %1 not found").arg(spectrum_identifier));
    }
  return it->second;
}
 
 
void
MzcborMsRunReader::readSpectrumCollectionWithMsrunReadConfig(
  const MsRunReadConfig &config, SpectrumCollectionHandlerInterface &handler)
{
  // acquireDevice();
  try
    {
      acquireDevice();
      pappso::UiMonitorVoid monitor;
      pappso::cbor::mzcbor::MzcborSpectrumCollectionReader mzcbor_spectrum_collection_reader(
        config, handler);
      mzcbor_spectrum_collection_reader.setMsRunId(getMsRunId());
      mzcbor_spectrum_collection_reader.setNativeId2SpectrumIndexMapPtr(
        &m_nativeId2SpectrumIndexMap);
      mzcbor_spectrum_collection_reader.readCbor(mpa_mzcborFileDevice, monitor);
 
 
      qDebug();
      mzcbor_spectrum_collection_reader.close();
      qDebug();
    }
  catch(const pappso::PappsoException &pappso_err)
    {
      throw pappso::PappsoException(
        QObject::tr("ERROR in MzcborMsRunReader::readSpectrumCollectionWithMsrunReadConfig:\n%1")
          .arg(pappso_err.qwhat()));
    }
  // End of
  // for(std::size_t iter = 0; iter < spectrum_list_size; iter++)
 
  // Now let the loading handler know that the loading of the data has ended.
  // The handler might need this "signal" to perform additional tasks or to
  // cleanup cruft.
}
 
 
void
MzcborMsRunReader::fillMzcborSpectrum(std::size_t spectrum_index,
                                      pappso::cbor::mzcbor::Spectrum &spectrum,
                                      bool want_binary_data) const
{
 
  if(spectrum_index >= m_spectrumIndexPositionInFile.size())
    {
      throw pappso::ExceptionNotFound(
        QObject::tr("spectrum index %1 not found").arg(spectrum_index));
    }
  if(mpa_mzcborFileDevice == nullptr)
    {
 
      throw pappso::PappsoException(
        QObject::tr("mzCBOR file device is not ready, use acquireDevice() before access"));
    }
  qDebug() << spectrum_index << " " << m_spectrumIndexPositionInFile[spectrum_index];
 
 
  // qDebug() << spectrum_index << " " << m_spectrumIndexPositionInFile[spectrum_index-1];
 
 
  if(mpa_mzcborFileDevice->seek(m_spectrumIndexPositionInFile[spectrum_index]))
    {
      pappso::cbor::CborStreamReader cbor_stream_reader(mpa_mzcborFileDevice);
      cbor_stream_reader.setDevice(mpa_mzcborFileDevice);
 
 
      try
        {
          spectrum.fromCbor(cbor_stream_reader, want_binary_data);
        }
      catch(const pappso::PappsoException &error)
        {
 
          throw pappso::PappsoException(QObject::tr("ERROR reading spectrum %1 cbor:\n%2")
                                          .arg(spectrum_index)
                                          .arg(error.qwhat()));
        }
      // cbor_stream_reader.leaveContainer();
 
      if(spectrum.id.isEmpty())
        {
          throw pappso::PappsoException(
            QObject::tr("ERROR reading spectrum %1 cbor:\nspectrum not found").arg(spectrum_index));
        }
    }
  else
    {
      throw pappso::PappsoException(
        QObject::tr("ERROR reading spectrum %1 cbor:\nseek failed").arg(spectrum_index));
    }
  // mzcbor_file.close();
}
 
 
Trace
MzcborMsRunReader::getTicChromatogram()
{
  qDebug() << m_spectrumMsLevelList.size() << " " << m_spectrumTotalIonCountList.size();
  Trace chromatogram;
  if((m_spectrumTotalIonCountList.size() > 0) &&
     (m_spectrumMsLevelList.size() == m_spectrumTotalIonCountList.size()))
    {
      for(std::size_t i = 0; i < m_spectrumTotalIonCountList.size(); i++)
        {
          if(m_spectrumMsLevelList.at(i) == 1)
            {
              chromatogram.push_back(
                {m_spectrumRtList.at(i), (double)m_spectrumTotalIonCountList.at(i)});
            }
        }
    }
  else
    {
      return MsRunReader::getTicChromatogram();
    }
  return chromatogram;
}
std::vector<double>
MzcborMsRunReader::getRetentionTimeLine()
{
  qDebug();
  std::vector<double> time_line;
  if((m_spectrumRtList.size() > 0) && (m_spectrumMsLevelList.size() == m_spectrumRtList.size()))
    {
      for(std::size_t i = 0; i < m_spectrumRtList.size(); i++)
        {
          if(m_spectrumMsLevelList.at(i) == 1)
            {
              time_line.push_back(m_spectrumRtList.at(i));
            }
        }
    }
  else
    {
      return MsRunReader::getRetentionTimeLine();
    }
  return time_line;
}
 
std::shared_ptr<pappso::cbor::mzcbor::Spectrum>
MzcborMsRunReader::getMzcborSpectrumSp(std::size_t spectrum_index, bool want_binary_data) const
{
  std::shared_ptr<pappso::cbor::mzcbor::Spectrum> cbor_spectrum_sp =
    std::make_shared<pappso::cbor::mzcbor::Spectrum>();
 
  fillMzcborSpectrum(spectrum_index, *cbor_spectrum_sp, want_binary_data);
 
  return cbor_spectrum_sp;
}
 
} // namespace pappso