libpappsomspp
Library for mass spectrometry
timsframebase.cpp
Go to the documentation of this file.
1/**
2 * \file pappsomspp/vendors/tims/timsframebase.cpp
3 * \date 16/12/2019
4 * \author Olivier Langella
5 * \brief handle a single Bruker's TimsTof frame without binary data
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#include "timsframebase.h"
28#include "../../../pappsomspp/pappsoexception.h"
29#include "../../../pappsomspp/exception/exceptionoutofrange.h"
31#include <QDebug>
32#include <QObject>
33#include <cmath>
34#include <algorithm>
35
36namespace pappso
37{
38
39TimsFrameBase::TimsFrameBase(std::size_t timsId, quint32 scanNum)
40{
41 qDebug() << timsId;
42 m_timsId = timsId;
43
44 m_scanNumber = scanNum;
45}
46
47TimsFrameBase::TimsFrameBase([[maybe_unused]] const TimsFrameBase &other)
48{
49}
50
52{
53}
54
55void
56TimsFrameBase::setAccumulationTime(double accumulation_time_ms)
57{
58 m_accumulationTime = accumulation_time_ms;
59}
60
61
62void
64 double T2_frame,
65 double digitizerTimebase,
66 double digitizerDelay,
67 double C0,
68 double C1,
69 double C2,
70 double C3,
71 double C4,
72 double T1_ref,
73 double T2_ref,
74 double dC1,
75 double dC2)
76{
77
78 /* MzCalibrationModel1 mzCalibration(temperature_correction,
79 digitizerTimebase,
80 digitizerDelay,
81 C0,
82 C1,
83 C2,
84 C3,
85 C4);
86 */
87 msp_mzCalibration = std::make_shared<MzCalibrationModel1>(T1_frame,
88 T2_frame,
89 digitizerTimebase,
90 digitizerDelay,
91 C0,
92 C1,
93 C2,
94 C3,
95 C4,
96 T1_ref,
97 T2_ref,
98 dC1,
99 dC2);
100}
101
102bool
103TimsFrameBase::checkScanNum(std::size_t scanNum) const
104{
105 if(scanNum >= m_scanNumber)
106 {
108 QObject::tr("Invalid scan number : scanNum %1 > m_scanNumber %2")
109 .arg(scanNum)
110 .arg(m_scanNumber));
111 }
112
113 return true;
114}
115
116std::size_t
117TimsFrameBase::getNbrPeaks(std::size_t scanNum) const
118{
119 throw PappsoException(
120 QObject::tr(
121 "ERROR unable to get number of peaks in TimsFrameBase for scan number %1")
122 .arg(scanNum));
123}
124
125std::size_t
127{
128 return m_scanNumber;
129}
130
132TimsFrameBase::getMassSpectrumSPtr(std::size_t scanNum) const
133{
134 throw PappsoException(
135 QObject::tr(
136 "ERROR unable to getMassSpectrumSPtr in TimsFrameBase for scan number %1")
137 .arg(scanNum));
138}
139Trace
140TimsFrameBase::cumulateScanToTrace(std::size_t scanNumBegin,
141 std::size_t scanNumEnd) const
142{
143 throw PappsoException(
144 QObject::tr("ERROR unable to cumulateScanToTrace in TimsFrameBase for scan "
145 "number begin %1 end %2")
146 .arg(scanNumBegin)
147 .arg(scanNumEnd));
148}
149void
150TimsFrameBase::cumulateScansInRawMap(std::map<quint32, quint32> &rawSpectrum
151 [[maybe_unused]],
152 std::size_t scanNumBegin,
153 std::size_t scanNumEnd) const
154{
155 throw PappsoException(
156 QObject::tr(
157 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
158 "number begin %1 end %2")
159 .arg(scanNumBegin)
160 .arg(scanNumEnd));
161}
162
163
164quint64
166{
167 throw PappsoException(
168 QObject::tr(
169 "ERROR unable to cumulateSingleScanIntensities in TimsFrameBase for scan "
170 "number %1.").arg(scanNum));
171
172 return 0;
173}
174
175
176quint64
178 std::size_t scanNumEnd) const
179{
180 throw PappsoException(
181 QObject::tr(
182 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
183 "number begin %1 end %2")
184 .arg(scanNumBegin)
185 .arg(scanNumEnd));
186
187 return 0;
188}
189
190void
192{
193 m_time = time;
194}
195
196void
198{
199
200 qDebug() << " m_msMsType=" << type;
201 m_msMsType = type;
202}
203
204unsigned int
206{
207 if(m_msMsType == 0)
208 return 1;
209 return 2;
210}
211
212double
214{
215 return m_time;
216}
217
218std::size_t
220{
221 return m_timsId;
222}
223void
225 double C0,
226 double C1,
227 double C2,
228 double C3,
229 double C4,
230 [[maybe_unused]] double C5,
231 double C6,
232 double C7,
233 double C8,
234 double C9)
235{
236 if(tims_model_type != 2)
237 {
238 throw pappso::PappsoException(QObject::tr(
239 "ERROR in TimsFrame::setTimsCalibration tims_model_type != 2"));
240 }
241 m_timsDvStart = C2; // C2 from TimsCalibration
242 m_timsTtrans = C4; // C4 from TimsCalibration
243 m_timsNdelay = C0; // C0 from TimsCalibration
244 m_timsVmin = C8; // C8 from TimsCalibration
245 m_timsVmax = C9; // C9 from TimsCalibration
246 m_timsC6 = C6;
247 m_timsC7 = C7;
248
249
251 (C3 - m_timsDvStart) / C1; // //C3 from TimsCalibration // C2 from
252 // TimsCalibration // C1 from TimsCalibration
253}
254double
256{
257 double v = m_timsDvStart +
258 m_timsSlope * ((double)scanNum - m_timsTtrans - m_timsNdelay);
259
260 if(v < m_timsVmin)
261 {
263 QObject::tr("ERROR in TimsFrame::getVoltageTransformation invalid tims "
264 "calibration, v < m_timsVmin"));
265 }
266
267
268 if(v > m_timsVmax)
269 {
271 QObject::tr("ERROR in TimsFrame::getVoltageTransformation invalid tims "
272 "calibration, v > m_timsVmax"));
273 }
274 return v;
275}
276double
277TimsFrameBase::getDriftTime(std::size_t scanNum) const
278{
279 return (m_accumulationTime / (double)m_scanNumber) * ((double)scanNum);
280}
281
282double
284{
285 return 1 / (m_timsC6 + (m_timsC7 / getVoltageTransformation(scanNum)));
286}
287
288
289std::size_t
291{
292 double temp = 1 / one_over_k0;
293 temp = temp - m_timsC6;
294 temp = temp / m_timsC7;
295 temp = 1 / temp;
296 temp = temp - m_timsDvStart;
297 temp = temp / m_timsSlope + m_timsTtrans + m_timsNdelay;
298 return (std::size_t)std::round(temp);
299}
300
301bool
303{
304 if((m_timsDvStart == other.m_timsDvStart) &&
305 (m_timsTtrans == other.m_timsTtrans) &&
306 (m_timsNdelay == other.m_timsNdelay) && (m_timsVmin == other.m_timsVmin) &&
307 (m_timsVmax == other.m_timsVmax) && (m_timsC6 == other.m_timsC6) &&
308 (m_timsC7 == other.m_timsC7) && (m_timsSlope == other.m_timsSlope))
309 {
310 return true;
311 }
312 return false;
313}
314
315
318 std::map<quint32, quint32> &accumulated_scans) const
319{
320 qDebug();
321 // qDebug();
322 // add flanking peaks
323 pappso::Trace local_trace;
324
325 MzCalibrationInterface *mz_calibration_p =
327
328
329 DataPoint element;
330 for(auto &scan_element : accumulated_scans)
331 {
332 // intensity normalization
333 element.y = ((double)scan_element.second) * 100.0 / m_accumulationTime;
334
335 // mz calibration
336 element.x = mz_calibration_p->getMzFromTofIndex(scan_element.first);
337
338 local_trace.push_back(element);
339 }
340 local_trace.sortX();
341
342 qDebug();
343 // qDebug();
344 return local_trace;
345}
346
349 std::map<quint32, quint32> &accumulated_scans) const
350{
351 qDebug();
352 // qDebug();
353 // add flanking peaks
354 std::vector<quint32> keys;
355 transform(begin(accumulated_scans),
356 end(accumulated_scans),
357 back_inserter(keys),
358 [](std::map<quint32, quint32>::value_type const &pair) {
359 return pair.first;
360 });
361 std::sort(keys.begin(), keys.end());
362 pappso::DataPoint data_point_cumul;
363 data_point_cumul.x = 0;
364 data_point_cumul.y = 0;
365
366 pappso::Trace local_trace;
367
368 MzCalibrationInterface *mz_calibration_p =
370
371
372 quint32 last_key = 0;
373
374 for(quint32 key : keys)
375 {
376 if(key == last_key + 1)
377 {
378 // cumulate
379 if(accumulated_scans[key] > accumulated_scans[last_key])
380 {
381 if(data_point_cumul.x == last_key)
382 {
383 // growing peak
384 data_point_cumul.x = key;
385 data_point_cumul.y += accumulated_scans[key];
386 }
387 else
388 {
389 // new peak
390 // flush
391 if(data_point_cumul.y > 0)
392 {
393 // intensity normalization
394 data_point_cumul.y *= 100.0 / m_accumulationTime;
395
396
397 // mz calibration
398 data_point_cumul.x =
399 mz_calibration_p->getMzFromTofIndex(data_point_cumul.x);
400 local_trace.push_back(data_point_cumul);
401 }
402
403 // new point
404 data_point_cumul.x = key;
405 data_point_cumul.y = accumulated_scans[key];
406 }
407 }
408 else
409 {
410 data_point_cumul.y += accumulated_scans[key];
411 }
412 }
413 else
414 {
415 // flush
416 if(data_point_cumul.y > 0)
417 {
418 // intensity normalization
419 data_point_cumul.y *= 100.0 / m_accumulationTime;
420
421
422 // qDebug() << "raw data x=" << data_point_cumul.x;
423 // mz calibration
424 data_point_cumul.x =
425 mz_calibration_p->getMzFromTofIndex(data_point_cumul.x);
426 // qDebug() << "mz=" << data_point_cumul.x;
427 local_trace.push_back(data_point_cumul);
428 }
429
430 // new point
431 data_point_cumul.x = key;
432 data_point_cumul.y = accumulated_scans[key];
433 }
434
435 last_key = key;
436 }
437 // flush
438 if(data_point_cumul.y > 0)
439 {
440 // intensity normalization
441 data_point_cumul.y *= 100.0 / m_accumulationTime;
442
443
444 // mz calibration
445 data_point_cumul.x =
446 mz_calibration_p->getMzFromTofIndex(data_point_cumul.x);
447 local_trace.push_back(data_point_cumul);
448 }
449
450 local_trace.sortX();
451 qDebug();
452 // qDebug();
453 return local_trace;
454}
455
458{
459 if(msp_mzCalibration == nullptr)
460 {
461
463 QObject::tr("ERROR in %1, %2, %3 msp_mzCalibration is null")
464 .arg(__FILE__)
465 .arg(__FUNCTION__)
466 .arg(__LINE__));
467 }
468 return msp_mzCalibration;
469}
470
471void
473 MzCalibrationInterfaceSPtr mzCalibration)
474{
475
476 if(mzCalibration == nullptr)
477 {
478
480 QObject::tr("ERROR in %1, %2, %3 msp_mzCalibration is null")
481 .arg(__FILE__)
482 .arg(__FUNCTION__)
483 .arg(__LINE__));
484 }
485 msp_mzCalibration = mzCalibration;
486}
487
488
489quint32
491{
492 quint32 max_value = 0;
493 for(quint32 i = 0; i < m_scanNumber; i++)
494 {
495 qDebug() << "m_scanNumber=" << m_scanNumber << " i=" << i;
496 std::vector<quint32> index_list = getScanIndexList(i);
497 auto it = std::max_element(index_list.begin(), index_list.end());
498 if(it != index_list.end())
499 {
500 max_value = std::max(max_value, *it);
501 }
502 }
503 return max_value;
504}
505
506std::vector<quint32>
507TimsFrameBase::getScanIndexList(std::size_t scanNum) const
508{
509 throw PappsoException(
510 QObject::tr(
511 "ERROR unable to getScanIndexList in TimsFrameBase for scan number %1")
512 .arg(scanNum));
513}
514
515
516std::vector<quint32>
517TimsFrameBase::getScanIntensities(std::size_t scanNum) const
518{
519 throw PappsoException(
520 QObject::tr(
521 "ERROR unable to getScanIntensities in TimsFrameBase for scan number %1")
522 .arg(scanNum));
523}
524
525Trace
527 std::size_t mz_index_lower_bound,
528 std::size_t mz_index_upper_bound,
529 XicExtractMethod method) const
530{
531 Trace im_trace;
532 DataPoint data_point;
533 for(quint32 i = 0; i < m_scanNumber; i++)
534 {
535 data_point.x = i;
536 data_point.y = 0;
537 qDebug() << "m_scanNumber=" << m_scanNumber << " i=" << i;
538 std::vector<quint32> index_list = getScanIndexList(i);
539 auto it_lower = std::find_if(index_list.begin(),
540 index_list.end(),
541 [mz_index_lower_bound](quint32 to_compare) {
542 if(to_compare < mz_index_lower_bound)
543 {
544 return false;
545 }
546 return true;
547 });
548
549
550 if(it_lower == index_list.end())
551 {
552 }
553 else
554 {
555
556
557 auto it_upper =
558 std::find_if(index_list.begin(),
559 index_list.end(),
560 [mz_index_upper_bound](quint32 to_compare) {
561 if(mz_index_upper_bound >= to_compare)
562 {
563 return false;
564 }
565 return true;
566 });
567 std::vector<quint32> intensity_list = getScanIntensities(i);
568 for(int j = std::distance(index_list.begin(), it_lower);
569 j < std::distance(index_list.begin(), it_upper);
570 j++)
571 {
572 if(method == XicExtractMethod::sum)
573 {
574 data_point.y += intensity_list[j];
575 }
576 else
577 {
578 data_point.y =
579 std::max((double)intensity_list[j], data_point.y);
580 }
581 }
582 }
583 im_trace.push_back(data_point);
584 }
585 qDebug();
586 return im_trace;
587}
588// namespace pappso
589} // namespace pappso
virtual double getMzFromTofIndex(quint32 tof_index)=0
get m/z from time of flight raw index
double m_accumulationTime
accumulation time in milliseconds
double getTime() const
virtual quint64 cumulateSingleScanIntensities(std::size_t scanNum) const
double getVoltageTransformation(std::size_t scanNum) const
get voltage for a given scan number
virtual Trace cumulateScanToTrace(std::size_t scanNumBegin, std::size_t scanNumEnd) const
cumulate spectrum given a scan number range need the binary file The intensities are normalized with ...
virtual std::size_t getTotalNumberOfScans() const
get the number of scans contained in this frame each scan represents an ion mobility slice
virtual Trace getIonMobilityTraceByMzIndexRange(std::size_t mz_index_lower_bound, std::size_t mz_index_upper_bound, XicExtractMethod method) const
get a mobility trace cumulating intensities inside the given mass index range
virtual std::size_t getNbrPeaks(std::size_t scanNum) const
get the number of peaks in this spectrum need the binary file
MzCalibrationInterfaceSPtr msp_mzCalibration
virtual MassSpectrumSPtr getMassSpectrumSPtr(std::size_t scanNum) const
get Mass spectrum with peaks for this scan number need the binary file
TimsFrameBase(std::size_t timsId, quint32 scanNum)
constructor for binary independant tims frame
virtual std::vector< quint32 > getScanIndexList(std::size_t scanNum) const
get raw index list for one given scan index are not TOF nor m/z, just index on digitizer
double getDriftTime(std::size_t scanNum) const
get drift time of a scan number in milliseconds
pappso::Trace getTraceFromCumulatedScansBuiltinCentroid(std::map< quint32, quint32 > &accumulated_scans) const
transform accumulation of raw scans into a real mass spectrum with a simple centroid on raw integers
double m_time
retention time
virtual quint64 cumulateScansIntensities(std::size_t scanNumBegin, std::size_t scanNumEnd) const
void setAccumulationTime(double accumulation_time_ms)
quint32 m_scanNumber
total number of scans contained in this frame
virtual void cumulateScansInRawMap(std::map< quint32, quint32 > &rawSpectrum, std::size_t scanNumBegin, std::size_t scanNumEnd) const
cumulate scan list into a trace into a raw spectrum map The intensities are NOT normalized with respe...
void setTime(double time)
std::size_t m_timsId
Tims frame database id (the SQL identifier of this frame)
pappso::Trace getTraceFromCumulatedScans(std::map< quint32, quint32 > &accumulated_scans) const
transform accumulation of raw scans into a real mass spectrum
virtual bool hasSameCalibrationData(const TimsFrameBase &other) const
tells if 2 tims frame has the same calibration data Usefull to know if raw data can be handled betwee...
virtual quint32 getMaximumRawMassIndex() const
get the maximum raw mass index contained in this frame
unsigned int getMsLevel() const
void setTimsCalibration(int tims_model_type, double C0, double C1, double C2, double C3, double C4, double C5, double C6, double C7, double C8, double C9)
virtual const MzCalibrationInterfaceSPtr & getMzCalibrationInterfaceSPtr() const final
get the MzCalibration model to compute mz and TOF for this frame
std::size_t getScanNumFromOneOverK0(double one_over_k0) const
get the scan number from a given 1/Ko mobility value
void setMsMsType(quint8 type)
void setMzCalibration(double T1_frame, double T2_frame, double digitizerTimebase, double digitizerDelay, double C0, double C1, double C2, double C3, double C4, double T1_ref, double T2_ref, double dC1, double dC2)
double getOneOverK0Transformation(std::size_t scanNum) const
get 1/K0 value of a given scan (mobility value)
bool checkScanNum(std::size_t scanNum) const
check that this scan number exists
void setMzCalibrationInterfaceSPtr(MzCalibrationInterfaceSPtr mzCalibration)
std::size_t getId() const
virtual std::vector< quint32 > getScanIntensities(std::size_t scanNum) const
get raw intensities without transformation from one scan it needs intensity normalization
A simple container of DataPoint instances.
Definition: trace.h:38
Q_INVOKABLE void sortX()
Definition: trace.cpp:983
implement Bruker's model type 1 formula to compute m/z
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition: aa.cpp:39
std::shared_ptr< MzCalibrationInterface > MzCalibrationInterfaceSPtr
std::shared_ptr< MassSpectrum > MassSpectrumSPtr
Definition: massspectrum.h:54
XicExtractMethod
Definition: types.h:201
@ sum
sum of intensities
pappso_double x
Definition: datapoint.h:23
pappso_double y
Definition: datapoint.h:24
handle a single Bruker's TimsTof frame without binary data