libpappsomspp
Library for mass spectrometry
Loading...
Searching...
No Matches
timsframesmsrunreader.cpp
Go to the documentation of this file.
1/**
2 * \file pappsomspp/msrun/private/timsmsrunreader.h
3 * \date 05/09/2019
4 * \author Olivier Langella
5 * \brief MSrun file reader for native Bruker TimsTOF raw 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
29#include "../../exception/exceptionnotimplemented.h"
30#include "../../exception/exceptioninterrupted.h"
31#include <QDebug>
32
33using namespace pappso;
34
36 : MsRunReader(msrun_id_csp)
37{
38 qDebug() << "Now initializing the TimsFramesMsRunReader.";
39
40 initialize();
41}
42
47
48void
50{
51 msp_timsData = std::make_shared<TimsData>(mcsp_msRunId.get()->getFileName());
52 if(msp_timsData == nullptr)
53 {
54 throw PappsoException(
55 QObject::tr("ERROR in TimsFramesMsRunReader::initialize "
56 "msp_timsData is null for MsRunId %1")
57 .arg(mcsp_msRunId.get()->toString()));
58 }
59}
60
61
62bool
63TimsFramesMsRunReader::accept(const QString &file_name) const
64{
65 qDebug() << file_name;
66 return true;
67}
68
69
72 [[maybe_unused]] std::size_t spectrum_index)
73{
75 QObject::tr("Not yet implemented in TimsFramesMsRunReader %1.\n")
76 .arg(__LINE__));
77
79}
80
81
84{
85 return msp_timsData->getMassSpectrumCstSPtrByRawIndex(spectrum_index);
86}
87
88
91 bool want_binary_data) const
92{
93
94 QualifiedMassSpectrum mass_spectrum;
95
96 msp_timsData->getQualifiedMassSpectrumByRawIndex(
97 getMsRunId(), mass_spectrum, spectrum_index, want_binary_data);
98 return mass_spectrum;
99}
100
101
102void
105{
106 qDebug() << "Reading the spectrum collection with no specific configuration.";
107 MsRunReadConfig config;
108 readSpectrumCollection2(config, handler);
109}
110
111
112void
115{
116 qDebug().noquote() << "Reading the spectrum collection with this "
117 "specific configuration:"
118 << config.toString();
119
120 // We want to restrict the data reading process to the configuration provided
121 // as parameter.
122 // The frame record subset will include all frames fulfiling the config
123 // requirements (retention time and ms levels)
124
125 std::vector<std::size_t> subset_of_tims_frame_ids;
126
127 // The very first sorting of frames is perfomed on the basis of the
128 // retention time range that the user might have defined. Otherwise, all the
129 // frames are deemed suitable for next processing steps.
130
131 // Then, the MS level is checked against the config provided by the user.
132 for(auto const &frame_record : msp_timsData->getTimsFrameRecordList())
133 {
134 if(frame_record.frame_id == 0)
135 continue;
136
137 if(config.acceptRetentionTimeInSeconds(frame_record.frame_time))
138 {
139 std::size_t ms_level = 2;
140 if(frame_record.msms_type == 0)
141 ms_level = 1;
142 if(config.acceptMsLevel(ms_level))
143 {
144 subset_of_tims_frame_ids.push_back(frame_record.frame_id);
145 }
146 }
147 }
148
149 // At this point, we have a subset of frame records.
150 std::size_t frames_count = subset_of_tims_frame_ids.size();
151 qDebug() << "The number of retained RT range-matching frames:"
152 << frames_count;
153
154 // Inform the handler of the spectrum list so that it can handle feedback to
155 // the user.
156 handler.spectrumListHasSize(frames_count);
157
158 // Check for ion mobility range :
159 bool asked_ion_mobility_range = false;
160 quint32 mobility_scan_num_begin = 0;
161 quint32 mobility_scan_num_end = 0;
162
163 if(!config
164 .getParameterValue(
166 .isNull())
167 {
168 asked_ion_mobility_range = true;
169
170 mobility_scan_num_begin =
171 config
174 .toUInt();
175 mobility_scan_num_end =
176 config
179 .toUInt();
180 }
181
182 // Check for m/z resolution downgrading
183 std::size_t mz_index_merge_window = 0;
184 if(!config
185 .getParameterValue(
187 .isNull())
188 {
189 mz_index_merge_window =
190 config
193 .toUInt();
194 }
195
196 // Check for m/z range selection
197 bool asked_mz_range = false;
198 double mz_range_begin = -1;
199 double mz_range_end = -1;
200
201 if(!config
202 .getParameterValue(
204 .isNull())
205 {
206 asked_mz_range = true;
207
208 mz_range_begin =
209 config
212 .toDouble();
213 mz_range_end = config
216 .toDouble();
217
218 qDebug() << "The m/z range asked is: " << mz_range_begin << "--"
219 << mz_range_end;
220 }
221
222 // Now that we have gathered all the configuration bits from the user,
223 // iterate in the vector of frames records and for each frame ask that
224 // all its scans be combined into a single mass spectrum, by taking into
225 // account the configuration bits provided by the user.
226
227 const std::vector<FrameIdDescr> &frame_id_descr_list =
228 msp_timsData->getFrameIdDescrList();
229
230 // The scan index is the index of the scan in the *whole* mass data file, it
231 // is a sequential number of scans over all the frames.
232 std::size_t scan_index = 0; // iterate in each spectrum
233
234 for(std::size_t tims_frame_id : subset_of_tims_frame_ids)
235 {
236
237 if(handler.shouldStop())
238 {
239 qDebug() << "The operation was cancelled. Breaking the loop.";
241 QObject::tr("reading TimsTOF job cancelled by the user"));
242 }
243
244 const FrameIdDescr &current_frame_record =
245 frame_id_descr_list[tims_frame_id];
246
247 scan_index = current_frame_record.m_cumulSize;
248 TimsFrameCstSPtr tims_frame_csp =
249 msp_timsData->getTimsFrameCstSPtrCached(tims_frame_id);
250
251 // If the user wants to select specific 1/Ko values, we need to
252 // compute the ion mobility scan value starting from that 1/Ko value in
253 // each frame.
254
255 if(!config
256 .getParameterValue(
259 .isNull())
260 {
261 asked_ion_mobility_range = true;
262
263 mobility_scan_num_begin =
264 tims_frame_csp.get()->getScanNumFromOneOverK0(
265 config
266 .getParameterValue(
269 .toDouble());
270
271 mobility_scan_num_end = tims_frame_csp.get()->getScanNumFromOneOverK0(
272 config
273 .getParameterValue(
276 .toDouble());
277 }
278
279 // Get to know the size of the frame, that is,
280 // the number of mobility scans in it.
281 quint32 mobility_scans_count = tims_frame_csp->getTotalNumberOfScans();
282
283 // Provide these two variables for the function below to fill in the
284 // values.
285 quint32 min_mz_index_out = 0;
286 quint32 max_mz_index_out = 0;
287
288 Trace spectrum;
289
290 if(config.needPeakList())
291 {
292 if(asked_ion_mobility_range)
293 {
294 if(mobility_scan_num_end > (mobility_scans_count - 1))
295 {
296 mobility_scan_num_end = mobility_scans_count - 1;
297 }
298 }
299 else
300 {
301 mobility_scan_num_begin = 0;
302 mobility_scan_num_end = mobility_scans_count - 1;
303 }
304
305 if(asked_mz_range)
306 {
307 spectrum = tims_frame_csp->cumulateScansToTraceMzDownResolution2(
308 mz_index_merge_window,
309 mz_range_begin,
310 mz_range_end,
311 mobility_scan_num_begin,
312 mobility_scan_num_end,
313 min_mz_index_out,
314 max_mz_index_out);
315 }
316 else
317 {
318 spectrum = tims_frame_csp->cumulateScansToTraceMzDownResolution(
319 mz_index_merge_window,
320 mobility_scan_num_begin,
321 mobility_scan_num_end,
322 min_mz_index_out,
323 max_mz_index_out);
324 }
325 }
326
327 // qDebug() << "Got min_mz_index_out:" << min_mz_index_out;
328 // qDebug() << "Got max_mz_index_out:" << max_mz_index_out;
329
330 // Now, with or without the peak list, we have to craft a qualified mass
331 // spectrum that will hold all the data about the data in it.
332 QualifiedMassSpectrum mass_spectrum;
333
334 MassSpectrumId spectrum_id;
335
336 spectrum_id.setSpectrumIndex(tims_frame_id);
337 spectrum_id.setMsRunId(getMsRunId());
338
339 // Can be modified to add bits that might help our case
340 spectrum_id.setNativeId(
341 QString("frame id=%1 scan index=%2 im_begin=%3 im_end=%4")
342 .arg(tims_frame_id)
343 .arg(scan_index)
344 .arg(mobility_scan_num_begin)
345 .arg(mobility_scan_num_end));
346
347 mass_spectrum.setMassSpectrumId(spectrum_id);
348
349 // We want to document the retention time!
350 mass_spectrum.setRtInSeconds(tims_frame_csp.get()->getTime());
351
352 // We do want to document the ms level of the spectrum and possibly
353 // the precursor's m/z and charge.
354 unsigned int frame_ms_level = tims_frame_csp.get()->getMsLevel();
355 mass_spectrum.setMsLevel(frame_ms_level);
356
357 // We cannot set this value because the inverse mobility (1/Ko) measured
358 // by Bruker does not correspond formally to a drift time.
359 mass_spectrum.setDtInMilliSeconds(-1);
360
361 mass_spectrum.setParameterValue(
363 static_cast<qlonglong>(tims_frame_csp->getTotalNumberOfScans()));
364 mass_spectrum.setParameterValue(
366 tims_frame_csp->getOneOverK0Transformation(0));
367 mass_spectrum.setParameterValue(
369 tims_frame_csp->getOneOverK0Transformation(
370 tims_frame_csp->getTotalNumberOfScans() - 1));
371
372 // qDebug() << "Setting mass spectrum parameter value mz index begin to:"
373 // << min_mz_index_out;
374 mass_spectrum.setParameterValue(
376 min_mz_index_out);
377
378 // qDebug() << "Setting mass spectrum parameter value mz index end to:"
379 // << max_mz_index_out;
380 mass_spectrum.setParameterValue(
382
383 if(config.needPeakList())
384 {
385 mass_spectrum.setEmptyMassSpectrum(false);
386
387 mass_spectrum.setMassSpectrumSPtr(
388 std::make_shared<MassSpectrum>(spectrum));
389 }
390 else
391 {
392 mass_spectrum.setEmptyMassSpectrum(true);
393 }
394
395 handler.setQualifiedMassSpectrum(mass_spectrum);
396 }
397}
398
399
400void
402 [[maybe_unused]] SpectrumCollectionHandlerInterface &handler,
403 [[maybe_unused]] unsigned int ms_level)
404{
405 qDebug();
406}
407
408
409std::size_t
411{
412 return msp_timsData->getTotalNumberOfScans();
413}
414
415
416bool
418{
419 return false;
420}
421
422
423bool
425{
426 msp_timsData = nullptr;
427 return true;
428}
429
430bool
432{
433 if(msp_timsData == nullptr)
434 {
435 initialize();
436 }
437 return true;
438}
439
440
443 std::size_t spectrum_index [[maybe_unused]],
444 pappso::PrecisionPtr precision [[maybe_unused]]) const
445{
446 throw ExceptionNotImplemented(QObject::tr("Not implemented %1 %2 %3")
447 .arg(__FILE__)
448 .arg(__FUNCTION__)
449 .arg(__LINE__));
450}
451
454 const pappso::QualifiedMassSpectrum &mass_spectrum [[maybe_unused]],
455 pappso::PrecisionPtr precision [[maybe_unused]]) const
456{
457 throw ExceptionNotImplemented(QObject::tr("Not implemented %1 %2 %3")
458 .arg(__FILE__)
459 .arg(__FUNCTION__)
460 .arg(__LINE__));
461}
462
469
470
471Trace
473{
474 // Use the Sqlite database to fetch the total ion current chromatogram (TIC
475 // chromatogram).
476
478
479 return msp_timsData->getTicChromatogram();
480}
481
482
483Trace
485{
486
487 // We want to compute the TIC chromatogram, not load the chromatogram that
488 // is located in the SQL database.
489 //
490 // For this, we need to iterated into the frames and ask for MS1 spectra
491 // only. msp_timsData has that information:
492 //
493 // std::vector<FrameIdDescr> m_frameIdDescrList;
494 //
495 // and
496
497 // struct FrameIdDescr
498 // {
499 // std::size_t m_frameId; // frame id
500 // std::size_t m_size; // frame size (number of TOF scans in frame)
501 // std::size_t m_cumulSize; // cumulative size
502 // };
503
504 Trace tic_chromatogram;
505
506 const std::vector<FrameIdDescr> frame_descr_list =
507 msp_timsData->getFrameIdDescrList();
508
509 for(FrameIdDescr frame_id_descr : frame_descr_list)
510 {
511 TimsFrameCstSPtr tims_frame_csp =
512 msp_timsData->getTimsFrameCstSPtrCached(frame_id_descr.m_frameId);
513 std::size_t scan_begin = 0;
514 std::size_t scan_end = tims_frame_csp->getTotalNumberOfScans() - 1;
515
516 // By convention, a TIC chromatogram is only performed using MS1
517 // spectra.
518 if(tims_frame_csp->getMsLevel() == 1)
519 {
520
521 // Retention times are in seconds in the Bruker world.
522 double rt = tims_frame_csp->getTime();
523
524 tic_chromatogram.append(DataPoint(
525 rt,
526 tims_frame_csp->cumulateScansIntensities(scan_begin, scan_end)));
527 }
528 else
529 continue;
530 }
531
532 return tic_chromatogram;
533}
void setNativeId(const QString &native_id)
void setMsRunId(MsRunIdCstSPtr other)
void setSpectrumIndex(std::size_t index)
const QVariant getParameterValue(MsRunReadConfigParameter parameter) const
bool acceptMsLevel(std::size_t ms_level) const
bool acceptRetentionTimeInSeconds(double retention_time_in_seconds) const
base class to read MSrun the only way to build a MsRunReader object is to use the MsRunReaderFactory
Definition msrunreader.h:63
MsRunIdCstSPtr mcsp_msRunId
const MsRunIdCstSPtr & getMsRunId() const
Class representing a fully specified mass spectrum.
void setDtInMilliSeconds(pappso_double rt)
Set the drift time in milliseconds.
void setMassSpectrumId(const MassSpectrumId &iD)
Set the MassSpectrumId.
void setMsLevel(uint ms_level)
Set the mass spectrum level.
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
virtual void setQualifiedMassSpectrum(const QualifiedMassSpectrum &spectrum)=0
virtual std::size_t spectrumListSize() const override
get the totat number of spectrum conained in the MSrun data file
virtual void readSpectrumCollection(SpectrumCollectionHandlerInterface &handler) override
function to visit an MsRunReader and get each Spectrum in a spectrum collection handler
virtual MassSpectrumCstSPtr massSpectrumCstSPtr(std::size_t spectrum_index) override
virtual pappso::XicCoordSPtr newXicCoordSPtrFromQualifiedMassSpectrum(const pappso::QualifiedMassSpectrum &mass_spectrum, pappso::PrecisionPtr precision) const override
get a xic coordinate object from a given spectrum
virtual QualifiedMassSpectrum qualifiedMassSpectrum(std::size_t spectrum_index, bool want_binary_data=true) const override
get a QualifiedMassSpectrum class given its scan number
virtual MassSpectrumSPtr massSpectrumSPtr(std::size_t spectrum_index) override
get a MassSpectrumSPtr class given its spectrum index
virtual Trace getTicChromatogram() override
get a TIC chromatogram
virtual bool releaseDevice() override
release data back end device if a the data back end is released, the developper has to use acquireDev...
virtual void readSpectrumCollection2(const MsRunReadConfig &config, SpectrumCollectionHandlerInterface &handler) override
virtual TimsDataSp getTimsDataSPtr()
give an access to the underlying raw data pointer
virtual bool hasScanNumbers() const override
tells if spectra can be accessed using scan numbers by default, it returns false. Only overrided func...
virtual bool acquireDevice() override
acquire data back end device
virtual bool accept(const QString &file_name) const override
tells if the reader is able to handle this file must be implemented by private MS run reader,...
TimsFramesMsRunReader(MsRunIdCstSPtr &msrun_id_csp)
virtual pappso::XicCoordSPtr newXicCoordSPtrFromSpectrumIndex(std::size_t spectrum_index, pappso::PrecisionPtr precision) const override
get a xic coordinate object from a given spectrum index
virtual void readSpectrumCollectionByMsLevel(SpectrumCollectionHandlerInterface &handler, unsigned int ms_level) override
function to visit an MsRunReader and get each Spectrum in a spectrum collection handler by Ms Levels
A simple container of DataPoint instances.
Definition trace.h:148
size_t append(const DataPoint &data_point)
appends a datapoint and return new size
Definition trace.cpp:648
tries to keep as much as possible monoisotopes, removing any possible C13 peaks and changes multichar...
Definition aa.cpp:39
std::shared_ptr< const MsRunId > MsRunIdCstSPtr
Definition msrunid.h:46
std::shared_ptr< TimsData > TimsDataSp
shared pointer on a TimsData object
Definition timsdata.h:50
std::shared_ptr< const MassSpectrum > MassSpectrumCstSPtr
@ TimsFrameInvKoBegin
Bruker's Tims tof 1/K0 frame start range.
@ TimsFrameMzIndexBegin
Bruker's Tims tof mz index frame start range.
@ TimsFrameMzIndexEnd
Bruker's Tims tof mz index frame end range.
@ TimsFrameInvKoEnd
Bruker's Tims tof 1/K0 frame end range.
@ rt
Retention time.
std::shared_ptr< MassSpectrum > MassSpectrumSPtr
std::shared_ptr< const TimsFrame > TimsFrameCstSPtr
Definition timsframe.h:42
MsRunReadConfigParameter
std::shared_ptr< XicCoord > XicCoordSPtr
Definition xiccoord.h:43
std::size_t m_cumulSize
Definition timsdata.h:57