Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages
FeatureFinderAlgorithmPickedHelperStructs.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-2013.
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: Stephan Aiche $
32 // $Authors: Marc Sturm, Stephan Aiche $
33 // --------------------------------------------------------------------------
34 
35 #ifndef OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMPICKEDHELPERSTRUCTS_H
36 #define OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMPICKEDHELPERSTRUCTS_H
37 
38 #include <OpenMS/CONCEPT/Types.h>
40 
44 
45 #include <vector>
46 
47 namespace OpenMS
48 {
49 
57  {
58 
62  struct OPENMS_DLLAPI Seed
63  {
70 
72  bool operator<(const Seed & rhs) const
73  {
74  return intensity < rhs.intensity;
75  }
76 
77  };
78 
82  template <class PeakType>
83  struct MassTrace
84  {
86  const PeakType * max_peak;
89 
92 
94  std::vector<std::pair<DoubleReal, const PeakType *> > peaks;
95 
98  {
99  ConvexHull2D::PointArrayType hull_points(peaks.size());
100  for (Size i = 0; i < peaks.size(); ++i)
101  {
102  hull_points[i][0] = peaks[i].first;
103  hull_points[i][1] = peaks[i].second->getMZ();
104  }
105  ConvexHull2D hull;
106  hull.addPoints(hull_points);
107  return hull;
108  }
109 
112  {
113  if (peaks.empty()) return;
114 
115  max_rt = peaks.begin()->first;
116  max_peak = peaks.begin()->second;
117 
118  for (Size i = 1; i < peaks.size(); ++i)
119  {
120  if (peaks[i].second->getIntensity() > max_peak->getIntensity())
121  {
122  max_rt = peaks[i].first;
123  max_peak = peaks[i].second;
124  }
125  }
126  }
127 
130  {
131  DoubleReal sum = 0.0;
132  DoubleReal intensities = 0.0;
133  for (Size i = 0; i < peaks.size(); ++i)
134  {
135  sum += peaks[i].second->getMZ() * peaks[i].second->getIntensity();
136  intensities += peaks[i].second->getIntensity();
137  }
138  return sum / intensities;
139  }
140 
142  bool isValid() const
143  {
144  return peaks.size() >= 3;
145  }
146 
147  };
148 
152  template <class PeakType>
153  struct MassTraces :
154  public std::vector<MassTrace<PeakType> >
155  {
158  max_trace(0)
159  {
160  }
161 
164  {
165  Size sum = 0;
166  for (Size i = 0; i < this->size(); ++i)
167  {
168  sum += this->at(i).peaks.size();
169  }
170  return sum;
171  }
172 
174  bool isValid(DoubleReal seed_mz, DoubleReal trace_tolerance)
175  {
176  //Abort if too few traces were found
177  if (this->size() < 2) return false;
178 
179  //Abort if the seed was removed
180  for (Size j = 0; j < this->size(); ++j)
181  {
182  if (std::fabs(seed_mz - this->at(j).getAvgMZ()) <= trace_tolerance)
183  {
184  return true;
185  }
186  }
187  return false;
188  }
189 
196  {
197  if (!this->size())
198  {
199  throw Exception::Precondition(__FILE__, __LINE__, __PRETTY_FUNCTION__, "There must be at least one trace to determine the theoretical maximum trace!");
200  }
201 
202  Size max = 0;
203  DoubleReal max_int = this->at(0).theoretical_int;
204  for (Size i = 1; i < this->size(); ++i)
205  {
206  if (this->at(i).theoretical_int > max_int)
207  {
208  max_int = this->at(i).theoretical_int;
209  max = i;
210  }
211  }
212  return max;
213  }
214 
217  {
218  if (this->size() == 0)
219  {
220  baseline = 0.0;
221  return;
222  }
223  bool first = true;
224  for (Size i = 0; i < this->size(); ++i)
225  {
226  for (Size j = 0; j < this->at(i).peaks.size(); ++j)
227  {
228  if (first)
229  {
230  baseline = this->at(i).peaks[j].second->getIntensity();
231  first = false;
232  }
233  if (this->at(i).peaks[j].second->getIntensity() < baseline)
234  {
235  baseline = this->at(i).peaks[j].second->getIntensity();
236  }
237  }
238  }
239  }
240 
246  std::pair<DoubleReal, DoubleReal> getRTBounds() const
247  {
248  if (!this->size())
249  {
250  throw Exception::Precondition(__FILE__, __LINE__, __PRETTY_FUNCTION__, "There must be at least one trace to determine the RT boundaries!");
251  }
252 
253  DoubleReal min = std::numeric_limits<DoubleReal>::max();
254  DoubleReal max = -std::numeric_limits<DoubleReal>::max();
255  //Abort if the seed was removed
256  for (Size i = 0; i < this->size(); ++i)
257  {
258  for (Size j = 0; j < this->at(i).peaks.size(); ++j)
259  {
260  DoubleReal rt = this->at(i).peaks[j].first;
261  if (rt > max) max = rt;
262  if (rt < min) min = rt;
263  }
264  }
265  return std::make_pair(min, max);
266  }
267 
272  };
273 
277  struct OPENMS_DLLAPI TheoreticalIsotopePattern
278  {
280  std::vector<DoubleReal> intensity;
290  Size size() const
291  {
292  return intensity.size();
293  }
294 
295  };
296 
300  struct OPENMS_DLLAPI IsotopePattern
301  {
303  std::vector<SignedSize> peak;
305  std::vector<Size> spectrum;
307  std::vector<DoubleReal> intensity;
309  std::vector<DoubleReal> mz_score;
311  std::vector<DoubleReal> theoretical_mz;
314 
317  peak(size, -1),
318  spectrum(size),
319  intensity(size),
320  mz_score(size),
321  theoretical_mz(size)
322  {
323  }
324 
325  };
326 
327  };
328 }
329 
330 #endif // #ifndef OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMPICKEDHELPERSTRUCTS_H
MassTraces()
Constructor.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:157
float Real
Real type.
Definition: Types.h:109
Size spectrum
Spectrum index.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:65
void updateMaximum()
Sets the maximum to the highest contained peak of the trace.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:111
A 2-dimensional raw data point or peak.
Definition: Peak2D.h:55
static DoubleReal sum(IteratorType begin, IteratorType end)
Calculates the sum of a range of values.
Definition: StatisticFunctions.h:56
Size max_trace
Maximum intensity trace.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:269
Helper struct for a collection of mass traces used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:153
bool isValid(DoubleReal seed_mz, DoubleReal trace_tolerance)
Checks if still valid (seed still contained and enough traces)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:174
std::vector< PointType > PointArrayType
Definition: ConvexHull2D.h:79
Wrapper struct for all the classes needed by the FeatureFinderAlgorithmPicked and the associated clas...
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:56
Size optional_begin
Number of optional peaks at the beginning of the pattern.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:282
Helper structure for a theoretical isotope pattern used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:277
std::vector< DoubleReal > mz_score
m/z score of peak (0 if peak index is -1 or -2)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:309
Size size() const
Returns the size.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:290
std::vector< DoubleReal > intensity
Peak intensity (0 if peak index is -1 or -2)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:307
Size getPeakCount() const
Returns the peak count of all traces.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:163
Precondition failed exception.
Definition: Exception.h:167
IsotopePattern(Size size)
Constructor that resizes the internal vectors.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:316
std::vector< DoubleReal > theoretical_mz
Theoretical m/z value of the isotope peak.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:311
A 2-dimensional hull representation in [counter]clockwise direction - depending on axis labelling...
Definition: ConvexHull2D.h:75
TheoreticalIsotopePattern theoretical_pattern
Theoretical isotope pattern.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:313
const PeakType * max_peak
Maximum peak pointer.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:86
std::vector< Size > spectrum
Spectrum index (undefined if peak index is -1 or -2)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:305
bool isValid() const
Checks if this Trace is valid (has more than 2 points)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:142
DoubleReal max
The maximum intensity contribution before scaling the pattern to 1.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:286
DoubleReal theoretical_int
Theoretical intensity value (scaled to [0,1])
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:91
ConvexHull2D getConvexhull() const
determines the convex hull of the trace
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:97
bool operator<(const Seed &rhs) const
Comparison operator.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:72
DoubleReal baseline
Estimated baseline in the region of the feature (used for the fit)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:271
Helper structure for a found isotope pattern used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:300
DoubleReal max_rt
RT of maximum peak.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:88
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:144
std::vector< std::pair< DoubleReal, const PeakType * > > peaks
Contained peaks (pair of RT and pointer to peak)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:94
Size optional_end
Number of optional peaks at the end of the pattern.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:284
DoubleReal getAvgMZ() const
Returns the average m/z of all peaks in this trace (weighted by intensity)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:129
Size trimmed_left
The number of isotopes trimmed on the left side. This is needed to reconstruct the monoisotopic peak...
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:288
void addPoints(const PointArrayType &points)
Size getTheoreticalmaxPosition() const
Returns the theoretical maximum trace index.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:195
void updateBaseline()
Sets the baseline to the lowest contained peak of the trace.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:216
std::vector< DoubleReal > intensity
Vector of intensity contributions.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:280
std::pair< DoubleReal, DoubleReal > getRTBounds() const
Returns the RT boundaries of the mass traces.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:246
Size peak
Peak index.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:67
std::vector< SignedSize > peak
Peak index (-1 if peak was not found, -2 if it was removed to improve the isotope fit) ...
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:303
Real intensity
Intensity.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:69
Helper struct for mass traces used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:83
Helper structure for seeds used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:62

OpenMS / TOPP release 1.11.1 Documentation generated on Thu Nov 14 2013 11:19:13 using doxygen 1.8.5