libpappsomspp
Library for mass spectrometry
Loading...
Searching...
No Matches
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"
31#include <QDebug>
32#include <QObject>
33#include <cmath>
34#include <algorithm>
35
36namespace pappso
37{
38
39TimsFrameBase::TimsFrameBase(std::size_t frame_id, quint32 scanCount)
40{
41 qDebug() << frame_id;
42 m_frameId = frame_id;
43
44 m_scanCount = scanCount;
45}
46
47TimsFrameBase::TimsFrameBase([[maybe_unused]] const TimsFrameBase &other)
48{
49}
50
54
55void
56TimsFrameBase::setAcqDurationInMilliseconds(double acquisition_duration_ms)
57{
58 m_acqDurationInMilliseconds = acquisition_duration_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_scanCount)
106 {
108 QObject::tr("Invalid scan number : scanNum %1 > m_scanNumber %2")
109 .arg(scanNum)
110 .arg(m_scanCount));
111 }
112
113 return true;
114}
115
116std::size_t
117TimsFrameBase::getScanPeakCount(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
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}
139
140
143{
144 // qDebug();
145
146 return getMassSpectrumSPtr(scanNum);
147}
148Trace
149TimsFrameBase::cumulateScansToTrace(std::size_t scanNumBegin,
150 std::size_t scanNumEnd) const
151{
152 throw PappsoException(
153 QObject::tr("ERROR unable to cumulateScanToTrace in TimsFrameBase for scan "
154 "number begin %1 end %2")
155 .arg(scanNumBegin)
156 .arg(scanNumEnd));
157}
158
159Trace
161 std::size_t mzindex_merge_window [[maybe_unused]],
162 std::size_t scanNumBegin [[maybe_unused]],
163 std::size_t scanNumEnd [[maybe_unused]],
164 quint32 &minimum_tof_index_out [[maybe_unused]],
165 quint32 &maximum_tof_index_out [[maybe_unused]]) const
166{
167 throw PappsoException(QObject::tr("Non implemented function %1 %2 %3")
168 .arg(__FILE__)
169 .arg(__FUNCTION__)
170 .arg(__LINE__));
171}
172
173Trace
175 std::size_t mz_index_merge_window [[maybe_unused]],
176 double mz_range_begin [[maybe_unused]],
177 double mz_range_end [[maybe_unused]],
178 std::size_t mobility_scan_begin [[maybe_unused]],
179 std::size_t mobility_scan_end [[maybe_unused]],
180 quint32 &mz_minimum_index_out [[maybe_unused]],
181 quint32 &mz_maximum_index_out [[maybe_unused]]) const
182{
183 throw PappsoException(QObject::tr("Non implemented function %1 %2 %3")
184 .arg(__FILE__)
185 .arg(__FUNCTION__)
186 .arg(__LINE__));
187}
188
189Trace
190TimsFrameBase::getMobilityScan(std::size_t scanNum [[maybe_unused]],
191 std::size_t mz_index_merge_window
192 [[maybe_unused]],
193 double mz_range_begin [[maybe_unused]],
194 double mz_range_end [[maybe_unused]],
195 quint32 &mz_minimum_index_out [[maybe_unused]],
196 quint32 &mz_maximum_index_out
197 [[maybe_unused]]) const
198{
199 throw PappsoException(QObject::tr("Non implemented function %1 %2 %3")
200 .arg(__FILE__)
201 .arg(__FUNCTION__)
202 .arg(__LINE__));
203}
204
205void
207 [[maybe_unused]],
208 std::size_t scanNumBegin,
209 std::size_t scanNumEnd) const
210{
211 throw PappsoException(
212 QObject::tr(
213 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
214 "number begin %1 end %2")
215 .arg(scanNumBegin)
216 .arg(scanNumEnd));
217}
218
219
220void
222 [[maybe_unused]],
223 std::size_t scanNumBegin,
224 std::size_t scanNumEnd,
225 quint32 tof_index_begin,
226 quint32 tof_index_end) const
227{
228 throw PappsoException(
229 QObject::tr(
230 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
231 "number begin %1 end %2, tof index begin %3, tof index end %4")
232 .arg(scanNumBegin)
233 .arg(scanNumEnd)
234 .arg(tof_index_begin)
235 .arg(tof_index_end));
236}
237
238quint64
240{
241 throw PappsoException(
242 QObject::tr(
243 "ERROR unable to cumulateSingleScanIntensities in TimsFrameBase for scan "
244 "number %1.")
245 .arg(scanNum));
246
247 return 0;
248}
249
250
251quint64
253 std::size_t scanNumEnd) const
254{
255 throw PappsoException(
256 QObject::tr(
257 "ERROR unable to cumulateScansInRawMap in TimsFrameBase for scan "
258 "number begin %1 end %2")
259 .arg(scanNumBegin)
260 .arg(scanNumEnd));
261
262 return 0;
263}
264
265void
267{
268 m_rtInSeconds = time;
269}
270
271void
273{
274
275 qDebug() << " m_msMsType=" << type;
276 m_msMsType = type;
277}
278
279unsigned int
281{
282 if(m_msMsType == 0)
283 return 1;
284 return 2;
285}
286
287double
289{
290 return m_rtInSeconds;
291}
292
293std::size_t
295{
296 return m_frameId;
297}
298void
300 double C0,
301 double C1,
302 double C2,
303 double C3,
304 double C4,
305 [[maybe_unused]] double C5,
306 double C6,
307 double C7,
308 double C8,
309 double C9)
310{
311 if(tims_model_type != 2)
312 {
313 throw pappso::PappsoException(QObject::tr(
314 "ERROR in TimsFrame::setTimsCalibration tims_model_type != 2"));
315 }
316 m_timsDvStart = C2; // C2 from TimsCalibration
317 m_timsTtrans = C4; // C4 from TimsCalibration
318 m_timsNdelay = C0; // C0 from TimsCalibration
319 m_timsVmin = C8; // C8 from TimsCalibration
320 m_timsVmax = C9; // C9 from TimsCalibration
321 m_timsC6 = C6;
322 m_timsC7 = C7;
323
324
326 (C3 - m_timsDvStart) / C1; // //C3 from TimsCalibration // C2 from
327 // TimsCalibration // C1 from TimsCalibration
328}
329double
331{
332 double v = m_timsDvStart +
333 m_timsSlope * ((double)scanNum - m_timsTtrans - m_timsNdelay);
334
335 if(v < m_timsVmin)
336 {
338 QObject::tr("ERROR in TimsFrame::getVoltageTransformation invalid tims "
339 "calibration, v < m_timsVmin %1 < %2")
340 .arg(v)
341 .arg(m_timsVmin));
342 }
343
344
345 if(v > m_timsVmax)
346 {
348 QObject::tr("ERROR in TimsFrame::getVoltageTransformation invalid tims "
349 "calibration, v > m_timsVmax %1 > %2")
350 .arg(v)
351 .arg(m_timsVmax));
352 }
353 return v;
354}
355double
357{
358 return (m_acqDurationInMilliseconds / (double)m_scanCount) *
359 ((double)scanNum);
360}
361
362double
364{
365 return 1 / (m_timsC6 + (m_timsC7 / getVoltageTransformation(scanNum)));
366}
367
368
369std::size_t
371{
372 double temp = 1 / one_over_k0;
373 temp = temp - m_timsC6;
374 temp = temp / m_timsC7;
375 temp = 1 / temp;
376 temp = temp - m_timsDvStart;
377 temp = temp / m_timsSlope + m_timsTtrans + m_timsNdelay;
378 return (std::size_t)std::round(temp);
379}
380
381bool
383{
384 if((m_timsDvStart == other.m_timsDvStart) &&
385 (m_timsTtrans == other.m_timsTtrans) &&
386 (m_timsNdelay == other.m_timsNdelay) && (m_timsVmin == other.m_timsVmin) &&
387 (m_timsVmax == other.m_timsVmax) && (m_timsC6 == other.m_timsC6) &&
388 (m_timsC7 == other.m_timsC7) && (m_timsSlope == other.m_timsSlope))
389 {
390 return true;
391 }
392 return false;
393}
394
395
398 TimsDataFastMap &accumulated_scans) const
399{
400 qDebug();
401 // qDebug();
402 // add flanking peaks
403 pappso::Trace local_trace;
404
405 MzCalibrationInterface *mz_calibration_p =
407
408
409 DataPoint element;
410 for(quint32 tof_index : accumulated_scans.getTofIndexList())
411 {
412 // intensity normalization
413 element.y = ((double)accumulated_scans.readIntensity(tof_index)) * 100.0 /
415
416 // mz calibration
417 element.x = mz_calibration_p->getMzFromTofIndex(tof_index);
418
419 local_trace.push_back(element);
420 }
421 local_trace.sortX();
422
423 qDebug();
424 // qDebug();
425 return local_trace;
426}
427
430{
431 if(msp_mzCalibration == nullptr)
432 {
433
435 QObject::tr("ERROR in %1, %2, %3 msp_mzCalibration is null")
436 .arg(__FILE__)
437 .arg(__FUNCTION__)
438 .arg(__LINE__));
439 }
440 return msp_mzCalibration;
441}
442
443void
445 MzCalibrationInterfaceSPtr mzCalibration)
446{
447
448 if(mzCalibration == nullptr)
449 {
450
452 QObject::tr("ERROR in %1, %2, %3 msp_mzCalibration is null")
453 .arg(__FILE__)
454 .arg(__FUNCTION__)
455 .arg(__LINE__));
456 }
457 msp_mzCalibration = mzCalibration;
458}
459
460
461quint32
463{
464 quint32 max_value = 0;
465 for(quint32 i = 0; i < m_scanCount; i++)
466 {
467 qDebug() << "m_scanNumber=" << m_scanCount << " i=" << i;
468 std::vector<quint32> index_list = getScanTofIndexList(i);
469 auto it = std::max_element(index_list.begin(), index_list.end());
470 if(it != index_list.end())
471 {
472 max_value = std::max(max_value, *it);
473 }
474 }
475 return max_value;
476}
477
478std::vector<quint32>
479TimsFrameBase::getScanTofIndexList(std::size_t scanNum) const
480{
481 throw PappsoException(
482 QObject::tr(
483 "ERROR unable to getScanIndexList in TimsFrameBase for scan number %1")
484 .arg(scanNum));
485}
486
487
488std::vector<quint32>
489TimsFrameBase::getScanIntensityList(std::size_t scanNum) const
490{
491 throw PappsoException(
492 QObject::tr(
493 "ERROR unable to getScanIntensities in TimsFrameBase for scan number %1")
494 .arg(scanNum));
495}
496
497Trace
499 std::size_t tof_index_begin,
500 std::size_t tof_index_end,
501 XicExtractMethod method,
502 std::size_t scan_index_begin,
503 std::size_t scan_index_end) const
504{
505 Trace im_trace;
506 DataPoint data_point;
507 if(scan_index_end > m_scanCount)
508 scan_index_end = m_scanCount;
509 for(quint32 i = scan_index_begin; i < (scan_index_end + 1); i++)
510 {
511 data_point.x = i;
512 data_point.y = 0;
513 qDebug() << "m_scanNumber=" << m_scanCount << " i=" << i;
514 std::vector<quint32> index_list = getScanTofIndexList(i);
515 auto it_lower = std::find_if(index_list.begin(),
516 index_list.end(),
517 [tof_index_begin](quint32 to_compare) {
518 if(to_compare < tof_index_begin)
519 {
520 return false;
521 }
522 return true;
523 });
524
525
526 if(it_lower == index_list.end())
527 {
528 }
529 else
530 {
531
532
533 auto it_upper = std::find_if(index_list.begin(),
534 index_list.end(),
535 [tof_index_end](quint32 to_compare) {
536 if(tof_index_end >= to_compare)
537 {
538 return false;
539 }
540 return true;
541 });
542 std::vector<quint32> intensity_list = getScanIntensityList(i);
543 for(int j = std::distance(index_list.begin(), it_lower);
544 j < std::distance(index_list.begin(), it_upper);
545 j++)
546 {
547 if(method == XicExtractMethod::sum)
548 {
549 data_point.y += intensity_list[j];
550 }
551 else
552 {
553 data_point.y =
554 std::max((double)intensity_list[j], data_point.y);
555 }
556 }
557 }
558 im_trace.push_back(data_point);
559 }
560 qDebug();
561 return im_trace;
562}
563
564
565std::vector<TimsFrameBase::TofIndexIntensityPair> &
566TimsFrameBase::downgradeResolutionOfTofIndexIntensityPairList(
567 std::size_t mzindex_merge_window,
568 std::vector<TimsFrameBase::TofIndexIntensityPair> &rawSpectrum) const
569{
570
571 qDebug() << rawSpectrum.size();
572 std::vector<TimsFrameBase::TofIndexIntensityPair> new_spectrum;
573
575 current_point.intensity_index = 0;
576 current_point.tof_index = 0;
577 for(auto &pair_mz_intensity : rawSpectrum)
578 {
579 quint32 mzkey = (pair_mz_intensity.tof_index / mzindex_merge_window);
580 mzkey = (mzkey * mzindex_merge_window) + (mzindex_merge_window / 2);
581
582 if(current_point.tof_index != mzkey)
583 {
584 if(current_point.tof_index > 0)
585 {
586 new_spectrum.push_back(current_point);
587 }
588
589 current_point.intensity_index = pair_mz_intensity.intensity_index;
590 current_point.tof_index = mzkey;
591 }
592 else
593 {
594 current_point.intensity_index += pair_mz_intensity.intensity_index;
595 }
596 }
597
598 if(current_point.tof_index > 0)
599 {
600 new_spectrum.push_back(current_point);
601 }
602 rawSpectrum = new_spectrum;
603 qDebug() << rawSpectrum.size();
604 return rawSpectrum;
605}
606// namespace pappso
607} // namespace pappso
virtual double getMzFromTofIndex(quint32 tof_index)=0
get m/z from time of flight raw index
replacement for std::map
const std::vector< quint32 > & getTofIndexList() const
std::size_t readIntensity(quint32)
reads intensity for a tof_index
double m_rtInSeconds
retention time
double getVoltageTransformation(std::size_t scanNum) const
get voltage for a given scan number
virtual std::size_t getTotalNumberOfScans() const
get the number of scans contained in this frame each scan represents an ion mobility slice
virtual pappso::MassSpectrumCstSPtr getMassSpectrumCstSPtr(std::size_t scan_index) const final
get the mass spectrum corresponding to a scan index
MzCalibrationInterfaceSPtr msp_mzCalibration
virtual std::size_t getScanPeakCount(std::size_t scanIndex) const
get the number of peaks in this spectrum need the binary file
double getDriftTimeInMilliseconds(std::size_t scan_index) const
get drift time of a scan number in milliseconds
virtual Trace cumulateScansToTrace(std::size_t scanIndexBegin, std::size_t scanIndexEnd) const
cumulate spectrum given a scan number range need the binary file The intensities are normalized with ...
std::size_t getScanIndexFromOneOverK0(double one_over_k0) const
get the scan number from a given 1/Ko mobility value
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 Trace combineScansToTraceWithDowngradedMzResolution(std::size_t tof_index_merge_window, std::size_t scanIndexBegin, std::size_t scanIndexEnd, quint32 &minimum_tof_index_out, quint32 &maximum_tof_index_out) const
cumulate spectrum given a scan index range need the binary file The intensities are normalized with r...
virtual Trace getIonMobilityTraceByTofIndexRange(std::size_t tof_index_begin, std::size_t tof_index_end, XicExtractMethod method, std::size_t scan_index_begin, std::size_t scan_index_end) const
get a mobility trace cumulating intensities inside the given mass index range
double m_acqDurationInMilliseconds
acquisition duration in milliseconds
double getRtInSeconds() const
virtual quint32 getMaximumRawMassIndex() const
get the maximum raw mass index contained in this frame
TimsFrameBase(std::size_t frameId, quint32 scanCount)
constructor for binary independant tims frame
pappso::Trace getTraceFromTofIndexIntensityMap(TimsDataFastMap &accumulated_scans) const
transform accumulation of raw scans into a real mass spectrum
virtual Trace combineScansToTraceWithDowngradedMzResolution2(std::size_t mz_index_merge_window, double mz_range_begin, double mz_range_end, std::size_t mobility_scan_begin, std::size_t mobility_scan_end, quint32 &mz_minimum_index_out, quint32 &mz_maximum_index_out) const
cumulate spectrum given a scan number range need the binary file The intensities are normalized with ...
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 MassSpectrumSPtr getMassSpectrumSPtr(std::size_t scan_index) const
get Mass spectrum with peaks for this scan index need the binary file
virtual const MzCalibrationInterfaceSPtr & getMzCalibrationInterfaceSPtr() const final
get the MzCalibration model to compute mz and TOF for this frame
virtual quint64 cumulateScanRangeIntensities(std::size_t scan_index_begin, std::size_t scan_index_end) const
void setRtInSeconds(double time)
double getOneOverK0Transformation(std::size_t scan_index) const
get 1/K0 value of a given scan (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)
quint32 m_scanCount
total number of scans contained in this frame
std::size_t m_frameId
Tims frame database id (the SQL identifier of this frame)
void setAcqDurationInMilliseconds(double acquisition_duration_ms)
virtual Trace getMobilityScan(std::size_t scan_index, std::size_t tof_index_merge_window, double mz_range_begin, double mz_range_end, quint32 &mz_minimum_index_out, quint32 &mz_maximum_index_out) const
get a single mobility scan m/z + intensities
bool checkScanNum(std::size_t scanNum) const
check that this scan number exists
virtual void combineScansInTofIndexIntensityMap(TimsDataFastMap &tof_index_intensity_map, std::size_t scan_index_begin, std::size_t scan_index_end) const
cumulate scan list into a trace into a raw spectrum map The intensities are NOT normalized with respe...
void setMzCalibrationInterfaceSPtr(MzCalibrationInterfaceSPtr mzCalibration)
virtual std::vector< quint32 > getScanTofIndexList(std::size_t scan_index) const
get raw index list for one given scan index are not TOF nor m/z, just index on digitizer
virtual std::vector< quint32 > getScanIntensityList(std::size_t scan_index) const
get raw intensities without transformation from one scan it needs intensity normalization
std::size_t getId() const
virtual quint64 cumulateScanIntensities(std::size_t scan_index) const
A simple container of DataPoint instances.
Definition trace.h:148
void sortX(SortOrder sort_order=SortOrder::ascending)
Definition trace.cpp:1086
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< const MassSpectrum > MassSpectrumCstSPtr
std::shared_ptr< MassSpectrum > MassSpectrumSPtr
XicExtractMethod
Definition types.h:247
@ 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