OpenMS  2.4.0
WindowMower.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-2018.
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: Mathias Walzer $
32 // $Authors: Mathias Walzer, Timo Sachsenberg$
33 // --------------------------------------------------------------------------
34 //
35 #pragma once
36 
41 
42 #include <set>
43 
44 namespace OpenMS
45 {
46 
54  class OPENMS_DLLAPI WindowMower :
55  public DefaultParamHandler
56  {
57 public:
58 
59  // @name Constructors, destructors and assignment operators
60  // @{
62  WindowMower();
64  ~WindowMower() override;
65 
67  WindowMower(const WindowMower& source);
69  WindowMower& operator=(const WindowMower& source);
70  // @}
71 
73  template <typename SpectrumType>
75  {
76  typedef typename SpectrumType::ConstIterator ConstIterator;
77 
78  windowsize_ = (double)param_.getValue("windowsize");
79  peakcount_ = (UInt)param_.getValue("peakcount");
80 
81  //copy spectrum
82  SpectrumType old_spectrum = spectrum;
83  old_spectrum.sortByPosition();
84 
85  //find high peak positions
86  bool end = false;
87  std::set<double> positions;
88  for (ConstIterator it = old_spectrum.begin(); it != old_spectrum.end(); ++it)
89  {
90  // copy the window from the spectrum
91  SpectrumType window;
92  for (ConstIterator it2 = it; (it2->getPosition() - it->getPosition() < windowsize_); )
93  {
94  window.push_back(*it2);
95  if (++it2 == old_spectrum.end())
96  {
97  end = true;
98  break;
99  }
100  }
101 
102  //extract peakcount most intense peaks
103  window.sortByIntensity(true);
104  for (Size i = 0; i < peakcount_; ++i)
105  {
106  if (i < window.size())
107  {
108  positions.insert(window[i].getMZ());
109  }
110  }
111  //abort at the end of the spectrum
112  if (end) break;
113  }
114 
115  // select peaks that were retained
116  std::vector<Size> indices;
117  for (ConstIterator it = spectrum.begin(); it != spectrum.end(); ++it)
118  {
119  if (positions.find(it->getMZ()) != positions.end())
120  {
121  Size index(it - spectrum.begin());
122  indices.push_back(index);
123  }
124  }
125  spectrum.select(indices);
126  }
127 
128  void filterPeakSpectrum(PeakSpectrum& spectrum);
129 
130  void filterPeakMap(PeakMap& exp);
131 
132  // jumping window version (faster)
133  template <typename SpectrumType>
135  {
136  if (spectrum.empty())
137  {
138  return;
139  }
140 
141  spectrum.sortByPosition();
142 
143  windowsize_ = static_cast<double>(param_.getValue("windowsize"));
144  peakcount_ = static_cast<UInt>(param_.getValue("peakcount"));
145 
146  // copy meta data
147  SpectrumType out = spectrum;
148  out.clear(false);
149 
150  SpectrumType peaks_in_window;
151  double window_start = spectrum[0].getMZ();
152  for (Size i = 0; i != spectrum.size(); ++i)
153  {
154  if (spectrum[i].getMZ() - window_start < windowsize_) // collect peaks in window
155  {
156  peaks_in_window.push_back(spectrum[i]);
157  }
158  else // step over window boundaries
159  {
160  window_start = spectrum[i].getMZ(); // as there might be large gaps between peaks resulting in empty windows, set new window start to next peak
161 
162  // copy N highest peaks to out
163  if (peaks_in_window.size() > peakcount_)
164  {
165  std::partial_sort(peaks_in_window.begin(), peaks_in_window.begin() + peakcount_, peaks_in_window.end(), reverseComparator(typename SpectrumType::PeakType::IntensityLess()));
166  copy(peaks_in_window.begin(), peaks_in_window.begin() + peakcount_, back_inserter(out));
167  }
168  else
169  {
170  std::sort(peaks_in_window.begin(), peaks_in_window.end(), reverseComparator(typename SpectrumType::PeakType::IntensityLess()));
171  copy(peaks_in_window.begin(), peaks_in_window.end(), back_inserter(out));
172  }
173 
174  peaks_in_window.clear(false);
175  peaks_in_window.push_back(spectrum[i]);
176  }
177  }
178 
179  if (peaks_in_window.empty()) // last window is empty -> no special handling needed
180  {
181  // select peaks that were retained
182  std::vector<Size> indices;
183  for (typename SpectrumType::ConstIterator it = spectrum.begin(); it != spectrum.end(); ++it)
184  {
185  if (std::find(out.begin(), out.end(), *it) != out.end())
186  {
187  Size index(it - spectrum.begin());
188  indices.push_back(index);
189  }
190  }
191  spectrum.select(indices);
192  return;
193  }
194 
195  // Note that the last window might be much smaller than windowsize.
196  // Therefor the number of peaks copied from this window should be adapted accordingly.
197  // Otherwise a lot of noise peaks are copied from each end of a spectrum.
198 
199  double last_window_size = peaks_in_window.back().getMZ() - window_start;
200  double last_window_size_fraction = last_window_size / windowsize_;
201  Size last_window_peakcount = last_window_size_fraction * peakcount_;
202 
203  if (last_window_peakcount) // handle single peak in last window (will produce no proper fraction)
204  {
205  last_window_peakcount = 1;
206  }
207 
208  // sort for last_window_peakcount highest peaks
209  std::partial_sort(peaks_in_window.begin(), peaks_in_window.begin() + last_window_peakcount, peaks_in_window.end(), reverseComparator(typename SpectrumType::PeakType::IntensityLess()));
210 
211  if (peaks_in_window.size() > last_window_peakcount)
212  {
213  std::copy(peaks_in_window.begin(), peaks_in_window.begin() + last_window_peakcount, back_inserter(out));
214  }
215  else
216  {
217  std::copy(peaks_in_window.begin(), peaks_in_window.end(), std::back_inserter(out));
218  }
219 
220  // select peaks that were retained
221  std::vector<Size> indices;
222  for (typename SpectrumType::ConstIterator it = spectrum.begin(); it != spectrum.end(); ++it)
223  {
224  if (std::find(out.begin(), out.end(), *it) != out.end())
225  {
226  Size index(it - spectrum.begin());
227  indices.push_back(index);
228  }
229  }
230  spectrum.select(indices);
231 
232  return;
233  }
234 
235  //TODO reimplement DefaultParamHandler::updateMembers_()
236 
237 private:
238  double windowsize_;
240  };
241 
242 }
243 
244 
WindowMower augments the highest peaks in a sliding or jumping window.
Definition: WindowMower.h:54
void sortByPosition()
Lexicographically sorts the peaks by their position.
unsigned int UInt
Unsigned integer type.
Definition: Types.h:94
ContainerType::const_iterator ConstIterator
Non-mutable iterator.
Definition: MSSpectrum.h:103
ReverseComparator< Cmp > reverseComparator(Cmp const &cmp)
Make-function to create a ReverseComparator from another comparator without the need to specify the t...
Definition: ComparatorUtils.h:260
void filterPeakSpectrumForTopNInSlidingWindow(SpectrumType &spectrum)
sliding window version (slower)
Definition: WindowMower.h:74
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:46
bool find(TFinder &finder, const Pattern< TNeedle, FuzzyAC > &me, PatternAuxData< TNeedle > &dh)
Definition: AhoCorasickAmbiguous.h:884
The representation of a 1D spectrum.
Definition: MSSpectrum.h:66
void filterPeakSpectrumForTopNInJumpingWindow(SpectrumType &spectrum)
Definition: WindowMower.h:134
void sortByIntensity(bool reverse=false)
Lexicographically sorts the peaks by their intensity.
UInt peakcount_
Definition: WindowMower.h:239
void clear(bool clear_meta_data)
Clears all data and meta data.
In-Memory representation of a mass spectrometry experiment.
Definition: MSExperiment.h:77
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:127
double windowsize_
Definition: WindowMower.h:238
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:91
MSSpectrum & select(const std::vector< Size > &indices)