OpenMS  3.0.0
MSExperiment.h
Go to the documentation of this file.
1 // --------------------------------------------------------------------------
2 // OpenMS -- Open-Source Mass Spectrometry
3 // --------------------------------------------------------------------------
4 // Copyright The OpenMS Team -- Eberhard Karls University Tuebingen,
5 // ETH Zurich, and Freie Universitaet Berlin 2002-2022.
6 //
7 // This software is released under a three-clause BSD license:
8 // * Redistributions of source code must retain the above copyright
9 // notice, this list of conditions and the following disclaimer.
10 // * Redistributions in binary form must reproduce the above copyright
11 // notice, this list of conditions and the following disclaimer in the
12 // documentation and/or other materials provided with the distribution.
13 // * Neither the name of any author or any participating institution
14 // may be used to endorse or promote products derived from this software
15 // without specific prior written permission.
16 // For a full list of authors, refer to the file AUTHORS.
17 // --------------------------------------------------------------------------
18 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
19 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
20 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
21 // ARE DISCLAIMED. IN NO EVENT SHALL ANY OF THE AUTHORS OR THE CONTRIBUTING
22 // INSTITUTIONS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
23 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
24 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
25 // OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
26 // WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
27 // OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
28 // ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 //
30 // --------------------------------------------------------------------------
31 // $Maintainer: Timo Sachsenberg $
32 // $Authors: Marc Sturm, Tom Waschischeck $
33 // --------------------------------------------------------------------------
34 
35 #pragma once
36 
43 
44 #include <vector>
45 
46 
47 namespace OpenMS
48 {
49  class Peak1D;
50  class ChromatogramPeak;
51 
70  class OPENMS_DLLAPI MSExperiment final :
71  public RangeManagerContainer<RangeRT, RangeMZ, RangeIntensity>,
73  {
74 
75 public:
76  typedef Peak1D PeakT;
78 
80 
81  typedef PeakT PeakType;
98  typedef std::vector<SpectrumType> Base;
100 
102 
103  typedef std::vector<SpectrumType>::iterator Iterator;
106  typedef std::vector<SpectrumType>::const_iterator ConstIterator;
112 
114  // Attention: these refer to the spectra vector only!
116  typedef Base::value_type value_type;
117  typedef Base::iterator iterator;
118  typedef Base::const_iterator const_iterator;
119 
120  inline Size size() const
121  {
122  return spectra_.size();
123  }
124 
125  inline void resize(Size s)
126  {
127  spectra_.resize(s);
128  }
129 
130  inline bool empty() const
131  {
132  return spectra_.empty();
133  }
134 
135  inline void reserve(Size s)
136  {
137  spectra_.reserve(s);
138  }
139 
140  inline SpectrumType& operator[] (Size n)
141  {
142  return spectra_[n];
143  }
144 
145  inline const SpectrumType& operator[] (Size n) const
146  {
147  return spectra_[n];
148  }
149 
150  inline Iterator begin()
151  {
152  return spectra_.begin();
153  }
154 
155  inline ConstIterator begin() const
156  {
157  return spectra_.begin();
158  }
159 
160  inline Iterator end()
161  {
162  return spectra_.end();
163  }
164 
165  inline ConstIterator end() const
166  {
167  return spectra_.end();
168  }
170 
171  // Aliases / chromatograms
172  void reserveSpaceSpectra(Size s);
173  void reserveSpaceChromatograms(Size s);
174 
176  MSExperiment();
177 
179  MSExperiment(const MSExperiment & source);
180 
182  MSExperiment(MSExperiment&&) = default;
183 
185  MSExperiment & operator=(const MSExperiment & source);
186 
188  MSExperiment& operator=(MSExperiment&&) & = default;
189 
191  MSExperiment & operator=(const ExperimentalSettings & source);
192 
194  bool operator==(const MSExperiment & rhs) const;
195 
197  bool operator!=(const MSExperiment & rhs) const;
198 
200 
201 
207  template <class Container>
208  void get2DData(Container& cont) const
209  {
210  for (typename Base::const_iterator spec = spectra_.begin(); spec != spectra_.end(); ++spec)
211  {
212  if (spec->getMSLevel() != 1)
213  {
214  continue;
215  }
216  typename Container::value_type s; // explicit object here, since instantiation within push_back() fails on VS<12
217  for (typename SpectrumType::const_iterator it = spec->begin(); it != spec->end(); ++it)
218  {
219  cont.push_back(s);
220  cont.back().setRT(spec->getRT());
221  cont.back().setMZ(it->getMZ());
222  cont.back().setIntensity(it->getIntensity());
223  }
224  }
225  }
226 
238  template <class Container>
239  void set2DData(const Container& container)
240  {
241  set2DData<false, Container>(container);
242  }
243 
258  template <class Container>
259  void set2DData(const Container& container, const StringList& store_metadata_names)
260  {
261  // clean up the container first
262  clear(true);
263  SpectrumType* spectrum = nullptr;
264  typename PeakType::CoordinateType current_rt = -std::numeric_limits<typename PeakType::CoordinateType>::max();
265  for (typename Container::const_iterator iter = container.begin(); iter != container.end(); ++iter)
266  {
267  // check if the retention time has changed
268  if (current_rt != iter->getRT() || spectrum == nullptr)
269  {
270  // append new spectrum
271  if (current_rt > iter->getRT())
272  {
273  throw Exception::Precondition(__FILE__, __LINE__, OPENMS_PRETTY_FUNCTION, "Input container is not sorted!");
274  }
275  current_rt = iter->getRT();
276  spectrum = createSpec_(current_rt, store_metadata_names);
277  }
278 
279  // add either data point or mass traces (depending on template argument value)
280  ContainerAdd_<typename Container::value_type, false>::addData_(spectrum, &(*iter), store_metadata_names);
281  }
282  }
283 
301  template <bool add_mass_traces, class Container>
302  void set2DData(const Container& container)
303  {
304  // clean up the container first
305  clear(true);
306  SpectrumType* spectrum = nullptr;
307  typename PeakType::CoordinateType current_rt = -std::numeric_limits<typename PeakType::CoordinateType>::max();
308  for (typename Container::const_iterator iter = container.begin(); iter != container.end(); ++iter)
309  {
310  // check if the retention time has changed
311  if (current_rt != iter->getRT() || spectrum == nullptr)
312  {
313  // append new spectrum
314  if (current_rt > iter->getRT())
315  {
316  throw Exception::Precondition(__FILE__, __LINE__, OPENMS_PRETTY_FUNCTION, "Input container is not sorted!");
317  }
318  current_rt = iter->getRT();
319  spectrum = createSpec_(current_rt);
320  }
321 
322  // add either data point or mass traces (depending on template argument value)
324  }
325  }
326 
328 
329 
331 
332  AreaIterator areaBegin(CoordinateType min_rt, CoordinateType max_rt, CoordinateType min_mz, CoordinateType max_mz);
334 
336  AreaIterator areaEnd();
337 
339  ConstAreaIterator areaBeginConst(CoordinateType min_rt, CoordinateType max_rt, CoordinateType min_mz, CoordinateType max_mz) const;
340 
342  ConstAreaIterator areaEndConst() const;
343 
344  // for fast pyOpenMS access to MS1 peak data in format: [rt, [mz, intensity]]
346  std::vector<float>& rt,
347  std::vector<std::vector<float>>& mz,
348  std::vector<std::vector<float>>& intensity) const
349  {
350  float t = -1.0;
351  for (auto it = areaBeginConst(min_rt, max_rt, min_mz, max_mz); it != areaEndConst(); ++it)
352  {
353  if (it.getRT() != t)
354  {
355  t = it.getRT();
356  rt.emplace_back(t);
357  mz.resize(mz.size() + 1);
358  rt.resize(rt.size() + 1);
359  intensity.resize(intensity.size() + 1);
360  }
361  mz.back().emplace_back(it->getMZ());
362  intensity.back().emplace_back(it->getIntensity());
363  }
364  }
365 
366  // for fast pyOpenMS access to MS1 peak data in format: [rt, mz, intensity]
368  std::vector<float>& rt,
369  std::vector<float>& mz,
370  std::vector<float>& intensity) const
371  {
372  for (auto it = areaBeginConst(min_rt, max_rt, min_mz, max_mz); it != areaEndConst(); ++it)
373  {
374  rt.emplace_back(it.getRT());
375  mz.emplace_back(it->getMZ());
376  intensity.emplace_back(it->getIntensity());
377  }
378  }
379 
380 
388  ConstIterator RTBegin(CoordinateType rt) const;
389 
397  ConstIterator RTEnd(CoordinateType rt) const;
398 
404  Iterator RTBegin(CoordinateType rt);
405 
411  Iterator RTEnd(CoordinateType rt);
412 
414 
420  // Docu in base class
422  void updateRanges() override;
423 
429  void updateRanges(Int ms_level);
430 
432  CoordinateType getMinMZ() const;
433 
435  CoordinateType getMaxMZ() const;
436 
438  CoordinateType getMinRT() const;
439 
441  CoordinateType getMaxRT() const;
442 
444  UInt64 getSize() const;
445 
447  const std::vector<UInt>& getMSLevels() const;
448 
450 
453  UInt64 getSqlRunID() const;
454 
456  void setSqlRunID(UInt64 id);
457 
460 
465  void sortSpectra(bool sort_mz = true);
466 
472  void sortChromatograms(bool sort_rt = true);
473 
479  bool isSorted(bool check_mz = true) const;
480 
482 
484  void reset();
485 
491  bool clearMetaDataArrays();
492 
494  const ExperimentalSettings& getExperimentalSettings() const;
495 
497  ExperimentalSettings& getExperimentalSettings();
498 
500  void getPrimaryMSRunPath(StringList& toFill) const;
501 
507  ConstIterator getPrecursorSpectrum(ConstIterator iterator) const;
508 
514  int getPrecursorSpectrum(int zero_based_index) const;
515 
517  void swap(MSExperiment& from);
518 
520  void setSpectra(const std::vector<MSSpectrum>& spectra);
521  void setSpectra(std::vector<MSSpectrum>&& spectra);
522 
524  void addSpectrum(const MSSpectrum& spectrum);
525  void addSpectrum(MSSpectrum&& spectrum);
526 
528  const std::vector<MSSpectrum>& getSpectra() const;
529 
531  std::vector<MSSpectrum>& getSpectra();
532 
534  void setChromatograms(const std::vector<MSChromatogram>& chromatograms);
535  void setChromatograms(std::vector<MSChromatogram>&& chromatograms);
536 
538  void addChromatogram(const MSChromatogram& chromatogram);
539  void addChromatogram(MSChromatogram&& chrom);
540 
542  const std::vector<MSChromatogram>& getChromatograms() const;
543 
545  std::vector<MSChromatogram>& getChromatograms();
546 
548 
549  MSChromatogram& getChromatogram(Size id);
551 
553  MSSpectrum& getSpectrum(Size id);
554 
556  Size getNrSpectra() const;
557 
559  Size getNrChromatograms() const;
561 
572  const MSChromatogram calculateTIC(float rt_bin_size = 0, UInt ms_level = 1) const;
573 
579  void clear(bool clear_meta_data);
580 
582  bool containsScanOfLevel(size_t ms_level) const;
583 
585  bool hasZeroIntensities(size_t ms_level) const;
586 
588  bool hasPeptideIdentifications() const;
589 
590  protected:
592  std::vector<UInt> ms_levels_;
596  std::vector<MSChromatogram > chromatograms_;
598  std::vector<SpectrumType> spectra_;
599 
600 private:
601 
603  template<typename ContainerValueType, bool addMassTraces>
605  {
606  static void addData_(SpectrumType* spectrum, const ContainerValueType* item);
607  static void addData_(SpectrumType* spectrum, const ContainerValueType* item, const StringList& store_metadata_names);
608  };
609 
610  template<typename ContainerValueType>
611  struct ContainerAdd_<ContainerValueType, false>
612  {
614  static void addData_(SpectrumType* spectrum, const ContainerValueType* item)
615  {
616  // create temporary peak and insert it into spectrum
617  spectrum->insert(spectrum->end(), PeakType());
618  spectrum->back().setIntensity(item->getIntensity());
619  spectrum->back().setPosition(item->getMZ());
620  }
622  static void addData_(SpectrumType* spectrum, const ContainerValueType* item, const StringList& store_metadata_names)
623  {
624  addData_(spectrum, item);
625  for (StringList::const_iterator itm = store_metadata_names.begin(); itm != store_metadata_names.end(); ++itm)
626  {
627  float val = std::numeric_limits<float>::quiet_NaN();
628  if (item->metaValueExists(*itm)) val = item->getMetaValue(*itm);
629  spectrum->getFloatDataArrays()[itm - store_metadata_names.begin()].push_back(val);
630  }
631  }
632  };
633 
634  template<typename ContainerValueType>
635  struct ContainerAdd_<ContainerValueType, true>
636  {
638  static void addData_(SpectrumType* spectrum, const ContainerValueType* item)
639  {
640  if (item->metaValueExists("num_of_masstraces"))
641  {
642  Size mts = item->getMetaValue("num_of_masstraces");
643  int charge = (item->getCharge()==0 ? 1 : item->getCharge()); // set to 1 if charge is 0, otherwise div/0 below
644  for (Size i = 0; i < mts; ++i)
645  {
646  String meta_name = String("masstrace_intensity_") + i;
647  if (!item->metaValueExists(meta_name))
648  {
649  throw Exception::Precondition(__FILE__, __LINE__, OPENMS_PRETTY_FUNCTION, String("Meta value '") + meta_name + "' expected but not found in container.");
650  }
651  ContainerValueType p;
652  p.setIntensity(item->getMetaValue(meta_name));
653  p.setPosition(item->getMZ() + Constants::C13C12_MASSDIFF_U / charge * i);
655  }
656  }
658  }
659  };
660 
661  /*
662  @brief Append a spectrum to current MSExperiment
663 
664  @param rt RT of new spectrum
665  @return Pointer to newly created spectrum
666  */
667  SpectrumType* createSpec_(PeakType::CoordinateType rt);
668 
669  /*
670  @brief Append a spectrum including floatdata arrays to current MSExperiment
671 
672  @param rt RT of new spectrum
673  @param metadata_names Names of floatdata arrays attached to this spectrum
674  @return Pointer to newly created spectrum
675  */
676  SpectrumType* createSpec_(PeakType::CoordinateType rt, const StringList& metadata_names);
677 
678  };
679 
681  OPENMS_DLLAPI std::ostream& operator<<(std::ostream& os, const MSExperiment& exp);
682 
683 } // namespace OpenMS
684 
686 
687 
OpenMS::PeakType
Peak2D PeakType
Definition: MassTrace.h:47
OpenMS::MSExperiment::set2DData
void set2DData(const Container &container)
Assignment of a data container with RT and MZ to an MSExperiment.
Definition: MSExperiment.h:239
OpenMS::UInt64
OPENMS_UINT64_TYPE UInt64
Unsigned integer type (64bit)
Definition: Types.h:77
AreaIterator.h
double
OpenMS::Internal::AreaIterator
Forward iterator for an area of peaks in an experiment.
Definition: AreaIterator.h:57
OpenMS::MSExperiment::Base
std::vector< SpectrumType > Base
STL base class type.
Definition: MSExperiment.h:98
OpenMS::MSExperiment::set2DData
void set2DData(const Container &container)
Assignment of a data container with RT and MZ to an MSExperiment.
Definition: MSExperiment.h:302
OpenMS::MSExperiment::value_type
Base::value_type value_type
Definition: MSExperiment.h:116
OpenMS::ExperimentalSettings
Description of the experimental settings.
Definition: ExperimentalSettings.h:59
OpenMS::String
A more convenient string class.
Definition: String.h:58
OpenMS::MSExperiment::begin
Iterator begin()
Definition: MSExperiment.h:150
OpenMS::Exception::Precondition
Precondition failed exception.
Definition: Exception.h:157
KDTree::operator!=
bool operator!=(_Iterator< _Val, _Ref, _Ptr > const &, _Iterator< _Val, _Ref, _Ptr > const &)
Definition: KDTree.h:824
OpenMS::MSExperiment
In-Memory representation of a mass spectrometry run.
Definition: MSExperiment.h:70
OpenMS::Size
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:127
OpenMS::MSExperiment::ChromatogramPeakType
ChromatogramPeakT ChromatogramPeakType
Chromatogram peak type.
Definition: MSExperiment.h:84
OpenMS::MSExperiment::IntensityType
PeakType::IntensityType IntensityType
Intensity type of peaks.
Definition: MSExperiment.h:88
OpenMS::MSExperiment::set2DData
void set2DData(const Container &container, const StringList &store_metadata_names)
Assignment of a data container with RT and MZ to an MSExperiment.
Definition: MSExperiment.h:259
OpenMS::MSExperiment::CoordinateType
PeakType::CoordinateType CoordinateType
Coordinate type of peak positions.
Definition: MSExperiment.h:86
OpenMS::MSExperiment::ContainerAdd_< ContainerValueType, false >::addData_
static void addData_(SpectrumType *spectrum, const ContainerValueType *item, const StringList &store_metadata_names)
general method for adding data points, including metadata arrays (populated from metainfointerface)
Definition: MSExperiment.h:622
OpenMS::MSExperiment::const_iterator
Base::const_iterator const_iterator
Definition: MSExperiment.h:118
OpenMS::MSExperiment::size
Size size() const
Definition: MSExperiment.h:120
OpenMS::MSExperiment::SpectrumType
MSSpectrum SpectrumType
Spectrum Type.
Definition: MSExperiment.h:94
OpenMS::MSExperiment::reserve
void reserve(Size s)
Definition: MSExperiment.h:135
OpenMS::MSExperiment::ChromatogramType
MSChromatogram ChromatogramType
Chromatogram type.
Definition: MSExperiment.h:96
OpenMS::Int
int Int
Signed integer type.
Definition: Types.h:102
OpenMS
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:47
OpenMS::MSExperiment::Iterator
std::vector< SpectrumType >::iterator Iterator
Mutable iterator.
Definition: MSExperiment.h:104
OpenMS::RangeManager
Handles the management of a multidimensional range, e.g. RangeMZ and RangeIntensity for spectra.
Definition: RangeManager.h:588
OpenMS::MSExperiment::RangeManagerContainerType
RangeManagerContainer< RangeRT, RangeMZ, RangeIntensity > RangeManagerContainerType
RangeManager type.
Definition: MSExperiment.h:92
OpenMS::MSExperiment::RangeManagerType
RangeManager< RangeRT, RangeMZ, RangeIntensity > RangeManagerType
RangeManager type.
Definition: MSExperiment.h:90
Exception.h
OpenMS::MSExperiment::ContainerAdd_
Helper class to add either general data points in set2DData or use mass traces from meta values.
Definition: MSExperiment.h:604
OpenMS::MSExperiment::iterator
Base::iterator iterator
Definition: MSExperiment.h:117
OpenMS::MSExperiment::end
ConstIterator end() const
Definition: MSExperiment.h:165
OpenMS::ChromatogramPeak
A 1-dimensional raw data point or peak for chromatograms.
Definition: ChromatogramPeak.h:53
OpenMS::MSExperiment::ChromatogramPeakT
ChromatogramPeak ChromatogramPeakT
Definition: MSExperiment.h:77
OpenMS::MSExperiment::begin
ConstIterator begin() const
Definition: MSExperiment.h:155
OpenMS::MSExperiment::empty
bool empty() const
Definition: MSExperiment.h:130
OpenMS::Peak1D
A 1-dimensional raw data point or peak.
Definition: Peak1D.h:53
OpenMS::operator<<
std::ostream & operator<<(std::ostream &os, const AccurateMassSearchResult &amsr)
OpenMS::StringList
std::vector< String > StringList
Vector of String.
Definition: ListUtils.h:70
OpenMS::UInt
unsigned int UInt
Unsigned integer type.
Definition: Types.h:94
OpenMS::MSExperiment::resize
void resize(Size s)
Definition: MSExperiment.h:125
OpenMS::RangeManagerContainer
Definition: RangeManager.h:919
OpenMS::MSSpectrum::getFloatDataArrays
const FloatDataArrays & getFloatDataArrays() const
Returns a const reference to the float meta data arrays.
Iterator
ExperimentalSettings.h
OpenMS::Internal::operator==
bool operator==(const IDBoostGraph::ProteinGroup &lhs, const IDBoostGraph::ProteinGroup &rhs)
OpenMS::MSExperiment::end
Iterator end()
Definition: MSExperiment.h:160
MSChromatogram.h
float
OpenMS::MSExperiment::get2DPeakData
void get2DPeakData(CoordinateType min_rt, CoordinateType max_rt, CoordinateType min_mz, CoordinateType max_mz, std::vector< float > &rt, std::vector< std::vector< float >> &mz, std::vector< std::vector< float >> &intensity) const
Definition: MSExperiment.h:345
OpenMS::MSExperiment::PeakT
Peak1D PeakT
Definition: MSExperiment.h:76
OpenMS::MSChromatogram
The representation of a chromatogram.
Definition: MSChromatogram.h:53
OpenMS::MSExperiment::get2DPeakData
void get2DPeakData(CoordinateType min_rt, CoordinateType max_rt, CoordinateType min_mz, CoordinateType max_mz, std::vector< float > &rt, std::vector< float > &mz, std::vector< float > &intensity) const
Definition: MSExperiment.h:367
StandardDeclarations.h
OpenMS::MSExperiment::ContainerAdd_< ContainerValueType, true >::addData_
static void addData_(SpectrumType *spectrum, const ContainerValueType *item)
specialization for adding feature mass traces (does not support metadata_names currently)
Definition: MSExperiment.h:638
OpenMS::MSExperiment::ms_levels_
std::vector< UInt > ms_levels_
MS levels of the data.
Definition: MSExperiment.h:592
OpenMS::MSExperiment::chromatograms_
std::vector< MSChromatogram > chromatograms_
chromatograms
Definition: MSExperiment.h:596
OpenMS::MSExperiment::get2DData
void get2DData(Container &cont) const
Reads out a 2D Spectrum.
Definition: MSExperiment.h:208
OpenMS::MSExperiment::ConstIterator
std::vector< SpectrumType >::const_iterator ConstIterator
Non-mutable iterator.
Definition: MSExperiment.h:106
OpenMS::MSExperiment::AreaIterator
Internal::AreaIterator< PeakT, PeakT &, PeakT *, Iterator, SpectrumType::Iterator > AreaIterator
Mutable area iterator type (for traversal of a rectangular subset of the peaks)
Definition: MSExperiment.h:108
OpenMS::MSExperiment::ConstAreaIterator
Internal::AreaIterator< const PeakT, const PeakT &, const PeakT *, ConstIterator, SpectrumType::ConstIterator > ConstAreaIterator
Immutable area iterator type (for traversal of a rectangular subset of the peaks)
Definition: MSExperiment.h:110
OpenMS::Constants::C13C12_MASSDIFF_U
const double C13C12_MASSDIFF_U
Definition: Constants.h:121
OpenMS::MSExperiment::ContainerAdd_< ContainerValueType, false >::addData_
static void addData_(SpectrumType *spectrum, const ContainerValueType *item)
general method for adding data points
Definition: MSExperiment.h:614
OpenMS::MSSpectrum
The representation of a 1D spectrum.
Definition: MSSpectrum.h:66
StandardTypes.h
OpenMS::MSExperiment::total_size_
UInt64 total_size_
Number of all data points.
Definition: MSExperiment.h:594
MSSpectrum.h
OpenMS::MSExperiment::spectra_
std::vector< SpectrumType > spectra_
spectra
Definition: MSExperiment.h:598