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