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}
149
150void
151TimsFrameBase::cumulateScansInRawMap(std::map<quint32, quint32> &rawSpectrum
152 [[maybe_unused]],
153 std::size_t scanNumBegin,
154 std::size_t scanNumEnd) const
155{
156 throw PappsoException(
157 QObject::tr(
158 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
159 "number begin %1 end %2")
160 .arg(scanNumBegin)
161 .arg(scanNumEnd));
162}
163
164
165quint64
167{
168 throw PappsoException(
169 QObject::tr(
170 "ERROR unable to cumulateSingleScanIntensities in TimsFrameBase for scan "
171 "number %1.").arg(scanNum));
172
173 return 0;
174}
175
176
177quint64
179 std::size_t scanNumEnd) const
180{
181 throw PappsoException(
182 QObject::tr(
183 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
184 "number begin %1 end %2")
185 .arg(scanNumBegin)
186 .arg(scanNumEnd));
187
188 return 0;
189}
190
191void
193{
194 m_time = time;
195}
196
197void
199{
200
201 qDebug() << " m_msMsType=" << type;
202 m_msMsType = type;
203}
204
205unsigned int
207{
208 if(m_msMsType == 0)
209 return 1;
210 return 2;
211}
212
213double
215{
216 return m_time;
217}
218
219std::size_t
221{
222 return m_timsId;
223}
224void
226 double C0,
227 double C1,
228 double C2,
229 double C3,
230 double C4,
231 [[maybe_unused]] double C5,
232 double C6,
233 double C7,
234 double C8,
235 double C9)
236{
237 if(tims_model_type != 2)
238 {
239 throw pappso::PappsoException(QObject::tr(
240 "ERROR in TimsFrame::setTimsCalibration tims_model_type != 2"));
241 }
242 m_timsDvStart = C2; // C2 from TimsCalibration
243 m_timsTtrans = C4; // C4 from TimsCalibration
244 m_timsNdelay = C0; // C0 from TimsCalibration
245 m_timsVmin = C8; // C8 from TimsCalibration
246 m_timsVmax = C9; // C9 from TimsCalibration
247 m_timsC6 = C6;
248 m_timsC7 = C7;
249
250
252 (C3 - m_timsDvStart) / C1; // //C3 from TimsCalibration // C2 from
253 // TimsCalibration // C1 from TimsCalibration
254}
255double
257{
258 double v = m_timsDvStart +
259 m_timsSlope * ((double)scanNum - m_timsTtrans - m_timsNdelay);
260
261 if(v < m_timsVmin)
262 {
264 QObject::tr("ERROR in TimsFrame::getVoltageTransformation invalid tims "
265 "calibration, v < m_timsVmin"));
266 }
267
268
269 if(v > m_timsVmax)
270 {
272 QObject::tr("ERROR in TimsFrame::getVoltageTransformation invalid tims "
273 "calibration, v > m_timsVmax"));
274 }
275 return v;
276}
277double
278TimsFrameBase::getDriftTime(std::size_t scanNum) const
279{
280 return (m_accumulationTime / (double)m_scanNumber) * ((double)scanNum);
281}
282
283double
285{
286 return 1 / (m_timsC6 + (m_timsC7 / getVoltageTransformation(scanNum)));
287}
288
289
290std::size_t
292{
293 double temp = 1 / one_over_k0;
294 temp = temp - m_timsC6;
295 temp = temp / m_timsC7;
296 temp = 1 / temp;
297 temp = temp - m_timsDvStart;
298 temp = temp / m_timsSlope + m_timsTtrans + m_timsNdelay;
299 return (std::size_t)std::round(temp);
300}
301
302bool
304{
305 if((m_timsDvStart == other.m_timsDvStart) &&
306 (m_timsTtrans == other.m_timsTtrans) &&
307 (m_timsNdelay == other.m_timsNdelay) && (m_timsVmin == other.m_timsVmin) &&
308 (m_timsVmax == other.m_timsVmax) && (m_timsC6 == other.m_timsC6) &&
309 (m_timsC7 == other.m_timsC7) && (m_timsSlope == other.m_timsSlope))
310 {
311 return true;
312 }
313 return false;
314}
315
316
319 std::map<quint32, quint32> &accumulated_scans) const
320{
321 qDebug();
322 // qDebug();
323 // add flanking peaks
324 pappso::Trace local_trace;
325
326 MzCalibrationInterface *mz_calibration_p =
328
329
330 DataPoint element;
331 for(auto &scan_element : accumulated_scans)
332 {
333 // intensity normalization
334 element.y = ((double)scan_element.second) * 100.0 / m_accumulationTime;
335
336 // mz calibration
337 element.x = mz_calibration_p->getMzFromTofIndex(scan_element.first);
338
339 local_trace.push_back(element);
340 }
341 local_trace.sortX();
342
343 qDebug();
344 // qDebug();
345 return local_trace;
346}
347
350 std::map<quint32, quint32> &accumulated_scans) const
351{
352 qDebug();
353 // qDebug();
354 // add flanking peaks
355 std::vector<quint32> keys;
356 transform(begin(accumulated_scans),
357 end(accumulated_scans),
358 back_inserter(keys),
359 [](std::map<quint32, quint32>::value_type const &pair) {
360 return pair.first;
361 });
362 std::sort(keys.begin(), keys.end());
363 pappso::DataPoint data_point_cumul;
364 data_point_cumul.x = 0;
365 data_point_cumul.y = 0;
366
367 pappso::Trace local_trace;
368
369 MzCalibrationInterface *mz_calibration_p =
371
372
373 quint32 last_key = 0;
374
375 for(quint32 key : keys)
376 {
377 if(key == last_key + 1)
378 {
379 // cumulate
380 if(accumulated_scans[key] > accumulated_scans[last_key])
381 {
382 if(data_point_cumul.x == last_key)
383 {
384 // growing peak
385 data_point_cumul.x = key;
386 data_point_cumul.y += accumulated_scans[key];
387 }
388 else
389 {
390 // new peak
391 // flush
392 if(data_point_cumul.y > 0)
393 {
394 // intensity normalization
395 data_point_cumul.y *= 100.0 / m_accumulationTime;
396
397
398 // mz calibration
399 data_point_cumul.x =
400 mz_calibration_p->getMzFromTofIndex(data_point_cumul.x);
401 local_trace.push_back(data_point_cumul);
402 }
403
404 // new point
405 data_point_cumul.x = key;
406 data_point_cumul.y = accumulated_scans[key];
407 }
408 }
409 else
410 {
411 data_point_cumul.y += accumulated_scans[key];
412 }
413 }
414 else
415 {
416 // flush
417 if(data_point_cumul.y > 0)
418 {
419 // intensity normalization
420 data_point_cumul.y *= 100.0 / m_accumulationTime;
421
422
423 // qDebug() << "raw data x=" << data_point_cumul.x;
424 // mz calibration
425 data_point_cumul.x =
426 mz_calibration_p->getMzFromTofIndex(data_point_cumul.x);
427 // qDebug() << "mz=" << data_point_cumul.x;
428 local_trace.push_back(data_point_cumul);
429 }
430
431 // new point
432 data_point_cumul.x = key;
433 data_point_cumul.y = accumulated_scans[key];
434 }
435
436 last_key = key;
437 }
438 // flush
439 if(data_point_cumul.y > 0)
440 {
441 // intensity normalization
442 data_point_cumul.y *= 100.0 / m_accumulationTime;
443
444
445 // mz calibration
446 data_point_cumul.x =
447 mz_calibration_p->getMzFromTofIndex(data_point_cumul.x);
448 local_trace.push_back(data_point_cumul);
449 }
450
451 local_trace.sortX();
452 qDebug();
453 // qDebug();
454 return local_trace;
455}
456
459{
460 if(msp_mzCalibration == nullptr)
461 {
462
464 QObject::tr("ERROR in %1, %2, %3 msp_mzCalibration is null")
465 .arg(__FILE__)
466 .arg(__FUNCTION__)
467 .arg(__LINE__));
468 }
469 return msp_mzCalibration;
470}
471
472void
474 MzCalibrationInterfaceSPtr mzCalibration)
475{
476
477 if(mzCalibration == nullptr)
478 {
479
481 QObject::tr("ERROR in %1, %2, %3 msp_mzCalibration is null")
482 .arg(__FILE__)
483 .arg(__FUNCTION__)
484 .arg(__LINE__));
485 }
486 msp_mzCalibration = mzCalibration;
487}
488
489
490quint32
492{
493 quint32 max_value = 0;
494 for(quint32 i = 0; i < m_scanNumber; i++)
495 {
496 qDebug() << "m_scanNumber=" << m_scanNumber << " i=" << i;
497 std::vector<quint32> index_list = getScanIndexList(i);
498 auto it = std::max_element(index_list.begin(), index_list.end());
499 if(it != index_list.end())
500 {
501 max_value = std::max(max_value, *it);
502 }
503 }
504 return max_value;
505}
506
507std::vector<quint32>
508TimsFrameBase::getScanIndexList(std::size_t scanNum) const
509{
510 throw PappsoException(
511 QObject::tr(
512 "ERROR unable to getScanIndexList in TimsFrameBase for scan number %1")
513 .arg(scanNum));
514}
515
516
517std::vector<quint32>
518TimsFrameBase::getScanIntensities(std::size_t scanNum) const
519{
520 throw PappsoException(
521 QObject::tr(
522 "ERROR unable to getScanIntensities in TimsFrameBase for scan number %1")
523 .arg(scanNum));
524}
525
526Trace
528 std::size_t mz_index_lower_bound,
529 std::size_t mz_index_upper_bound,
530 XicExtractMethod method) const
531{
532 Trace im_trace;
533 DataPoint data_point;
534 for(quint32 i = 0; i < m_scanNumber; i++)
535 {
536 data_point.x = i;
537 data_point.y = 0;
538 qDebug() << "m_scanNumber=" << m_scanNumber << " i=" << i;
539 std::vector<quint32> index_list = getScanIndexList(i);
540 auto it_lower = std::find_if(index_list.begin(),
541 index_list.end(),
542 [mz_index_lower_bound](quint32 to_compare) {
543 if(to_compare < mz_index_lower_bound)
544 {
545 return false;
546 }
547 return true;
548 });
549
550
551 if(it_lower == index_list.end())
552 {
553 }
554 else
555 {
556
557
558 auto it_upper =
559 std::find_if(index_list.begin(),
560 index_list.end(),
561 [mz_index_upper_bound](quint32 to_compare) {
562 if(mz_index_upper_bound >= to_compare)
563 {
564 return false;
565 }
566 return true;
567 });
568 std::vector<quint32> intensity_list = getScanIntensities(i);
569 for(int j = std::distance(index_list.begin(), it_lower);
570 j < std::distance(index_list.begin(), it_upper);
571 j++)
572 {
573 if(method == XicExtractMethod::sum)
574 {
575 data_point.y += intensity_list[j];
576 }
577 else
578 {
579 data_point.y =
580 std::max((double)intensity_list[j], data_point.y);
581 }
582 }
583 }
584 im_trace.push_back(data_point);
585 }
586 qDebug();
587 return im_trace;
588}
589// namespace pappso
590} // 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:148
void sortX()
Definition: trace.cpp:987
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