libpappsomspp
Library for mass spectrometry
timsdata.cpp
Go to the documentation of this file.
1 /**
2  * \file pappsomspp/vendors/tims/timsdata.cpp
3  * \date 27/08/2019
4  * \author Olivier Langella
5  * \brief main Tims data handler
6  */
7 
8 /*******************************************************************************
9  * Copyright (c) 2019 Olivier Langella <Olivier.Langella@u-psud.fr>.
10  *
11  * This file is part of the PAPPSOms++ library.
12  *
13  * PAPPSOms++ is free software: you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License as published by
15  * the Free Software Foundation, either version 3 of the License, or
16  * (at your option) any later version.
17  *
18  * PAPPSOms++ is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with PAPPSOms++. If not, see <http://www.gnu.org/licenses/>.
25  *
26  ******************************************************************************/
27 
28 #include "timsdata.h"
29 #include "../../exception/exceptionnotfound.h"
30 #include "../../exception/exceptioninterrupted.h"
31 #include "../../processing/combiners/tracepluscombiner.h"
32 #include "../../processing/filters/filtertriangle.h"
33 #include "../../processing/filters/filterpass.h"
34 #include "../../processing/filters/filtersuitestring.h"
36 #include <QDebug>
37 #include <solvers.h>
38 #include <QSqlError>
39 #include <QSqlQuery>
40 #include <QSqlRecord>
41 #include <QMutexLocker>
42 #include <QThread>
43 #include <set>
44 #include <QtConcurrent>
45 
46 namespace pappso
47 {
48 
49 TimsData::TimsData(QDir timsDataDirectory)
50  : m_timsDataDirectory(timsDataDirectory)
51 {
52 
53  qDebug() << "Start of construction of TimsData";
55  if(!m_timsDataDirectory.exists())
56  {
57  throw PappsoException(
58  QObject::tr("ERROR TIMS data directory %1 not found")
59  .arg(m_timsDataDirectory.absolutePath()));
60  }
61 
62  if(!QFileInfo(m_timsDataDirectory.absoluteFilePath("analysis.tdf")).exists())
63  {
64 
65  throw PappsoException(
66  QObject::tr("ERROR TIMS data directory, %1 sqlite file not found")
67  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf")));
68  }
69 
70  // Open the database
71  QSqlDatabase qdb = openDatabaseConnection();
72 
73 
74  QSqlQuery q(qdb);
75  if(!q.exec("select Key, Value from GlobalMetadata where "
76  "Key='TimsCompressionType';"))
77  {
78 
79  qDebug();
80  throw PappsoException(
81  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
82  "command %2:\n%3\n%4\n%5")
83  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
84  .arg(q.lastQuery())
85  .arg(q.lastError().databaseText())
86  .arg(q.lastError().driverText())
87  .arg(q.lastError().nativeErrorCode()));
88  }
89 
90 
91  int compression_type = 0;
92  if(q.next())
93  {
94  compression_type = q.value(1).toInt();
95  }
96  qDebug() << " compression_type=" << compression_type;
98  QFileInfo(m_timsDataDirectory.absoluteFilePath("analysis.tdf_bin")),
99  compression_type);
100 
101  qDebug();
102 
103  // get number of precursors
105  if(!q.exec("SELECT COUNT( DISTINCT Id) FROM Precursors;"))
106  {
107  m_hasPrecursorTable = false;
108  }
109  else
110  {
111  m_hasPrecursorTable = true;
112  if(q.next())
113  {
114  m_totalNumberOfPrecursors = q.value(0).toLongLong();
115  }
116  }
117 
118 
120 
121  // get number of scans
122  if(!q.exec("SELECT SUM(NumScans),COUNT(Id) FROM Frames"))
123  {
124  qDebug();
125  throw PappsoException(
126  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
127  "command %2:\n%3\n%4\n%5")
128  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
129  .arg(q.lastQuery())
130  .arg(qdb.lastError().databaseText())
131  .arg(qdb.lastError().driverText())
132  .arg(qdb.lastError().nativeErrorCode()));
133  }
134  if(q.next())
135  {
136  m_totalNumberOfScans = q.value(0).toLongLong();
137  m_totalNumberOfFrames = q.value(1).toLongLong();
138  }
139 
140  if(!q.exec("select * from MzCalibration;"))
141  {
142  qDebug();
143  throw PappsoException(
144  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
145  "command %2:\n%3\n%4\n%5")
146  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
147  .arg(q.lastQuery())
148  .arg(q.lastError().databaseText())
149  .arg(q.lastError().driverText())
150  .arg(q.lastError().nativeErrorCode()));
151  }
152 
153  while(q.next())
154  {
155  QSqlRecord record = q.record();
157  std::pair<int, QSqlRecord>(record.value(0).toInt(), record));
158  }
159 
160  // m_mapTimsCalibrationRecord
161 
162  if(!q.exec("select * from TimsCalibration;"))
163  {
164  qDebug();
165  throw PappsoException(
166  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
167  "command %2:\n%3\n%4\n%5")
168  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
169  .arg(q.lastQuery())
170  .arg(q.lastError().databaseText())
171  .arg(q.lastError().driverText())
172  .arg(q.lastError().nativeErrorCode()));
173  }
174  while(q.next())
175  {
176  QSqlRecord record = q.record();
178  std::pair<int, QSqlRecord>(record.value(0).toInt(), record));
179  }
180 
181 
182  // store frames
183  if(!q.exec("select Frames.TimsId, Frames.AccumulationTime, " // 1
184  "Frames.MzCalibration, " // 2
185  "Frames.T1, Frames.T2, " // 4
186  "Frames.Time, Frames.MsMsType, Frames.TimsCalibration, " // 7
187  "Frames.Id " // 8
188  " FROM Frames;"))
189  {
190  qDebug();
191  throw PappsoException(
192  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
193  "command %2:\n%3\n%4\n%5")
194  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
195  .arg(q.lastQuery())
196  .arg(q.lastError().databaseText())
197  .arg(q.lastError().driverText())
198  .arg(q.lastError().nativeErrorCode()));
199  }
200 
202  while(q.next())
203  {
204  QSqlRecord record = q.record();
205  TimsFrameRecord &frame_record =
206  m_mapFramesRecord[record.value(8).toULongLong()];
207 
208  frame_record.tims_offset = record.value(0).toULongLong();
209  frame_record.accumulation_time = record.value(1).toDouble();
210  frame_record.mz_calibration_id = record.value(2).toULongLong();
211  frame_record.frame_t1 = record.value(3).toDouble();
212  frame_record.frame_t2 = record.value(4).toDouble();
213  frame_record.frame_time = record.value(5).toDouble();
214  frame_record.msms_type = record.value(6).toInt();
215  frame_record.tims_calibration_id = record.value(7).toULongLong();
216  }
217 
218  mcsp_ms2Filter = std::make_shared<pappso::FilterSuiteString>(
219  "chargeDeconvolution|0.02dalton mzExclusion|0.01dalton");
220 
221 
222  std::shared_ptr<FilterTriangle> ms1filter =
223  std::make_shared<FilterTriangle>();
224  ms1filter.get()->setTriangleSlope(50, 0.01);
225  mcsp_ms1Filter = ms1filter;
226  qDebug();
227 }
228 
229 void
230 TimsData::setMonoThread(bool is_mono_thread)
231 {
232  m_isMonoThread = is_mono_thread;
233 }
234 
235 QSqlDatabase
237 {
238  QString database_connection_name = QString("%1_%2")
239  .arg(m_timsDataDirectory.absolutePath())
240  .arg((quintptr)QThread::currentThread());
241  // Open the database
242  QSqlDatabase qdb = QSqlDatabase::database(database_connection_name);
243  if(!qdb.isValid())
244  {
245  qDebug() << database_connection_name;
246  qdb = QSqlDatabase::addDatabase("QSQLITE", database_connection_name);
247  qdb.setDatabaseName(m_timsDataDirectory.absoluteFilePath("analysis.tdf"));
248  }
249 
250 
251  if(!qdb.open())
252  {
253  qDebug();
254  throw PappsoException(
255  QObject::tr("ERROR opening TIMS sqlite database file %1, database name "
256  "%2 :\n%3\n%4\n%5")
257  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
258  .arg(database_connection_name)
259  .arg(qdb.lastError().databaseText())
260  .arg(qdb.lastError().driverText())
261  .arg(qdb.lastError().nativeErrorCode()));
262  }
263  return qdb;
264 }
265 
266 TimsData::TimsData([[maybe_unused]] const pappso::TimsData &other)
267 {
268  qDebug();
269 }
270 
272 {
273  // m_qdb.close();
274  if(mpa_timsBinDec != nullptr)
275  {
276  delete mpa_timsBinDec;
277  }
278  if(mpa_mzCalibrationStore != nullptr)
279  {
280  delete mpa_mzCalibrationStore;
281  }
282 }
283 
284 void
286 {
287  m_builtinMs2Centroid = centroid;
288 }
289 
290 bool
292 {
293  return m_builtinMs2Centroid;
294 }
295 
296 void
298 {
299  qDebug();
300  QSqlDatabase qdb = openDatabaseConnection();
301 
302  QSqlQuery q =
303  qdb.exec(QString("SELECT Id, NumScans FROM "
304  "Frames ORDER BY Id"));
305  if(q.lastError().isValid())
306  {
307 
308  throw PappsoException(
309  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
310  "command %2:\n%3\n%4\n%5")
311  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
312  .arg(q.lastQuery())
313  .arg(qdb.lastError().databaseText())
314  .arg(qdb.lastError().driverText())
315  .arg(qdb.lastError().nativeErrorCode()));
316  }
317  TimsFrameSPtr tims_frame;
318  bool index_found = false;
319  std::size_t timsId;
320  /** @brief number of scans in mobility dimension (number of TOF scans)
321  */
322  std::size_t numberScans;
323  std::size_t cumulScans = 0;
324  while(q.next() && (!index_found))
325  {
326  timsId = q.value(0).toULongLong();
327  numberScans = q.value(1).toULongLong();
328 
329  // qDebug() << timsId;
330 
332  std::pair<std::size_t, std::size_t>((cumulScans / 1000),
333  m_frameIdDescrList.size()));
334 
335  m_frameIdDescrList.push_back({timsId, numberScans, cumulScans});
336  cumulScans += numberScans;
337  }
338  qDebug();
339 }
340 
341 std::pair<std::size_t, std::size_t>
342 TimsData::getScanCoordinateFromRawIndex(std::size_t raw_index) const
343 {
344 
345  std::size_t fast_access = raw_index / 1000;
346  qDebug() << " fast_access=" << fast_access;
347  auto map_it = m_thousandIndexToFrameIdDescrListIndex.find(fast_access);
348  if(map_it == m_thousandIndexToFrameIdDescrListIndex.end())
349  {
350  throw ExceptionNotFound(
351  QObject::tr("ERROR raw index %1 not found (fast_access)")
352  .arg(raw_index));
353  }
354  std::size_t start_point_index = map_it->second;
355  while((start_point_index > 0) &&
356  (m_frameIdDescrList[start_point_index].m_cumulSize > raw_index))
357  {
358  start_point_index--;
359  }
360  for(std::size_t i = start_point_index; i < m_frameIdDescrList.size(); i++)
361  {
362 
363  if(raw_index <
364  (m_frameIdDescrList[i].m_cumulSize + m_frameIdDescrList[i].m_size))
365  {
366  return std::pair<std::size_t, std::size_t>(
367  m_frameIdDescrList[i].m_frameId,
368  raw_index - m_frameIdDescrList[i].m_cumulSize);
369  }
370  }
371 
372  throw ExceptionNotFound(
373  QObject::tr("ERROR raw index %1 not found").arg(raw_index));
374 }
375 
376 
377 std::size_t
379  std::size_t scan_num) const
380 {
381 
382  for(auto frameDescr : m_frameIdDescrList)
383  {
384  if(frameDescr.m_frameId == frame_id)
385  {
386  return frameDescr.m_cumulSize + scan_num;
387  }
388  }
389 
390  throw ExceptionNotFound(
391  QObject::tr("ERROR raw index with frame=%1 scan=%2 not found")
392  .arg(frame_id)
393  .arg(scan_num));
394 }
395 
396 /** @brief get a mass spectrum given its spectrum index
397  * @param raw_index a number begining at 0, corresponding to a Tims Scan in
398  * the order they lies in the binary data file
399  */
402 {
403 
404  qDebug() << " raw_index=" << raw_index;
405  try
406  {
407  auto coordinate = getScanCoordinateFromRawIndex(raw_index);
408  return getMassSpectrumCstSPtr(coordinate.first, coordinate.second);
409  }
410  catch(PappsoException &error)
411  {
412  throw PappsoException(
413  QObject::tr("Error TimsData::getMassSpectrumCstSPtrByRawIndex "
414  "raw_index=%1 :\n%2")
415  .arg(raw_index)
416  .arg(error.qwhat()));
417  }
418 }
419 
420 
423 {
424 
425  qDebug() << " timsId=" << timsId;
426 
427  const TimsFrameRecord &frame_record = m_mapFramesRecord[timsId];
428  if(timsId > m_totalNumberOfScans)
429  {
430  throw ExceptionNotFound(
431  QObject::tr("ERROR Frames database id %1 not found").arg(timsId));
432  }
433  TimsFrameBaseSPtr tims_frame;
434 
435 
436  tims_frame = std::make_shared<TimsFrameBase>(
437  TimsFrameBase(timsId, frame_record.tims_offset));
438 
439  auto it_map_record =
440  m_mapMzCalibrationRecord.find(frame_record.mz_calibration_id);
441  if(it_map_record != m_mapMzCalibrationRecord.end())
442  {
443 
444  double T1_frame = frame_record.frame_t1; // Frames.T1
445  double T2_frame = frame_record.frame_t2; // Frames.T2
446 
447 
448  tims_frame.get()->setMzCalibrationInterfaceSPtr(
450  T1_frame, T2_frame, it_map_record->second));
451  }
452  else
453  {
454  throw ExceptionNotFound(
455  QObject::tr("ERROR MzCalibration database id %1 not found")
456  .arg(frame_record.mz_calibration_id));
457  }
458 
459  tims_frame.get()->setAccumulationTime(frame_record.accumulation_time);
460 
461  tims_frame.get()->setTime(frame_record.frame_time);
462  tims_frame.get()->setMsMsType(frame_record.msms_type);
463 
464 
465  auto it_map_record_tims_calibration =
467  if(it_map_record_tims_calibration != m_mapTimsCalibrationRecord.end())
468  {
469 
470  tims_frame.get()->setTimsCalibration(
471  it_map_record_tims_calibration->second.value(1).toInt(),
472  it_map_record_tims_calibration->second.value(2).toDouble(),
473  it_map_record_tims_calibration->second.value(3).toDouble(),
474  it_map_record_tims_calibration->second.value(4).toDouble(),
475  it_map_record_tims_calibration->second.value(5).toDouble(),
476  it_map_record_tims_calibration->second.value(6).toDouble(),
477  it_map_record_tims_calibration->second.value(7).toDouble(),
478  it_map_record_tims_calibration->second.value(8).toDouble(),
479  it_map_record_tims_calibration->second.value(9).toDouble(),
480  it_map_record_tims_calibration->second.value(10).toDouble(),
481  it_map_record_tims_calibration->second.value(11).toDouble());
482  }
483  else
484  {
485  throw ExceptionNotFound(
486  QObject::tr("ERROR TimsCalibration database id %1 not found")
487  .arg(frame_record.tims_calibration_id));
488  }
489 
490  return tims_frame;
491 }
492 
493 std::vector<std::size_t>
494 TimsData::getTimsMS1FrameIdRange(double rt_begin, double rt_end) const
495 {
496 
497  qDebug() << " rt_begin=" << rt_begin << " rt_end=" << rt_end;
498  if(rt_begin < 0)
499  rt_begin = 0;
500  std::vector<std::size_t> tims_frameid_list;
501  QSqlDatabase qdb = openDatabaseConnection();
502  QSqlQuery q = qdb.exec(QString("SELECT Frames.Id FROM Frames WHERE "
503  "Frames.MsMsType=0 AND (Frames.Time>=%1) AND "
504  "(Frames.Time<=%2) ORDER BY Frames.Time;")
505  .arg(rt_begin)
506  .arg(rt_end));
507  if(q.lastError().isValid())
508  {
509 
510  throw PappsoException(
511  QObject::tr("ERROR in TIMS sqlite database file %1, database name %2, "
512  "executing SQL "
513  "command %3:\n%4\n%5\n%6")
514  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
515  .arg(qdb.databaseName())
516  .arg(q.lastQuery())
517  .arg(qdb.lastError().databaseText())
518  .arg(qdb.lastError().driverText())
519  .arg(qdb.lastError().nativeErrorCode()));
520  }
521  while(q.next())
522  {
523 
524  tims_frameid_list.push_back(q.value(0).toULongLong());
525  }
526  return tims_frameid_list;
527 }
528 
530 TimsData::getTimsFrameCstSPtr(std::size_t timsId)
531 {
532 
533  qDebug() << " timsId=" << timsId
534  << " m_mapFramesRecord.size()=" << m_mapFramesRecord.size();
535 
536  /*
537  for(auto pair_i : m_mapFramesRecord)
538  {
539 
540  qDebug() << " pair_i=" << pair_i.first;
541  }
542  */
543 
544  const TimsFrameRecord &frame_record = m_mapFramesRecord[timsId];
545  if(timsId > m_totalNumberOfScans)
546  {
547  throw ExceptionNotFound(
548  QObject::tr("ERROR Frames database id %1 not found").arg(timsId));
549  }
550 
551  TimsFrameSPtr tims_frame;
552 
553 
554  // QMutexLocker lock(&m_mutex);
555  tims_frame =
557  // lock.unlock();
558 
559  qDebug();
560  auto it_map_record =
561  m_mapMzCalibrationRecord.find(frame_record.mz_calibration_id);
562  if(it_map_record != m_mapMzCalibrationRecord.end())
563  {
564 
565  double T1_frame = frame_record.frame_t1; // Frames.T1
566  double T2_frame = frame_record.frame_t2; // Frames.T2
567 
568 
569  tims_frame.get()->setMzCalibrationInterfaceSPtr(
571  T1_frame, T2_frame, it_map_record->second));
572  }
573  else
574  {
575  throw ExceptionNotFound(
576  QObject::tr("ERROR MzCalibration database id %1 not found")
577  .arg(frame_record.mz_calibration_id));
578  }
579 
580  tims_frame.get()->setAccumulationTime(frame_record.accumulation_time);
581 
582  tims_frame.get()->setTime(frame_record.frame_time);
583  tims_frame.get()->setMsMsType(frame_record.msms_type);
584 
585  qDebug();
586  auto it_map_record_tims_calibration =
588  if(it_map_record_tims_calibration != m_mapTimsCalibrationRecord.end())
589  {
590 
591  tims_frame.get()->setTimsCalibration(
592  it_map_record_tims_calibration->second.value(1).toInt(),
593  it_map_record_tims_calibration->second.value(2).toDouble(),
594  it_map_record_tims_calibration->second.value(3).toDouble(),
595  it_map_record_tims_calibration->second.value(4).toDouble(),
596  it_map_record_tims_calibration->second.value(5).toDouble(),
597  it_map_record_tims_calibration->second.value(6).toDouble(),
598  it_map_record_tims_calibration->second.value(7).toDouble(),
599  it_map_record_tims_calibration->second.value(8).toDouble(),
600  it_map_record_tims_calibration->second.value(9).toDouble(),
601  it_map_record_tims_calibration->second.value(10).toDouble(),
602  it_map_record_tims_calibration->second.value(11).toDouble());
603  }
604  else
605  {
606  throw ExceptionNotFound(
607  QObject::tr("ERROR TimsCalibration database id %1 not found")
608  .arg(frame_record.tims_calibration_id));
609  }
610  qDebug();
611  return tims_frame;
612 }
613 
614 
616 TimsData::getMassSpectrumCstSPtr(std::size_t timsId, std::size_t scanNum)
617 {
618  qDebug() << " timsId=" << timsId << " scanNum=" << scanNum;
620 
621  return frame->getMassSpectrumCstSPtr(scanNum);
622 }
623 
624 std::size_t
626 {
627  return m_totalNumberOfScans;
628 }
629 
630 
631 std::size_t
633 {
635 }
636 
637 std::vector<std::size_t>
639  double mz_val,
640  double rt_sec,
641  double k0)
642 {
643  std::vector<std::size_t> precursor_ids;
644  std::vector<std::vector<double>> ids;
645 
646  QSqlDatabase qdb = openDatabaseConnection();
647  QSqlQuery q = qdb.exec(
648  QString(
649  "SELECT Frames.Time, Precursors.MonoisotopicMz, Precursors.Charge, "
650  "Precursors.Id, Frames.Id, PasefFrameMsMsInfo.ScanNumBegin, "
651  "PasefFrameMsMsInfo.scanNumEnd "
652  "FROM Frames "
653  "INNER JOIN PasefFrameMsMsInfo ON Frames.Id = PasefFrameMsMsInfo.Frame "
654  "INNER JOIN Precursors ON PasefFrameMsMsInfo.Precursor = Precursors.Id "
655  "WHERE Precursors.Charge == %1 "
656  "AND Precursors.MonoisotopicMz > %2 -0.01 "
657  "AND Precursors.MonoisotopicMz < %2 +0.01 "
658  "AND Frames.Time >= %3 -1 "
659  "AND Frames.Time < %3 +1; ")
660  .arg(charge)
661  .arg(mz_val)
662  .arg(rt_sec));
663  if(q.lastError().isValid())
664  {
665 
666  throw PappsoException(
667  QObject::tr("ERROR in TIMS sqlite database file %1, database name %2, "
668  "executing SQL "
669  "command %3:\n%4\n%5\n%6")
670  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
671  .arg(qdb.databaseName())
672  .arg(q.lastQuery())
673  .arg(qdb.lastError().databaseText())
674  .arg(qdb.lastError().driverText())
675  .arg(qdb.lastError().nativeErrorCode()));
676  }
677  while(q.next())
678  {
679  // qInfo() << q.value(0).toDouble() << q.value(1).toDouble()
680  // << q.value(2).toDouble() << q.value(3).toDouble();
681 
682  std::vector<double> sql_values;
683  sql_values.push_back(q.value(4).toDouble()); // frame id
684  sql_values.push_back(q.value(3).toDouble()); // precursor id
685  sql_values.push_back(q.value(5).toDouble()); // scan num begin
686  sql_values.push_back(q.value(6).toDouble()); // scan num end
687  sql_values.push_back(q.value(1).toDouble()); // mz_value
688 
689  ids.push_back(sql_values);
690 
691 
692  if(std::find(precursor_ids.begin(),
693  precursor_ids.end(),
694  q.value(3).toDouble()) == precursor_ids.end())
695  {
696  precursor_ids.push_back(q.value(3).toDouble());
697  }
698  }
699 
700  if(precursor_ids.size() > 1)
701  {
702  // std::vector<std::size_t> precursor_ids_ko =
703  // getMatchPrecursorIdByKo(ids, values[3]);
704  if(precursor_ids.size() > 1)
705  {
706  precursor_ids = getClosestPrecursorIdByMz(ids, k0);
707  }
708  return precursor_ids;
709  }
710  else
711  {
712  return precursor_ids;
713  }
714 }
715 
716 std::vector<std::size_t>
717 TimsData::getMatchPrecursorIdByKo(std::vector<std::vector<double>> ids,
718  double ko_value)
719 {
720  std::vector<std::size_t> precursor_id;
721  for(std::vector<double> index : ids)
722  {
723  auto coordinate = getScanCoordinateFromRawIndex(index[0]);
724 
725  TimsFrameBaseCstSPtr tims_frame;
726  tims_frame = getTimsFrameBaseCstSPtrCached(coordinate.first);
727 
728  double bko = tims_frame.get()->getOneOverK0Transformation(index[2]);
729  double eko = tims_frame.get()->getOneOverK0Transformation(index[3]);
730 
731  // qInfo() << "diff" << (bko + eko) / 2;
732  double mean_ko = (bko + eko) / 2;
733 
734  if(mean_ko > ko_value - 0.1 && mean_ko < ko_value + 0.1)
735  {
736  precursor_id.push_back(index[1]);
737  }
738  }
739  return precursor_id;
740 }
741 
742 std::vector<std::size_t>
743 TimsData::getClosestPrecursorIdByMz(std::vector<std::vector<double>> ids,
744  double mz_value)
745 {
746  std::vector<std::size_t> best_precursor;
747  double best_value = 1;
748  int count = 1;
749  int best_val_position = 0;
750 
751  for(std::vector<double> values : ids)
752  {
753  double new_val = abs(mz_value - values[4]);
754  if(new_val < best_value)
755  {
756  best_value = new_val;
757  best_val_position = count;
758  }
759  count++;
760  }
761  best_precursor.push_back(ids[best_val_position][1]);
762  return best_precursor;
763 }
764 
765 
766 unsigned int
767 TimsData::getMsLevelBySpectrumIndex(std::size_t spectrum_index)
768 {
769  auto coordinate = getScanCoordinateFromRawIndex(spectrum_index);
770  auto tims_frame = getTimsFrameCstSPtrCached(coordinate.first);
771  return tims_frame.get()->getMsLevel();
772 }
773 
774 
775 void
777  const MsRunIdCstSPtr &msrun_id,
778  QualifiedMassSpectrum &mass_spectrum,
779  std::size_t spectrum_index,
780  bool want_binary_data)
781 {
782  try
783  {
784  auto coordinate = getScanCoordinateFromRawIndex(spectrum_index);
785  TimsFrameBaseCstSPtr tims_frame;
786  if(want_binary_data)
787  {
788  tims_frame = getTimsFrameCstSPtrCached(coordinate.first);
789  }
790  else
791  {
792  tims_frame = getTimsFrameBaseCstSPtrCached(coordinate.first);
793  }
794  MassSpectrumId spectrum_id;
795 
796  spectrum_id.setSpectrumIndex(spectrum_index);
797  spectrum_id.setMsRunId(msrun_id);
798  spectrum_id.setNativeId(QString("frame=%1 scan=%2 index=%3")
799  .arg(coordinate.first)
800  .arg(coordinate.second)
801  .arg(spectrum_index));
802 
803  mass_spectrum.setMassSpectrumId(spectrum_id);
804 
805  mass_spectrum.setMsLevel(tims_frame.get()->getMsLevel());
806  mass_spectrum.setRtInSeconds(tims_frame.get()->getTime());
807 
808  mass_spectrum.setDtInMilliSeconds(
809  tims_frame.get()->getDriftTime(coordinate.second));
810  // 1/K0
811  mass_spectrum.setParameterValue(
813  tims_frame.get()->getOneOverK0Transformation(coordinate.second));
814 
815  mass_spectrum.setEmptyMassSpectrum(true);
816  if(want_binary_data)
817  {
818  mass_spectrum.setMassSpectrumSPtr(
819  tims_frame.get()->getMassSpectrumSPtr(coordinate.second));
820  if(mass_spectrum.size() > 0)
821  {
822  mass_spectrum.setEmptyMassSpectrum(false);
823  }
824  }
825  else
826  {
827  // if(tims_frame.get()->getNbrPeaks(coordinate.second) > 0)
828  //{
829  mass_spectrum.setEmptyMassSpectrum(false);
830  // }
831  }
832  if(tims_frame.get()->getMsLevel() > 1)
833  {
834 
835  auto spectrum_descr = getSpectrumDescrWithScanCoordinate(coordinate);
836  if(spectrum_descr.precursor_id > 0)
837  {
838 
839  mass_spectrum.appendPrecursorIonData(
840  spectrum_descr.precursor_ion_data);
841 
842 
843  MassSpectrumId spectrum_id;
844  std::size_t prec_spectrum_index = getRawIndexFromCoordinate(
845  spectrum_descr.parent_frame, coordinate.second);
846 
847  mass_spectrum.setPrecursorSpectrumIndex(prec_spectrum_index);
848  mass_spectrum.setPrecursorNativeId(
849  QString("frame=%1 scan=%2 index=%3")
850  .arg(spectrum_descr.parent_frame)
851  .arg(coordinate.second)
852  .arg(prec_spectrum_index));
853 
854  mass_spectrum.setParameterValue(
856  spectrum_descr.isolationMz);
857  mass_spectrum.setParameterValue(
859  spectrum_descr.isolationWidth);
860 
861  mass_spectrum.setParameterValue(
863  spectrum_descr.collisionEnergy);
864  mass_spectrum.setParameterValue(
866  (quint64)spectrum_descr.precursor_id);
867  }
868  }
869  }
870  catch(PappsoException &error)
871  {
873  QObject::tr("Error TimsData::getQualifiedMassSpectrumByRawIndex "
874  "spectrum_index=%1 :\n%2")
875  .arg(spectrum_index)
876  .arg(error.qwhat()));
877  }
878 }
879 
880 
881 Trace
883 {
884  // In the Frames table, each frame has a record describing the
885  // SummedIntensities for all the mobility spectra.
886 
887 
888  MapTrace rt_tic_map_trace;
889 
890  using Pair = std::pair<double, double>;
891  using Map = std::map<double, double>;
892  using Iterator = Map::iterator;
893 
894 
895  QSqlDatabase qdb = openDatabaseConnection();
896  QSqlQuery q =
897  qdb.exec(QString("SELECT Time, SummedIntensities "
898  "FROM Frames WHERE MsMsType = 0 "
899  "ORDER BY Time;"));
900 
901  if(q.lastError().isValid())
902  {
903 
904  throw PappsoException(
905  QObject::tr("ERROR in TIMS sqlite database file %1, database name %2, "
906  "executing SQL "
907  "command %3:\n%4\n%5\n%6")
908  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
909  .arg(qdb.databaseName())
910  .arg(q.lastQuery())
911  .arg(qdb.lastError().databaseText())
912  .arg(qdb.lastError().driverText())
913  .arg(qdb.lastError().nativeErrorCode()));
914  }
915 
916  while(q.next())
917  {
918 
919  bool ok = false;
920 
921  int cumulated_results = 2;
922 
923  double rt = q.value(0).toDouble(&ok);
924  cumulated_results -= ok;
925 
926  double sumY = q.value(1).toDouble(&ok);
927  cumulated_results -= ok;
928 
929  if(cumulated_results)
930  {
931  throw PappsoException(
932  QObject::tr(
933  "ERROR in TIMS sqlite database file: could not read either the "
934  "retention time or the summed intensities (%1, database name "
935  "%2, "
936  "executing SQL "
937  "command %3:\n%4\n%5\n%6")
938  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
939  .arg(qdb.databaseName())
940  .arg(q.lastQuery())
941  .arg(qdb.lastError().databaseText())
942  .arg(qdb.lastError().driverText())
943  .arg(qdb.lastError().nativeErrorCode()));
944  }
945 
946  // Try to insert value sumY at key rt.
947  std::pair<Iterator, bool> res = rt_tic_map_trace.insert(Pair(rt, sumY));
948 
949  if(!res.second)
950  {
951  // One other same rt value was seen already (like in ion mobility
952  // mass spectrometry, for example). Only increment the y value.
953 
954  res.first->second += sumY;
955  }
956  }
957 
958  // qDebug().noquote() << "The TIC chromatogram:\n"
959  //<< rt_tic_map_trace.toTrace().toString();
960 
961  return rt_tic_map_trace.toTrace();
962 }
963 
964 
965 void
967  const MsRunIdCstSPtr &msrun_id,
968  QualifiedMassSpectrum &mass_spectrum,
969  SpectrumDescr &spectrum_descr,
970  bool want_binary_data)
971 {
972 
973  qDebug() << " ms2_index=" << spectrum_descr.ms2_index
974  << " precursor_index=" << spectrum_descr.precursor_id;
975 
976  TracePlusCombiner combiner;
977  MapTrace combiner_result;
978 
979  try
980  {
981  mass_spectrum.setMsLevel(1);
982  mass_spectrum.setPrecursorSpectrumIndex(0);
983  mass_spectrum.setEmptyMassSpectrum(true);
984 
985  MassSpectrumId spectrum_id;
986  spectrum_id.setSpectrumIndex(spectrum_descr.ms1_index);
987  spectrum_id.setNativeId(
988  QString("frame=%1 begin=%2 end=%3 precursor=%4 idxms1=%5")
989  .arg(spectrum_descr.parent_frame)
990  .arg(spectrum_descr.scan_mobility_start)
991  .arg(spectrum_descr.scan_mobility_end)
992  .arg(spectrum_descr.precursor_id)
993  .arg(spectrum_descr.ms1_index));
994 
995  spectrum_id.setMsRunId(msrun_id);
996 
997  mass_spectrum.setMassSpectrumId(spectrum_id);
998 
999 
1000  TimsFrameBaseCstSPtr tims_frame;
1001  if(want_binary_data)
1002  {
1003  qDebug() << "bindec";
1004  tims_frame = getTimsFrameCstSPtrCached(spectrum_descr.parent_frame);
1005  }
1006  else
1007  {
1008  tims_frame =
1010  }
1011  mass_spectrum.setRtInSeconds(tims_frame.get()->getTime());
1012 
1013  mass_spectrum.setParameterValue(
1015  tims_frame.get()->getOneOverK0Transformation(
1016  spectrum_descr.scan_mobility_start));
1017 
1018  mass_spectrum.setParameterValue(
1020  tims_frame.get()->getOneOverK0Transformation(
1021  spectrum_descr.scan_mobility_end));
1022 
1023 
1024  if(want_binary_data)
1025  {
1026  combiner.combine(combiner_result,
1027  tims_frame.get()->cumulateScanToTrace(
1028  spectrum_descr.scan_mobility_start,
1029  spectrum_descr.scan_mobility_end));
1030 
1031  pappso::Trace trace(combiner_result);
1032  qDebug();
1033 
1034  if(trace.size() > 0)
1035  {
1036  if(mcsp_ms1Filter != nullptr)
1037  {
1038  mcsp_ms1Filter->filter(trace);
1039  }
1040 
1041  qDebug();
1042  mass_spectrum.setMassSpectrumSPtr(
1043  MassSpectrum(trace).makeMassSpectrumSPtr());
1044  mass_spectrum.setEmptyMassSpectrum(false);
1045  }
1046  else
1047  {
1048  mass_spectrum.setMassSpectrumSPtr(nullptr);
1049  mass_spectrum.setEmptyMassSpectrum(true);
1050  }
1051  }
1052  qDebug();
1053  }
1054 
1055  catch(PappsoException &error)
1056  {
1057  throw error;
1058  }
1059  catch(std::exception &error)
1060  {
1061  qDebug() << QString("Failure %1 ").arg(error.what());
1062  }
1063 }
1064 
1065 
1068 {
1069  QMutexLocker locker(&m_mutex);
1070  for(auto &tims_frame : m_timsFrameBaseCache)
1071  {
1072  if(tims_frame.get()->getId() == timsId)
1073  {
1074  m_timsFrameBaseCache.push_back(tims_frame);
1075  if(m_timsFrameBaseCache.size() > m_cacheSize)
1076  m_timsFrameBaseCache.pop_front();
1077  return tims_frame;
1078  }
1079  }
1080 
1081  m_timsFrameBaseCache.push_back(getTimsFrameBaseCstSPtr(timsId));
1082  if(m_timsFrameBaseCache.size() > m_cacheSize)
1083  m_timsFrameBaseCache.pop_front();
1084  return m_timsFrameBaseCache.back();
1085 }
1086 
1089 {
1090  qDebug();
1091  QMutexLocker locker(&m_mutex);
1092  for(auto &tims_frame : m_timsFrameCache)
1093  {
1094  if(tims_frame.get()->getId() == timsId)
1095  {
1096  m_timsFrameCache.push_back(tims_frame);
1097  if(m_timsFrameCache.size() > m_cacheSize)
1098  m_timsFrameCache.pop_front();
1099  return tims_frame;
1100  }
1101  }
1102  pappso::TimsFrameCstSPtr frame_sptr = getTimsFrameCstSPtr(timsId);
1103 
1104  // locker.relock();
1105  qDebug();
1106 
1107  m_timsFrameCache.push_back(frame_sptr);
1108  if(m_timsFrameCache.size() > m_cacheSize)
1109  m_timsFrameCache.pop_front();
1110  qDebug();
1111  return m_timsFrameCache.back();
1112 
1113 
1114  /*
1115 // the frame is not in the cache
1116 if(std::find(m_someoneIsLoadingFrameId.begin(),
1117  m_someoneIsLoadingFrameId.end(),
1118  timsId) == m_someoneIsLoadingFrameId.end())
1119  {
1120  // not found, we are alone on this frame
1121  m_someoneIsLoadingFrameId.push_back(timsId);
1122  qDebug();
1123  //locker.unlock();
1124  pappso::TimsFrameCstSPtr frame_sptr = getTimsFrameCstSPtr(timsId);
1125 
1126  // locker.relock();
1127  qDebug();
1128  m_someoneIsLoadingFrameId.erase(
1129  std::find(m_someoneIsLoadingFrameId.begin(),
1130  m_someoneIsLoadingFrameId.end(),
1131  timsId));
1132 
1133  m_timsFrameCache.push_back(frame_sptr);
1134  if(m_timsFrameCache.size() > m_cacheSize)
1135  m_timsFrameCache.pop_front();
1136  qDebug();
1137  return m_timsFrameCache.back();
1138  }
1139 else
1140  {
1141  // this frame is loading by someone else, we have to wait
1142  qDebug();
1143  // locker.unlock();
1144  // std::size_t another_frame_id = timsId;
1145  while(true)
1146  {
1147  QThread::usleep(1);
1148  // locker.relock();
1149 
1150  for(auto &tims_frame : m_timsFrameCache)
1151  {
1152  if(tims_frame.get()->getId() == timsId)
1153  {
1154  m_timsFrameCache.push_back(tims_frame);
1155  return tims_frame;
1156  }
1157  }
1158  // locker.unlock();
1159 }
1160 } // namespace pappso
1161 */
1162 }
1163 
1164 void
1166 {
1167  mcsp_ms2Filter = filter;
1168 }
1169 void
1171 {
1172  mcsp_ms1Filter = filter;
1173 }
1174 
1177  PrecisionPtr precision_ptr)
1178 {
1179 
1180  qDebug();
1181  XicCoordTims xic_coord_tims_struct;
1182 
1183  try
1184  {
1185  if(m_mapXicCoordRecord.size() == 0)
1186  {
1187  QMutexLocker lock(&m_mutex);
1188  // Go get records!
1189 
1190  // We proceed in this way:
1191 
1192  // 1. For each Precursor reference to the Precursors table's ID
1193  // found in the PasefFrameMsMsInfo table, store the precursor ID for
1194  // step 2.
1195 
1196  // 2. From the Precursors table's ID from step 1, get the
1197  // MonoisotopicMz.
1198 
1199  // 3. From the PasefFrameMsMsInfo table, for the Precursors table's
1200  // ID reference, get a reference to the Frames table's ID. Thanks to
1201  // the Frames ID, look for the Time value (acquisition retention
1202  // time) for the MS/MS spectrum. The Time value in the Frames tables
1203  // always corresponds to a Frame of MsMsType 8 (that is, MS/MS),
1204  // which is expected since we are looking into MS/MS data.
1205 
1206  // 4. From the PasefFrameMsMsInfo table, associate the values
1207  // ScanNumBegin and ScanNumEnd, the mobility bins in which the
1208  // precursor was found.
1209 
1210 
1211  QSqlDatabase qdb = openDatabaseConnection();
1212  QSqlQuery q = qdb.exec(
1213  QString("SELECT Precursors.id, "
1214  "min(Frames.Time), "
1215  "min(PasefFrameMsMsInfo.ScanNumBegin), "
1216  "max(PasefFrameMsMsInfo.ScanNumEnd), "
1217  "Precursors.MonoisotopicMz "
1218  "FROM "
1219  "PasefFrameMsMsInfo INNER JOIN Precursors ON "
1220  "PasefFrameMsMsInfo.Precursor=Precursors.Id INNER JOIN "
1221  "Frames ON PasefFrameMsMsInfo.Frame=Frames.Id "
1222  "GROUP BY Precursors.id;"));
1223  if(q.lastError().isValid())
1224  {
1225  qDebug();
1226  throw PappsoException(
1227  QObject::tr(
1228  "ERROR in TIMS sqlite database file %1, executing SQL "
1229  "command %2:\n%3\n%4\n%5")
1230  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
1231  .arg(q.lastQuery())
1232  .arg(qdb.lastError().databaseText())
1233  .arg(qdb.lastError().driverText())
1234  .arg(qdb.lastError().nativeErrorCode()));
1235  }
1236 
1237  q.last(); // strange bug : get the last sql record and get back,
1238  // otherwise it will not retrieve all records.
1239  q.first();
1240  // std::size_t i = 0;
1241  do
1242  {
1243  QSqlRecord record = q.record();
1244  m_mapXicCoordRecord.insert(std::pair<std::size_t, QSqlRecord>(
1245  (std::size_t)record.value(0).toULongLong(), record));
1246  }
1247  while(q.next());
1248  }
1249 
1250 
1251  auto it_map_xiccoord = m_mapXicCoordRecord.find(precursor_id);
1252  if(it_map_xiccoord == m_mapXicCoordRecord.end())
1253  {
1254 
1255  throw ExceptionNotFound(
1256  QObject::tr("ERROR Precursors database id %1 not found")
1257  .arg(precursor_id));
1258  }
1259 
1260  auto &q = it_map_xiccoord->second;
1261  xic_coord_tims_struct.mzRange =
1262  MzRange(q.value(4).toDouble(), precision_ptr);
1263  xic_coord_tims_struct.scanNumBegin = q.value(2).toUInt();
1264  xic_coord_tims_struct.scanNumEnd = q.value(3).toUInt();
1265  xic_coord_tims_struct.rtTarget = q.value(1).toDouble();
1266  // xic_structure.charge = q.value(5).toUInt();
1267  xic_coord_tims_struct.xicSptr = std::make_shared<Xic>();
1268  }
1269  catch(PappsoException &error)
1270  {
1271  throw error;
1272  }
1273  catch(std::exception &error)
1274  {
1275  qDebug() << QString("Failure %1 ").arg(error.what());
1276  }
1277  return xic_coord_tims_struct;
1278 }
1279 
1280 
1281 std::map<quint32, quint32>
1282 TimsData::getRawMs2ByPrecursorId(std::size_t precursor_index)
1283 {
1284  qDebug();
1285  std::map<quint32, quint32> raw_spectrum;
1286  try
1287  {
1288  QSqlDatabase qdb = openDatabaseConnection();
1289 
1290  qdb = openDatabaseConnection();
1291  QSqlQuery q =
1292  qdb.exec(QString("SELECT PasefFrameMsMsInfo.*, Precursors.* FROM "
1293  "PasefFrameMsMsInfo INNER JOIN Precursors ON "
1294  "PasefFrameMsMsInfo.Precursor=Precursors.Id where "
1295  "Precursors.Id=%1;")
1296  .arg(precursor_index));
1297  if(q.lastError().isValid())
1298  {
1299  qDebug();
1300  throw PappsoException(
1301  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
1302  "command %2:\n%3\n%4\n%5")
1303  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
1304  .arg(q.lastQuery())
1305  .arg(qdb.lastError().databaseText())
1306  .arg(qdb.lastError().driverText())
1307  .arg(qdb.lastError().nativeErrorCode()));
1308  }
1309  qDebug();
1310  // m_mutex.unlock();
1311  if(q.size() == 0)
1312  {
1313 
1314  throw ExceptionNotFound(
1315  QObject::tr(
1316  "ERROR in getQualifiedMassSpectrumByPrecursorId, precursor "
1317  "id=%1 not found")
1318  .arg(precursor_index));
1319  }
1320  else
1321  {
1322  // qDebug() << " q.size()="<< q.size();
1323  qDebug();
1324  bool first = true;
1325  std::size_t scan_mobility_start = 0;
1326  std::size_t scan_mobility_end = 0;
1327  std::vector<std::size_t> tims_frame_list;
1328 
1329  while(q.next())
1330  {
1331  tims_frame_list.push_back(q.value(0).toLongLong());
1332  if(first)
1333  {
1334 
1335  scan_mobility_start = q.value(1).toLongLong();
1336  scan_mobility_end = q.value(2).toLongLong();
1337 
1338  first = false;
1339  }
1340  }
1341  // QMutexLocker locker(&m_mutex_spectrum);
1342  qDebug();
1343  pappso::TimsFrameCstSPtr tims_frame, previous_frame;
1344  // TracePlusCombiner combiner;
1345  // MapTrace combiner_result;
1346  for(std::size_t tims_id : tims_frame_list)
1347  {
1348  tims_frame = getTimsFrameCstSPtrCached(tims_id);
1349  qDebug();
1350  /*combiner.combine(combiner_result,
1351  tims_frame.get()->cumulateScanToTrace(
1352  scan_mobility_start, scan_mobility_end));*/
1353  if(previous_frame.get() != nullptr)
1354  {
1355  if(previous_frame.get()->hasSameCalibrationData(
1356  *tims_frame.get()))
1357  {
1358  }
1359  else
1360  {
1361  throw ExceptionNotFound(
1362  QObject::tr(
1363  "ERROR in %1 %2, different calibration data "
1364  "between frame id %3 and frame id %4")
1365  .arg(__FILE__)
1366  .arg(__FUNCTION__)
1367  .arg(previous_frame.get()->getId())
1368  .arg(tims_frame.get()->getId()));
1369  }
1370  }
1371  tims_frame.get()->cumulateScansInRawMap(
1372  raw_spectrum, scan_mobility_start, scan_mobility_end);
1373  qDebug();
1374 
1375  previous_frame = tims_frame;
1376  }
1377  qDebug() << " precursor_index=" << precursor_index
1378  << " num_rows=" << tims_frame_list.size()
1379  << " sql=" << q.lastQuery() << " "
1380  << (std::size_t)QThread::currentThreadId();
1381  if(first == true)
1382  {
1383  throw ExceptionNotFound(
1384  QObject::tr(
1385  "ERROR in getQualifiedMassSpectrumByPrecursorId, precursor "
1386  "id=%1 not found")
1387  .arg(precursor_index));
1388  }
1389  qDebug();
1390  }
1391  }
1392 
1393  catch(PappsoException &error)
1394  {
1395  throw PappsoException(QObject::tr("ERROR in %1 (precursor_index=%2):\n%3")
1396  .arg(__FUNCTION__)
1397  .arg(precursor_index)
1398  .arg(error.qwhat()));
1399  }
1400  catch(std::exception &error)
1401  {
1402  qDebug() << QString("Failure %1 ").arg(error.what());
1403  }
1404  return raw_spectrum;
1405  qDebug();
1406 }
1407 
1408 
1409 void
1411  const MsRunIdCstSPtr &msrun_id,
1412  QualifiedMassSpectrum &mass_spectrum,
1413  SpectrumDescr &spectrum_descr,
1414  bool want_binary_data)
1415 {
1416  try
1417  {
1418  qDebug();
1419  MassSpectrumId spectrum_id;
1420 
1421  spectrum_id.setSpectrumIndex(spectrum_descr.ms2_index);
1422  spectrum_id.setNativeId(QString("precursor=%1 idxms2=%2")
1423  .arg(spectrum_descr.precursor_id)
1424  .arg(spectrum_descr.ms2_index));
1425  spectrum_id.setMsRunId(msrun_id);
1426 
1427  mass_spectrum.setMassSpectrumId(spectrum_id);
1428 
1429  mass_spectrum.setMsLevel(2);
1430  qDebug() << "spectrum_descr.precursor_id=" << spectrum_descr.precursor_id
1431  << " spectrum_descr.ms1_index=" << spectrum_descr.ms1_index
1432  << " spectrum_descr.ms2_index=" << spectrum_descr.ms2_index;
1433  mass_spectrum.setPrecursorSpectrumIndex(spectrum_descr.ms1_index);
1434 
1435  mass_spectrum.setEmptyMassSpectrum(true);
1436 
1437  qDebug();
1438 
1439 
1440  mass_spectrum.appendPrecursorIonData(spectrum_descr.precursor_ion_data);
1441 
1442  mass_spectrum.setPrecursorNativeId(
1443  QString("frame=%1 begin=%2 end=%3 precursor=%4 idxms1=%5")
1444  .arg(spectrum_descr.parent_frame)
1445  .arg(spectrum_descr.scan_mobility_start)
1446  .arg(spectrum_descr.scan_mobility_end)
1447  .arg(spectrum_descr.precursor_id)
1448  .arg(spectrum_descr.ms1_index));
1449 
1450  mass_spectrum.setParameterValue(
1452  spectrum_descr.isolationMz);
1453  mass_spectrum.setParameterValue(
1455  spectrum_descr.isolationWidth);
1456 
1457  mass_spectrum.setParameterValue(
1459  spectrum_descr.collisionEnergy);
1460  mass_spectrum.setParameterValue(
1462  (quint64)spectrum_descr.precursor_id);
1463 
1464  // QMutexLocker locker(&m_mutex_spectrum);
1465  qDebug();
1466  pappso::TimsFrameBaseCstSPtr tims_frame, previous_frame;
1467  // TracePlusCombiner combiner;
1468  // MapTrace combiner_result;
1469  std::map<quint32, quint32> raw_spectrum;
1470  bool first = true;
1471  for(std::size_t tims_id : spectrum_descr.tims_frame_list)
1472  {
1473  qDebug() << " precursor_index=" << spectrum_descr.precursor_id
1474  << " tims_id=" << tims_id
1475  << (std::size_t)QThread::currentThreadId();
1476  ;
1477  if(want_binary_data)
1478  {
1479  qDebug() << "bindec";
1480  tims_frame = getTimsFrameCstSPtrCached(tims_id);
1481  }
1482  else
1483  {
1484  tims_frame = getTimsFrameBaseCstSPtrCached(tims_id);
1485  }
1486  qDebug() << (std::size_t)QThread::currentThreadId();
1487  ;
1488  if(first)
1489  {
1490  mass_spectrum.setRtInSeconds(tims_frame.get()->getTime());
1491 
1492  mass_spectrum.setParameterValue(
1494  tims_frame.get()->getOneOverK0Transformation(
1495  spectrum_descr.scan_mobility_start));
1496 
1497  mass_spectrum.setParameterValue(
1499  tims_frame.get()->getOneOverK0Transformation(
1500  spectrum_descr.scan_mobility_end));
1501 
1502  first = false;
1503  }
1504 
1505 
1506  if(want_binary_data)
1507  {
1508  qDebug();
1509  /*combiner.combine(combiner_result,
1510  tims_frame.get()->cumulateScanToTrace(
1511  scan_mobility_start, scan_mobility_end));*/
1512  if(previous_frame.get() != nullptr)
1513  {
1514  if(previous_frame.get()->hasSameCalibrationData(
1515  *tims_frame.get()))
1516  {
1517  }
1518  else
1519  {
1520  throw ExceptionNotFound(
1521  QObject::tr(
1522  "ERROR in %1 %2, different calibration data "
1523  "between frame id %3 and frame id %4")
1524  .arg(__FILE__)
1525  .arg(__FUNCTION__)
1526  .arg(previous_frame.get()->getId())
1527  .arg(tims_frame.get()->getId()));
1528  }
1529  }
1530  qDebug() << (std::size_t)QThread::currentThreadId();
1531  ;
1532  tims_frame.get()->cumulateScansInRawMap(
1533  raw_spectrum,
1534  spectrum_descr.scan_mobility_start,
1535  spectrum_descr.scan_mobility_end);
1536  qDebug() << (std::size_t)QThread::currentThreadId();
1537  ;
1538  }
1539  previous_frame = tims_frame;
1540  }
1541  qDebug() << " precursor_index=" << spectrum_descr.precursor_id
1542  << " num_rows=" << spectrum_descr.tims_frame_list.size()
1543  << (std::size_t)QThread::currentThreadId();
1544  if(first == true)
1545  {
1546  throw ExceptionNotFound(
1547  QObject::tr(
1548  "ERROR in getQualifiedMassSpectrumByPrecursorId, precursor "
1549  "id=%1 not found")
1550  .arg(spectrum_descr.precursor_id));
1551  }
1552  if(want_binary_data)
1553  {
1554  qDebug() << " precursor_index=" << spectrum_descr.precursor_id;
1555  // peak_pick.filter(trace);
1556  pappso::Trace trace;
1558  {
1559  trace =
1560  tims_frame.get()->getTraceFromCumulatedScansBuiltinCentroid(
1561  raw_spectrum);
1562  }
1563  else
1564  {
1565  // no builtin centroid:
1566 
1567  trace =
1568  tims_frame.get()->getTraceFromCumulatedScans(raw_spectrum);
1569  }
1570 
1571  if(trace.size() > 0)
1572  {
1573  qDebug() << " precursor_index=" << spectrum_descr.precursor_id
1574  << " " << trace.size() << " "
1575  << (std::size_t)QThread::currentThreadId();
1576 
1577  if(mcsp_ms2Filter != nullptr)
1578  {
1579  // FilterTriangle filter;
1580  // filter.setTriangleSlope(50, 0.02);
1581  // filter.filter(trace);
1582  // trace.filter(pappso::FilterHighPass(10));
1583  mcsp_ms2Filter->filter(trace);
1584  }
1585 
1586  // FilterScaleFactorY filter_scale((double)1 /
1587  // (double)tims_frame_list.size());
1588  // filter_scale.filter(trace);
1589  qDebug() << " precursor_index=" << spectrum_descr.precursor_id;
1590  mass_spectrum.setMassSpectrumSPtr(
1591  MassSpectrum(trace).makeMassSpectrumSPtr());
1592  mass_spectrum.setEmptyMassSpectrum(false);
1593  }
1594  else
1595  {
1596  mass_spectrum.setMassSpectrumSPtr(nullptr);
1597  mass_spectrum.setEmptyMassSpectrum(true);
1598  }
1599 
1600  qDebug();
1601  }
1602  qDebug();
1603  }
1604 
1605  catch(PappsoException &error)
1606  {
1607  throw PappsoException(
1608  QObject::tr("ERROR in %1 (ms2_index=%2 precursor_index=%3):\n%4")
1609  .arg(__FUNCTION__)
1610  .arg(spectrum_descr.ms2_index)
1611  .arg(spectrum_descr.precursor_id)
1612  .arg(error.qwhat()));
1613  }
1614  catch(std::exception &error)
1615  {
1616  qDebug() << QString("Failure %1 ").arg(error.what());
1617  }
1618  qDebug();
1619 }
1620 
1621 void
1623  const MsRunIdCstSPtr &msrun_id,
1625  unsigned int ms_level)
1626 {
1627  qDebug() << " ms_level=" << ms_level;
1628  if(!m_hasPrecursorTable)
1629  {
1630  throw PappsoException(
1631  QObject::tr("unable to read spectrum list : this data file does not "
1632  "contain MS2 data, no precursor found."));
1633  }
1634 
1635  QSqlDatabase qdb = openDatabaseConnection();
1636  QSqlQuery qprecursor_list = qdb.exec(QString(
1637  "SELECT PasefFrameMsMsInfo.Frame, " // 0
1638  "PasefFrameMsMsInfo.ScanNumBegin, " // 1
1639  "PasefFrameMsMsInfo.ScanNumEnd, " // 2
1640  "PasefFrameMsMsInfo.IsolationMz, " // 3
1641  "PasefFrameMsMsInfo.IsolationWidth, " // 4
1642  "PasefFrameMsMsInfo.CollisionEnergy, " // 5
1643  "PasefFrameMsMsInfo.Precursor, " // 6
1644  "Precursors.Id, " // 7
1645  "Precursors.LargestPeakMz, " // 8
1646  "Precursors.AverageMz, " // 9
1647  "Precursors.MonoisotopicMz, " // 10
1648  "Precursors.Charge, " // 11
1649  "Precursors.ScanNumber, " // 12
1650  "Precursors.Intensity, " // 13
1651  "Precursors.Parent " // 14
1652  "FROM PasefFrameMsMsInfo "
1653  "INNER JOIN Precursors ON "
1654  "PasefFrameMsMsInfo.Precursor=Precursors.Id "
1655  "ORDER BY PasefFrameMsMsInfo.Precursor, PasefFrameMsMsInfo.Frame ;"));
1656  if(qprecursor_list.lastError().isValid())
1657  {
1658 
1659  throw PappsoException(
1660  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
1661  "command %2:\n%3\n%4\n%5")
1662  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
1663  .arg(qprecursor_list.lastQuery())
1664  .arg(qdb.lastError().databaseText())
1665  .arg(qdb.lastError().driverText())
1666  .arg(qdb.lastError().nativeErrorCode()));
1667  }
1668 
1669 
1670  qDebug() << "qprecursor_list.size()=" << qprecursor_list.size();
1671  qDebug() << QObject::tr(
1672  "TIMS sqlite database file %1, executing SQL "
1673  "command %2:\n%3\n%4\n%5")
1674  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
1675  .arg(qprecursor_list.lastQuery())
1676  .arg(qdb.lastError().databaseText())
1677  .arg(qdb.lastError().driverText())
1678  .arg(qdb.lastError().nativeErrorCode());
1679 
1680  qDebug() << "qprecursor_list.isActive()=" << qprecursor_list.isActive();
1681  qDebug() << "qprecursor_list.isSelect()=" << qprecursor_list.isSelect();
1682  bool first = true;
1683  SpectrumDescr spectrum_descr;
1684  std::size_t i = 0; /*
1685  while(qprecursor_list.next())
1686  {
1687  qDebug() << "i=" << i;
1688  i++;
1689  }*/
1690  qprecursor_list.last(); // strange bug : get the last sql record and get
1691  // back, unless it will not retrieve all records.
1692 
1693  qDebug() << "qprecursor_list.at()=" << qprecursor_list.at();
1694  qprecursor_list.first();
1695  std::vector<pappso::TimsData::SpectrumDescr> spectrum_description_list;
1696  spectrum_descr.precursor_id = 0;
1697  // std::size_t i = 0;
1698 
1699  do
1700  {
1701 
1702  if(spectrum_descr.precursor_id !=
1703  (std::size_t)qprecursor_list.value(6).toLongLong())
1704  {
1705  // new precursor
1706  if(spectrum_descr.precursor_id > 0)
1707  {
1708  spectrum_description_list.push_back(spectrum_descr);
1709  }
1710 
1711  spectrum_descr.tims_frame_list.clear();
1712  first = true;
1713  }
1714  qDebug() << " qprecursor_list.value(6).toLongLong() ="
1715  << qprecursor_list.value(6).toLongLong();
1716  spectrum_descr.precursor_id =
1717  (std::size_t)qprecursor_list.value(6).toLongLong();
1718  qDebug() << " spectrum_descr.precursor_id ="
1719  << spectrum_descr.precursor_id;
1720  qDebug() << " cumul tims frame:" << qprecursor_list.value(0).toLongLong();
1721  spectrum_descr.tims_frame_list.push_back(
1722  qprecursor_list.value(0).toLongLong());
1723  qDebug() << " first =" << first;
1724  if(first)
1725  {
1726  qDebug();
1727  // mass_spectrum.setPrecursorCharge(q.value(11).toInt());
1728  // mass_spectrum.setPrecursorMz(q.value(10).toDouble());
1729  // mass_spectrum.setPrecursorIntensity(q.value(13).toDouble());
1730  spectrum_descr.precursor_ion_data =
1731  PrecursorIonData(qprecursor_list.value(10).toDouble(),
1732  qprecursor_list.value(11).toInt(),
1733  qprecursor_list.value(13).toDouble());
1734 
1735  // spectrum_descr.precursor_id = q.value(6).toLongLong();
1736  spectrum_descr.ms2_index = (spectrum_descr.precursor_id * 2) - 1;
1737  spectrum_descr.ms1_index = (spectrum_descr.precursor_id * 2) - 2;
1738 
1739  spectrum_descr.scan_mobility_start =
1740  qprecursor_list.value(1).toLongLong();
1741  spectrum_descr.scan_mobility_end =
1742  qprecursor_list.value(2).toLongLong();
1743 
1744  spectrum_descr.isolationMz = qprecursor_list.value(3).toDouble();
1745  spectrum_descr.isolationWidth = qprecursor_list.value(4).toDouble();
1746  spectrum_descr.collisionEnergy = qprecursor_list.value(5).toFloat();
1747  spectrum_descr.parent_frame = qprecursor_list.value(14).toLongLong();
1748 
1749 
1750  first = false;
1751  }
1752  // qDebug() << "qprecursor_list.executedQuery()="
1753  // << qprecursor_list.executedQuery();
1754  // qDebug() << "qprecursor_list.last()=" << qprecursor_list.last();
1755  i++;
1756  }
1757  while(qprecursor_list.next());
1758 
1759  // last One
1760 
1761  // new precursor
1762  if(spectrum_descr.precursor_id > 0)
1763  {
1764  spectrum_description_list.push_back(spectrum_descr);
1765  }
1766 
1767 
1768  QString local_filepath = m_timsDataDirectory.absoluteFilePath("analysis.tdf");
1769 
1770  if(m_isMonoThread)
1771  {
1772  for(SpectrumDescr &spectrum_descr : spectrum_description_list)
1773  {
1774 
1775  std::vector<QualifiedMassSpectrum> mass_spectrum_list;
1776  ms2ReaderGenerateMS1MS2Spectrum(
1777  msrun_id, mass_spectrum_list, handler, spectrum_descr, ms_level);
1778 
1779  for(auto &qualified_spectrum : mass_spectrum_list)
1780  {
1781  handler.setQualifiedMassSpectrum(qualified_spectrum);
1782  }
1783 
1784  if(handler.shouldStop())
1785  {
1786  qDebug() << "The operation was cancelled. Breaking the loop.";
1787  throw ExceptionInterrupted(
1788  QObject::tr("reading TimsTOF job cancelled by the user :\n%1")
1789  .arg(local_filepath));
1790  }
1791  }
1792  }
1793  else
1794  {
1795 
1796 
1797  TimsData *itself = this;
1798  pappso::SpectrumCollectionHandlerInterface *pointer_handler = &handler;
1799 
1800 
1801  std::function<std::vector<QualifiedMassSpectrum>(
1803  generate_spectrum = [itself, msrun_id, pointer_handler, ms_level](
1804  pappso::TimsData::SpectrumDescr &spectrum_descr)
1805  -> std::vector<QualifiedMassSpectrum> {
1806  std::vector<QualifiedMassSpectrum> mass_spectrum_list;
1807  itself->ms2ReaderGenerateMS1MS2Spectrum(msrun_id,
1808  mass_spectrum_list,
1809  *pointer_handler,
1810  spectrum_descr,
1811  ms_level);
1812 
1813 
1814  return mass_spectrum_list;
1815  };
1816 
1817  QFuture<std::size_t> res;
1818  res = QtConcurrent::mappedReduced<std::size_t>(
1819  spectrum_description_list.begin(),
1820  spectrum_description_list.end(),
1821  generate_spectrum,
1822  [pointer_handler, res, local_filepath](
1823  std::size_t &result,
1824  std::vector<QualifiedMassSpectrum> qualified_spectrum_list) {
1825  for(auto &qualified_spectrum : qualified_spectrum_list)
1826  {
1827  pointer_handler->setQualifiedMassSpectrum(qualified_spectrum);
1828  }
1829 
1830  if(pointer_handler->shouldStop())
1831  {
1832  qDebug() << "The operation was cancelled. Breaking the loop.";
1833  throw ExceptionInterrupted(
1834  QObject::tr("reading TimsTOF job cancelled by the user :\n%1")
1835  .arg(local_filepath));
1836  }
1837  result++;
1838  },
1839  QtConcurrent::OrderedReduce);
1840  res.waitForFinished();
1841  }
1842  handler.loadingEnded();
1843  mpa_timsBinDec->closeLinearRead();
1844 }
1845 
1846 
1847 void
1849  const MsRunIdCstSPtr &msrun_id,
1850  std::vector<QualifiedMassSpectrum> &qualified_mass_spectrum_list,
1852  pappso::TimsData::SpectrumDescr &spectrum_descr,
1853  unsigned int ms_level)
1854 {
1855 
1856  qDebug() << " ms_level=" << ms_level;
1857  // The handler will receive the index of the mass spectrum in the
1858  // current run via the mass spectrum id member datum.
1859  if((ms_level == 0) || (ms_level == 1))
1860  {
1861  qualified_mass_spectrum_list.push_back(QualifiedMassSpectrum());
1863  msrun_id,
1864  qualified_mass_spectrum_list.back(),
1865  spectrum_descr,
1866  handler.needMsLevelPeakList(1));
1867  }
1868  if((ms_level == 0) || (ms_level == 2))
1869  {
1870  qualified_mass_spectrum_list.push_back(QualifiedMassSpectrum());
1872  msrun_id,
1873  qualified_mass_spectrum_list.back(),
1874  spectrum_descr,
1875  handler.needMsLevelPeakList(2));
1876  }
1877  qDebug();
1878 }
1879 
1880 
1883 {
1884 
1885  SpectrumDescr spectrum_descr;
1886  QSqlDatabase qdb = openDatabaseConnection();
1887  QSqlQuery q = qdb.exec(QString("SELECT PasefFrameMsMsInfo.Frame, " // 0
1888  "PasefFrameMsMsInfo.ScanNumBegin, " // 1
1889  "PasefFrameMsMsInfo.ScanNumEnd, " // 2
1890  "PasefFrameMsMsInfo.IsolationMz, " // 3
1891  "PasefFrameMsMsInfo.IsolationWidth, " // 4
1892  "PasefFrameMsMsInfo.CollisionEnergy, " // 5
1893  "PasefFrameMsMsInfo.Precursor, " // 6
1894  "Precursors.Id, " // 7
1895  "Precursors.LargestPeakMz, " // 8
1896  "Precursors.AverageMz, " // 9
1897  "Precursors.MonoisotopicMz, " // 10
1898  "Precursors.Charge, " // 11
1899  "Precursors.ScanNumber, " // 12
1900  "Precursors.Intensity, " // 13
1901  "Precursors.Parent " // 14
1902  "FROM PasefFrameMsMsInfo "
1903  "INNER JOIN Precursors ON "
1904  "PasefFrameMsMsInfo.Precursor=Precursors.Id "
1905  "WHERE Precursors.Id=%1;")
1906  .arg(precursor_id));
1907  if(q.lastError().isValid())
1908  {
1909 
1910  throw PappsoException(
1911  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
1912  "command %2:\n%3\n%4\n%5")
1913  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
1914  .arg(q.lastQuery())
1915  .arg(qdb.lastError().databaseText())
1916  .arg(qdb.lastError().driverText())
1917  .arg(qdb.lastError().nativeErrorCode()));
1918  }
1919 
1920 
1921  bool first = true;
1922  while(q.next())
1923  {
1924 
1925  qDebug() << " cumul tims frame:" << q.value(0).toLongLong();
1926  spectrum_descr.tims_frame_list.push_back(q.value(0).toLongLong());
1927  if(first)
1928  {
1929  // mass_spectrum.setPrecursorCharge(q.value(11).toInt());
1930  // mass_spectrum.setPrecursorMz(q.value(10).toDouble());
1931  // mass_spectrum.setPrecursorIntensity(q.value(13).toDouble());
1932  spectrum_descr.precursor_ion_data =
1933  PrecursorIonData(q.value(10).toDouble(),
1934  q.value(11).toInt(),
1935  q.value(13).toDouble());
1936 
1937  spectrum_descr.precursor_id = q.value(6).toLongLong();
1938  spectrum_descr.ms2_index = (spectrum_descr.precursor_id * 2) - 1;
1939  spectrum_descr.ms1_index = (spectrum_descr.precursor_id * 2) - 2;
1940 
1941  spectrum_descr.scan_mobility_start = q.value(1).toLongLong();
1942  spectrum_descr.scan_mobility_end = q.value(2).toLongLong();
1943 
1944  spectrum_descr.isolationMz = q.value(3).toDouble();
1945  spectrum_descr.isolationWidth = q.value(4).toDouble();
1946  spectrum_descr.collisionEnergy = q.value(5).toFloat();
1947  spectrum_descr.parent_frame = q.value(14).toLongLong();
1948 
1949 
1950  first = false;
1951  }
1952  }
1953  if(spectrum_descr.precursor_id == 0)
1954  {
1955  throw ExceptionNotFound(
1956  QObject::tr("ERROR in %1 %2 : precursor id (%3) NOT FOUND ")
1957  .arg(__FILE__)
1958  .arg(__FUNCTION__)
1959  .arg(precursor_id));
1960  }
1961  return spectrum_descr;
1962 }
1963 
1964 std::vector<double>
1966 {
1967  std::vector<double> timeline;
1968  timeline.reserve(m_mapFramesRecord.size());
1969  for(const TimsFrameRecord &frame_record : m_mapFramesRecord)
1970  {
1971  if(frame_record.mz_calibration_id != 0)
1972  {
1973  timeline.push_back(frame_record.frame_time);
1974  }
1975  }
1976  return timeline;
1977 }
1978 
1981  const std::pair<std::size_t, std::size_t> &scan_coordinate)
1982 {
1983 
1984  SpectrumDescr spectrum_descr;
1985  QSqlDatabase qdb = openDatabaseConnection();
1986  QSqlQuery q =
1987  qdb.exec(QString("SELECT PasefFrameMsMsInfo.Frame, " // 0
1988  "PasefFrameMsMsInfo.ScanNumBegin, " // 1
1989  "PasefFrameMsMsInfo.ScanNumEnd, " // 2
1990  "PasefFrameMsMsInfo.IsolationMz, " // 3
1991  "PasefFrameMsMsInfo.IsolationWidth, " // 4
1992  "PasefFrameMsMsInfo.CollisionEnergy, " // 5
1993  "PasefFrameMsMsInfo.Precursor, " // 6
1994  "Precursors.Id, " // 7
1995  "Precursors.LargestPeakMz, " // 8
1996  "Precursors.AverageMz, " // 9
1997  "Precursors.MonoisotopicMz, " // 10
1998  "Precursors.Charge, " // 11
1999  "Precursors.ScanNumber, " // 12
2000  "Precursors.Intensity, " // 13
2001  "Precursors.Parent " // 14
2002  "FROM PasefFrameMsMsInfo "
2003  "INNER JOIN Precursors ON "
2004  "PasefFrameMsMsInfo.Precursor=Precursors.Id "
2005  "WHERE "
2006  "PasefFrameMsMsInfo.Frame=%1 and "
2007  "(PasefFrameMsMsInfo.ScanNumBegin "
2008  "<= %2 and PasefFrameMsMsInfo.ScanNumEnd >= %2);")
2009  .arg(scan_coordinate.first)
2010  .arg(scan_coordinate.second));
2011  if(q.lastError().isValid())
2012  {
2013 
2014  throw PappsoException(
2015  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
2016  "command %2:\n%3\n%4\n%5")
2017  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
2018  .arg(q.lastQuery())
2019  .arg(qdb.lastError().databaseText())
2020  .arg(qdb.lastError().driverText())
2021  .arg(qdb.lastError().nativeErrorCode()));
2022  }
2023 
2024  if(q.next())
2025  {
2026 
2027  qDebug() << " cumul tims frame:" << q.value(0).toLongLong();
2028  spectrum_descr.tims_frame_list.push_back(q.value(0).toLongLong());
2029  // mass_spectrum.setPrecursorCharge(q.value(11).toInt());
2030  // mass_spectrum.setPrecursorMz(q.value(10).toDouble());
2031  // mass_spectrum.setPrecursorIntensity(q.value(13).toDouble());
2032  spectrum_descr.precursor_ion_data = PrecursorIonData(
2033  q.value(10).toDouble(), q.value(11).toInt(), q.value(13).toDouble());
2034 
2035  spectrum_descr.precursor_id = q.value(6).toLongLong();
2036  spectrum_descr.ms2_index = (spectrum_descr.precursor_id * 2) - 1;
2037  spectrum_descr.ms1_index = (spectrum_descr.precursor_id * 2) - 2;
2038 
2039  spectrum_descr.scan_mobility_start = q.value(1).toLongLong();
2040  spectrum_descr.scan_mobility_end = q.value(2).toLongLong();
2041 
2042  spectrum_descr.isolationMz = q.value(3).toDouble();
2043  spectrum_descr.isolationWidth = q.value(4).toDouble();
2044  spectrum_descr.collisionEnergy = q.value(5).toFloat();
2045  spectrum_descr.parent_frame = q.value(14).toLongLong();
2046  }
2047  return spectrum_descr;
2048 }
2049 
2050 
2051 void
2053  pappso::TimsData::SpectrumDescr &spectrum_descr, QSqlQuery &qprecursor_list)
2054 {
2055 
2056  spectrum_descr.tims_frame_list.clear();
2057  spectrum_descr.tims_frame_list.push_back(
2058  qprecursor_list.value(0).toLongLong());
2059  // mass_spectrum.setPrecursorCharge(q.value(11).toInt());
2060  // mass_spectrum.setPrecursorMz(q.value(10).toDouble());
2061  // mass_spectrum.setPrecursorIntensity(q.value(13).toDouble());
2062  spectrum_descr.precursor_ion_data =
2063  PrecursorIonData(qprecursor_list.value(10).toDouble(),
2064  qprecursor_list.value(11).toInt(),
2065  qprecursor_list.value(13).toDouble());
2066 
2067  spectrum_descr.precursor_id = qprecursor_list.value(6).toLongLong();
2068  spectrum_descr.ms2_index = (spectrum_descr.precursor_id * 2) - 1;
2069  spectrum_descr.ms1_index = (spectrum_descr.precursor_id * 2) - 2;
2070 
2071  spectrum_descr.scan_mobility_start = qprecursor_list.value(1).toLongLong();
2072  spectrum_descr.scan_mobility_end = qprecursor_list.value(2).toLongLong();
2073 
2074  spectrum_descr.isolationMz = qprecursor_list.value(3).toDouble();
2075  spectrum_descr.isolationWidth = qprecursor_list.value(4).toDouble();
2076  spectrum_descr.collisionEnergy = qprecursor_list.value(5).toFloat();
2077  spectrum_descr.parent_frame = qprecursor_list.value(14).toLongLong();
2078 }
2079 
2080 
2081 void
2083  const pappso::MsRunIdCstSPtr &msrun_id,
2085  unsigned int ms_level)
2086 {
2087 
2088  if(!m_hasPrecursorTable)
2089  {
2090  throw PappsoException(
2091  QObject::tr("unable to read spectrum list : this data file does not "
2092  "contain MS2 data, no precursor found."));
2093  }
2094 
2095  // We'll need it to perform the looping in the spectrum list.
2096  std::size_t spectrum_list_size = getTotalNumberOfScans();
2097 
2098  // qDebug() << "The spectrum list has size:" << spectrum_list_size;
2099 
2100  // Inform the handler of the spectrum list so that it can handle feedback to
2101  // the user.
2102  handler.spectrumListHasSize(spectrum_list_size);
2103 
2104  QSqlDatabase qdb = openDatabaseConnection();
2105  QSqlQuery qprecursor_list = qdb.exec(QString(
2106  "SELECT DISTINCT "
2107  "PasefFrameMsMsInfo.Frame, " // 0
2108  "PasefFrameMsMsInfo.ScanNumBegin, " // 1
2109  "PasefFrameMsMsInfo.ScanNumEnd, " // 2
2110  "PasefFrameMsMsInfo.IsolationMz, " // 3
2111  "PasefFrameMsMsInfo.IsolationWidth, " // 4
2112  "PasefFrameMsMsInfo.CollisionEnergy, " // 5
2113  "PasefFrameMsMsInfo.Precursor, " // 6
2114  "Precursors.Id, " // 7
2115  "Precursors.LargestPeakMz, " // 8
2116  "Precursors.AverageMz, " // 9
2117  "Precursors.MonoisotopicMz, " // 10
2118  "Precursors.Charge, " // 11
2119  "Precursors.ScanNumber, " // 12
2120  "Precursors.Intensity, " // 13
2121  "Precursors.Parent " // 14
2122  "FROM PasefFrameMsMsInfo "
2123  "INNER JOIN Precursors ON "
2124  "PasefFrameMsMsInfo.Precursor=Precursors.Id "
2125  "ORDER BY PasefFrameMsMsInfo.Frame, PasefFrameMsMsInfo.ScanNumBegin ;"));
2126  if(qprecursor_list.lastError().isValid())
2127  {
2128  throw PappsoException(
2129  QObject::tr("ERROR in TIMS sqlite database file %1, executing SQL "
2130  "command %2:\n%3\n%4\n%5")
2131  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf"))
2132  .arg(qprecursor_list.lastQuery())
2133  .arg(qdb.lastError().databaseText())
2134  .arg(qdb.lastError().driverText())
2135  .arg(qdb.lastError().nativeErrorCode()));
2136  }
2137 
2138 
2139  std::size_t i = 0; // iterate on each Spectrum
2140 
2141  qprecursor_list.last(); // strange bug : get the last sql record and get
2142  // back, unless it will not retrieve all records.
2143 
2144  qDebug() << "qprecursor_list.at()=" << qprecursor_list.at();
2145  qprecursor_list.first();
2146 
2147  TimsFrameBaseCstSPtr tims_frame;
2148  SpectrumDescr spectrum_descr;
2149 
2150  for(FrameIdDescr &current_frame : m_frameIdDescrList)
2151  {
2152 
2153  // If the user of this reader instance wants to stop reading the
2154  // spectra, then break this loop.
2155  if(handler.shouldStop())
2156  {
2157  qDebug() << "The operation was cancelled. Breaking the loop.";
2158  throw ExceptionInterrupted(
2159  QObject::tr("reading TimsTOF job cancelled by the user :\n%1")
2160  .arg(m_timsDataDirectory.absoluteFilePath("analysis.tdf")));
2161  }
2162 
2163  tims_frame = getTimsFrameBaseCstSPtrCached(current_frame.m_frameId);
2164  unsigned int tims_ms_level = tims_frame.get()->getMsLevel();
2165 
2166  if((ms_level != 0) && (ms_level != tims_ms_level))
2167  { // bypass
2168  i += current_frame.m_size;
2169  }
2170  else
2171  {
2172  bool want_binary_data = handler.needMsLevelPeakList(tims_ms_level);
2173  qDebug() << "want_binary_data=" << want_binary_data;
2174  if(want_binary_data)
2175  {
2176  qDebug() << "bindec";
2177  tims_frame = getTimsFrameCstSPtrCached(current_frame.m_frameId);
2178  }
2179 
2180  bool possible_precursor = false;
2181  if(tims_ms_level == 2)
2182  {
2183  // seek the precursor record:
2184  while(qprecursor_list.value(0).toULongLong() <
2185  current_frame.m_frameId)
2186  {
2187  qprecursor_list.next();
2188 
2189  if(qprecursor_list.value(0).toULongLong() ==
2190  current_frame.m_frameId)
2191  {
2192  possible_precursor = true;
2193  }
2194  fillSpectrumDescriptionWithSqlRecord(spectrum_descr,
2195  qprecursor_list);
2196  }
2197  }
2198 
2199  for(std::size_t scan_num = 0; scan_num < current_frame.m_size;
2200  scan_num++)
2201  {
2202  bool has_a_precursor = false;
2203  if(possible_precursor)
2204  {
2205  if(spectrum_descr.scan_mobility_end < scan_num)
2206  {
2207  // seek the precursor record:
2208  while(qprecursor_list.value(0).toULongLong() <
2209  current_frame.m_frameId)
2210  {
2211  qprecursor_list.next();
2212 
2213  if(qprecursor_list.value(0).toULongLong() !=
2214  current_frame.m_frameId)
2215  {
2216  possible_precursor = false;
2217  }
2218  fillSpectrumDescriptionWithSqlRecord(spectrum_descr,
2219  qprecursor_list);
2220  }
2221  }
2222 
2223  if(possible_precursor &&
2224  (spectrum_descr.scan_mobility_start < scan_num))
2225  {
2226  // we are in
2227  has_a_precursor = true;
2228  }
2229  } // end to determine if we are in a precursor for this
2230  // spectrum
2231 
2232  QualifiedMassSpectrum mass_spectrum;
2233 
2234 
2235  MassSpectrumId spectrum_id;
2236 
2237  spectrum_id.setSpectrumIndex(i);
2238  spectrum_id.setMsRunId(msrun_id);
2239  spectrum_id.setNativeId(QString("frame=%1 scan=%2 index=%3")
2240  .arg(current_frame.m_frameId)
2241  .arg(scan_num)
2242  .arg(i));
2243 
2244  mass_spectrum.setMassSpectrumId(spectrum_id);
2245 
2246  mass_spectrum.setMsLevel(tims_frame.get()->getMsLevel());
2247  mass_spectrum.setRtInSeconds(tims_frame.get()->getTime());
2248 
2249  mass_spectrum.setDtInMilliSeconds(
2250  tims_frame.get()->getDriftTime(scan_num));
2251  // 1/K0
2252  mass_spectrum.setParameterValue(
2254  tims_frame.get()->getOneOverK0Transformation(scan_num));
2255 
2256  mass_spectrum.setEmptyMassSpectrum(true);
2257  if(want_binary_data)
2258  {
2259  try
2260  {
2261  mass_spectrum.setMassSpectrumSPtr(
2262  tims_frame.get()->getMassSpectrumSPtr(scan_num));
2263  }
2264  catch(PappsoException &error)
2265  {
2266  throw PappsoException(
2267  QObject::tr(
2268  "ERROR in %1 (scan_num=%2 spectrum_index=%3):\n%4")
2269  .arg(__FUNCTION__)
2270  .arg(scan_num)
2271  .arg(spectrum_id.getSpectrumIndex())
2272  .arg(error.qwhat()));
2273  }
2274  if(mass_spectrum.size() > 0)
2275  {
2276  mass_spectrum.setEmptyMassSpectrum(false);
2277  }
2278  }
2279  else
2280  {
2281  // if(tims_frame.get()->getNbrPeaks(coordinate.second) > 0)
2282  //{
2283  mass_spectrum.setEmptyMassSpectrum(false);
2284  // }
2285  }
2286  if(has_a_precursor)
2287  {
2288  if(spectrum_descr.precursor_id > 0)
2289  {
2290 
2291  mass_spectrum.appendPrecursorIonData(
2292  spectrum_descr.precursor_ion_data);
2293 
2294  std::size_t prec_spectrum_index =
2295  getRawIndexFromCoordinate(spectrum_descr.parent_frame,
2296  scan_num);
2297 
2298  mass_spectrum.setPrecursorSpectrumIndex(
2299  prec_spectrum_index);
2300  mass_spectrum.setPrecursorNativeId(
2301  QString("frame=%1 scan=%2 index=%3")
2302  .arg(spectrum_descr.parent_frame)
2303  .arg(scan_num)
2304  .arg(prec_spectrum_index));
2305 
2306  mass_spectrum.setParameterValue(
2308  spectrum_descr.isolationMz);
2309  mass_spectrum.setParameterValue(
2311  spectrum_descr.isolationWidth);
2312 
2313  mass_spectrum.setParameterValue(
2315  spectrum_descr.collisionEnergy);
2316  mass_spectrum.setParameterValue(
2318  (quint64)spectrum_descr.precursor_id);
2319  }
2320  }
2321 
2322  handler.setQualifiedMassSpectrum(mass_spectrum);
2323  i++;
2324  }
2325  }
2326  }
2327 }
2328 
2329 std::map<quint32, quint32>
2331 {
2332 
2333  qDebug() << " spectrum_index=" << spectrum_index;
2334  auto coordinate = getScanCoordinateFromRawIndex(spectrum_index);
2335  TimsFrameBaseCstSPtr tims_frame;
2336  tims_frame = getTimsFrameCstSPtrCached(coordinate.first);
2337 
2338  std::map<quint32, quint32> raw_spectrum;
2339  tims_frame.get()->cumulateScansInRawMap(
2340  raw_spectrum, coordinate.second, coordinate.second);
2341  return raw_spectrum;
2342 }
2343 } // namespace pappso
Trace toTrace() const
Definition: maptrace.cpp:218
void setNativeId(const QString &native_id)
void setMsRunId(MsRunIdCstSPtr other)
std::size_t getSpectrumIndex() const
void setSpectrumIndex(std::size_t index)
Class to represent a mass spectrum.
Definition: massspectrum.h:71
MzCalibrationInterfaceSPtr getInstance(double T1_frame, double T2_frame, const QSqlRecord &mzcalibration_record)
virtual const QString & qwhat() const
const char * what() const noexcept override
Class representing a fully specified mass spectrum.
void setPrecursorNativeId(const QString &native_id)
Set the scan native id of the precursor ion.
void setDtInMilliSeconds(pappso_double rt)
Set the drift time in milliseconds.
void appendPrecursorIonData(const PrecursorIonData &precursor_ion_data)
void setMassSpectrumId(const MassSpectrumId &iD)
Set the MassSpectrumId.
void setMsLevel(uint ms_level)
Set the mass spectrum level.
void setPrecursorSpectrumIndex(std::size_t precursor_scan_num)
Set the scan number of the precursor ion.
void setParameterValue(QualifiedMassSpectrumParameter parameter, const QVariant &value)
void setMassSpectrumSPtr(MassSpectrumSPtr massSpectrum)
Set the MassSpectrumSPtr.
void setRtInSeconds(pappso_double rt)
Set the retention time in seconds.
void setEmptyMassSpectrum(bool is_empty_mass_spectrum)
interface to collect spectrums from the MsRunReader class
Definition: msrunreader.h:56
virtual bool needMsLevelPeakList(unsigned int ms_level) const final
tells if we need the peak list (if we want the binary data) for each spectrum, given an MS level
Definition: msrunreader.cpp:69
virtual void spectrumListHasSize(std::size_t size)
Definition: msrunreader.cpp:52
virtual void setQualifiedMassSpectrum(const QualifiedMassSpectrum &spectrum)=0
TimsFrameSPtr getTimsFrameSPtrByOffset(std::size_t frameId, const std::vector< pappso::TimsFrameRecord > &frame_record_list)
Definition: timsbindec.cpp:147
QSqlDatabase openDatabaseConnection() const
Definition: timsdata.cpp:236
TimsFrameCstSPtr getTimsFrameCstSPtr(std::size_t timsId)
get a Tims frame with his database ID
Definition: timsdata.cpp:530
std::vector< FrameIdDescr > m_frameIdDescrList
store every frame id and corresponding sizes
Definition: timsdata.h:328
void ms2ReaderSpectrumCollectionByMsLevel(const MsRunIdCstSPtr &msrun_id, SpectrumCollectionHandlerInterface &handler, unsigned int ms_level)
function to visit an MsRunReader and get each Spectrum in a spectrum collection handler by Ms Levels
Definition: timsdata.cpp:1622
TimsFrameCstSPtr getTimsFrameCstSPtrCached(std::size_t timsId)
get a Tims frame with his database ID but look in the cache first
Definition: timsdata.cpp:1088
pappso::MassSpectrumCstSPtr getMassSpectrumCstSPtrByRawIndex(std::size_t raw_index)
get a mass spectrum given its spectrum index
Definition: timsdata.cpp:401
virtual ~TimsData()
Definition: timsdata.cpp:271
SpectrumDescr getSpectrumDescrWithPrecursorId(std::size_t precursor_id)
get an intermediate structure describing a spectrum
Definition: timsdata.cpp:1882
TimsData(QDir timsDataDirectory)
build using the tims data directory
Definition: timsdata.cpp:49
std::map< quint32, quint32 > getRawMs2ByPrecursorId(std::size_t precursor_index)
get cumulated raw signal for a given precursor only to use to see the raw signal
Definition: timsdata.cpp:1282
std::size_t m_totalNumberOfFrames
Definition: timsdata.h:297
TimsFrameBaseCstSPtr getTimsFrameBaseCstSPtrCached(std::size_t timsId)
Definition: timsdata.cpp:1067
std::size_t m_totalNumberOfScans
Definition: timsdata.h:295
std::deque< TimsFrameCstSPtr > m_timsFrameCache
Definition: timsdata.h:299
pappso::MassSpectrumCstSPtr getMassSpectrumCstSPtr(std::size_t timsId, std::size_t scanNum)
get a mass spectrum given the tims frame database id and scan number within tims frame
Definition: timsdata.cpp:616
std::size_t m_cacheSize
Definition: timsdata.h:298
std::pair< std::size_t, std::size_t > getScanCoordinateFromRawIndex(std::size_t spectrum_index) const
Definition: timsdata.cpp:342
std::vector< TimsFrameRecord > m_mapFramesRecord
Definition: timsdata.h:312
std::vector< std::size_t > getClosestPrecursorIdByMz(std::vector< std::vector< double >> ids, double mz_value)
Definition: timsdata.cpp:743
std::map< int, QSqlRecord > m_mapMzCalibrationRecord
Definition: timsdata.h:310
void ms2ReaderGenerateMS1MS2Spectrum(const MsRunIdCstSPtr &msrun_id, std::vector< QualifiedMassSpectrum > &qualified_mass_spectrum_list, SpectrumCollectionHandlerInterface &handler, SpectrumDescr &spectrum_descr, unsigned int ms_level)
Definition: timsdata.cpp:1848
std::vector< std::size_t > getPrecursorsFromMzRtCharge(int charge, double mz_val, double rt_sec, double k0)
guess possible precursor ids given a charge, m/z, retention time and k0
Definition: timsdata.cpp:638
void fillSpectrumDescriptionWithSqlRecord(SpectrumDescr &spectrum_descr, QSqlQuery &qprecursor_list)
Definition: timsdata.cpp:2052
std::map< int, QSqlRecord > m_mapTimsCalibrationRecord
Definition: timsdata.h:311
QMutex m_mutex
Definition: timsdata.h:346
bool m_builtinMs2Centroid
enable builtin centroid on raw tims integers by default
Definition: timsdata.h:307
void setMs2BuiltinCentroid(bool centroid)
enable or disable simple centroid filter on raw tims data for MS2
Definition: timsdata.cpp:285
void fillFrameIdDescrList()
private function to fill m_frameIdDescrList
Definition: timsdata.cpp:297
TimsFrameBaseCstSPtr getTimsFrameBaseCstSPtr(std::size_t timsId)
get a Tims frame base (no binary data file access) with his database ID
Definition: timsdata.cpp:422
bool m_hasPrecursorTable
Definition: timsdata.h:344
void rawReaderSpectrumCollectionByMsLevel(const MsRunIdCstSPtr &msrun_id, SpectrumCollectionHandlerInterface &handler, unsigned int ms_level)
function to visit an MsRunReader and get each raw Spectrum in a spectrum collection handler by Ms Lev...
Definition: timsdata.cpp:2082
QDir m_timsDataDirectory
Definition: timsdata.h:292
std::size_t getTotalNumberOfPrecursors() const
get the number of precursors analyzes by PASEF
Definition: timsdata.cpp:632
MzCalibrationStore * mpa_mzCalibrationStore
Definition: timsdata.h:315
std::vector< std::size_t > getTimsMS1FrameIdRange(double rt_begin, double rt_end) const
Definition: timsdata.cpp:494
virtual std::vector< double > getRetentionTimeLine() const
retention timeline get retention times along the MSrun in seconds
Definition: timsdata.cpp:1965
unsigned int getMsLevelBySpectrumIndex(std::size_t spectrum_index)
Definition: timsdata.cpp:767
bool getMs2BuiltinCentroid() const
tells if simple centroid filter on raw tims data for MS2 is enabled or not
Definition: timsdata.cpp:291
void getQualifiedMs1MassSpectrumByPrecursorId(const MsRunIdCstSPtr &msrun_id, QualifiedMassSpectrum &mass_spectrum, SpectrumDescr &spectrum_descr, bool want_binary_data)
Definition: timsdata.cpp:966
SpectrumDescr getSpectrumDescrWithScanCoordinate(const std::pair< std::size_t, std::size_t > &scan_coordinate)
Definition: timsdata.cpp:1980
std::map< quint32, quint32 > getRawMsBySpectrumIndex(std::size_t spectrum_index)
get raw signal for a spectrum index only to use to see the raw signal
Definition: timsdata.cpp:2330
std::deque< TimsFrameBaseCstSPtr > m_timsFrameBaseCache
Definition: timsdata.h:300
std::map< std::size_t, std::size_t > m_thousandIndexToFrameIdDescrListIndex
index to find quickly a frameId in the description list with the raw index of spectrum modulo 1000 @k...
Definition: timsdata.h:335
TimsBinDec * mpa_timsBinDec
Definition: timsdata.h:293
void getQualifiedMassSpectrumByRawIndex(const MsRunIdCstSPtr &msrun_id, QualifiedMassSpectrum &mass_spectrum, std::size_t spectrum_index, bool want_binary_data)
Definition: timsdata.cpp:776
void setMs1FilterCstSPtr(pappso::FilterInterfaceCstSPtr &filter)
filter interface to apply just after raw MS1 specturm extraction the filter can be a list of filters ...
Definition: timsdata.cpp:1170
void setMs2FilterCstSPtr(pappso::FilterInterfaceCstSPtr &filter)
filter interface to apply just after raw MS2 specturm extraction the filter can be a list of filters ...
Definition: timsdata.cpp:1165
pappso::FilterInterfaceCstSPtr mcsp_ms1Filter
Definition: timsdata.h:303
bool m_isMonoThread
Definition: timsdata.h:342
std::size_t getTotalNumberOfScans() const
get the total number of scans
Definition: timsdata.cpp:625
std::vector< std::size_t > getMatchPrecursorIdByKo(std::vector< std::vector< double >> ids, double ko_value)
Definition: timsdata.cpp:717
std::size_t getRawIndexFromCoordinate(std::size_t frame_id, std::size_t scan_num) const
Definition: timsdata.cpp:378
pappso::FilterInterfaceCstSPtr mcsp_ms2Filter
Definition: timsdata.h:302
void getQualifiedMs2MassSpectrumByPrecursorId(const MsRunIdCstSPtr &msrun_id, QualifiedMassSpectrum &mass_spectrum, SpectrumDescr &spectrum_descr, bool want_binary_data)
Definition: timsdata.cpp:1410
void setMonoThread(bool is_mono_thread)
set only one is_mono_thread to true
Definition: timsdata.cpp:230
Trace getTicChromatogram()
Definition: timsdata.cpp:882
XicCoordTims getXicCoordTimsFromPrecursorId(std::size_t precursor_id, PrecisionPtr precision_ptr)
Definition: timsdata.cpp:1176
std::size_t m_totalNumberOfPrecursors
Definition: timsdata.h:296
virtual MapTrace & combine(MapTrace &map_trace, const Trace &trace) const override
A simple container of DataPoint instances.
Definition: trace.h:147
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition: aa.cpp:39
std::shared_ptr< const TimsFrameBase > TimsFrameBaseCstSPtr
Definition: timsframebase.h:41
std::shared_ptr< TimsFrame > TimsFrameSPtr
Definition: timsframe.h:41
std::shared_ptr< TimsFrameBase > TimsFrameBaseSPtr
Definition: timsframebase.h:39
std::shared_ptr< const MsRunId > MsRunIdCstSPtr
Definition: msrunid.h:44
std::shared_ptr< const MassSpectrum > MassSpectrumCstSPtr
Definition: massspectrum.h:55
@ CollisionEnergy
Bruker's Tims tof collision energy.
@ OneOverK0end
1/k0 of last acquisition for composite pasef MS/MS spectrum
@ IsolationWidth
isolation window width
@ BrukerPrecursorIndex
Bruker's Tims tof precursor index.
@ rt
Retention time.
std::shared_ptr< const FilterInterface > FilterInterfaceCstSPtr
std::shared_ptr< const TimsFrame > TimsFrameCstSPtr
Definition: timsframe.h:43
std::vector< std::size_t > tims_frame_list
Definition: timsdata.h:125
PrecursorIonData precursor_ion_data
Definition: timsdata.h:126
std::size_t tims_calibration_id
coordinates of the XIC to extract and the resulting XIC after extraction
Definition: xiccoordtims.h:51
std::size_t scanNumEnd
mobility index end
Definition: xiccoordtims.h:91
std::size_t scanNumBegin
mobility index begin
Definition: xiccoordtims.h:87
XicSPtr xicSptr
extracted xic
Definition: xiccoord.h:113
double rtTarget
the targeted retention time to extract around intended in seconds, and related to one msrun....
Definition: xiccoord.h:109
MzRange mzRange
the mass to extract
Definition: xiccoord.h:103
main Tims data handler
minimum functions to extract XICs from Tims Data