Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages
OfflinePrecursorIonSelection.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-2015.
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: Alexandra Zerck $
32 // $Authors: Alexandra Zerck $
33 // --------------------------------------------------------------------------
34 
35 #ifndef OPENMS_ANALYSIS_TARGETED_OFFLINEPRECURSORIONSELECTION_H
36 #define OPENMS_ANALYSIS_TARGETED_OFFLINEPRECURSORIONSELECTION_H
37 
38 
45 
46 namespace OpenMS
47 {
48  class PeptideIdentification;
49  class ProteinIdentification;
50  class String;
51 
52 
62  class OPENMS_DLLAPI OfflinePrecursorIonSelection :
63  public DefaultParamHandler
64  {
65 public:
67 
70 
80  template <typename InputPeakType>
81  void makePrecursorSelectionForKnownLCMSMap(const FeatureMap& features,
82  const MSExperiment<InputPeakType>& experiment,
84  std::set<Int>& charges_set,
85  bool feature_based);
86 
94  template <typename InputPeakType>
95  void getMassRanges(const FeatureMap& features,
96  const MSExperiment<InputPeakType>& experiment,
97  std::vector<std::vector<std::pair<Size, Size> > >& indices);
98 
99  void createProteinSequenceBasedLPInclusionList(String include, String rt_model_file, String pt_model_file, FeatureMap& precursors);
100 
102  {
103  solver_ = solver;
104  std::cout << " LPSolver set to " << solver_ << std::endl;
105  }
106 
108  {
109  return solver_;
110  }
111 
112 private:
116  template <typename InputPeakType>
117  void calculateXICs_(const FeatureMap& features,
118  const std::vector<std::vector<std::pair<Size, Size> > >& mass_ranges,
119  const MSExperiment<InputPeakType>& experiment,
120  const std::set<Int>& charges_set,
121  std::vector<std::vector<std::pair<Size, double> > >& xics);
122 
126  template <typename InputPeakType>
127  void checkMassRanges_(std::vector<std::vector<std::pair<Size, Size> > >& mass_ranges,
128  const MSExperiment<InputPeakType>& experiment);
129 
130  template <typename T>
131  void updateExclusionList_(std::vector<std::pair<T, Size> >& exclusion_list);
132 
133  void updateExclusionList_(std::map<std::pair<double, double>, Size, PairComparatorSecondElement<std::pair<double, double> > >& exclusion_list);
134 
136  };
137 
138  template <typename InputPeakType>
140  {
141  bool enclose_hit = false;
142  const std::vector<ConvexHull2D>& hulls = f.getConvexHulls();
143  for (Size i = 0; i < hulls.size(); ++i)
144  {
145  if (hulls[i].getBoundingBox().encloses(rt, mz))
146  {
147  enclose_hit = true;
148  return enclose_hit;
149  }
150  }
151  return enclose_hit;
152  }
153 
154  template <typename InputPeakType>
156  const MSExperiment<InputPeakType>& experiment,
157  std::vector<std::vector<std::pair<Size, Size> > >& indices)
158  {
159  if (experiment.empty())
160  throw Exception::InvalidSize(__FILE__, __LINE__, __PRETTY_FUNCTION__, 0);
161  for (Size f = 0; f < features.size(); ++f)
162  {
163  std::vector<std::pair<Size, Size> > vec;
164 
165  for (Size rt = 0; rt < experiment.size(); ++rt)
166  {
167  // is scan relevant?
168  if (!enclosesBoundingBox<InputPeakType>(features[f], experiment[rt].getRT(), features[f].getMZ()))
169  continue;
170 
171  std::pair<Size, Size> start;
172  std::pair<Size, Size> end;
173  bool start_found = false;
174  bool end_found = false;
175  typename MSSpectrum<InputPeakType>::ConstIterator mz_iter = experiment[rt].MZBegin(features[f].getMZ());
176  typename MSSpectrum<InputPeakType>::ConstIterator mz_end = mz_iter;
177  if (mz_iter == experiment[rt].end())
178  continue;
179  // check to the left
180  while (enclosesBoundingBox<InputPeakType>(features[f], experiment[rt].getRT(), mz_iter->getMZ()))
181  {
182  start_found = true;
183  start.first = rt;
184  start.second = distance(experiment[rt].begin(), mz_iter);
185  if (mz_iter == experiment[rt].begin())
186  break;
187  --mz_iter;
188  end_found = true;
189  end.first = rt;
190  end.second = start.second;
191  }
192  // and now to the right
193  while (mz_end != experiment[rt].end() && enclosesBoundingBox<InputPeakType>(features[f], experiment[rt].getRT(), mz_end->getMZ()))
194  {
195  end_found = true;
196  end.first = rt;
197  end.second = distance(experiment[rt].begin(), mz_end);
198  ++mz_end;
199  }
200  if (start_found && end_found)
201  {
202  vec.push_back(start);
203  vec.push_back(end);
204  }
205 #ifdef DEBUG_OPS
206  else if (start_found || end_found)
207  {
208  std::cout << "start " << start_found << " end " << end_found << std::endl;
209  std::cout << "feature: " << f << " rt: " << rt << std::endl;
210  }
211 #endif
212  }
213 #ifdef DEBUG_OPS
214  if (vec.size() > 0)
215  {
216  std::cout << vec.size() << " / 2 scans" << std::endl;
217  for (Size i = 0; i < vec.size(); i += 2)
218  {
219  std::cout << "Feature " << f << " RT : " << vec[i].first
220  << " MZ : " << experiment[vec[i].first][vec[i].second].getMZ() << " "
221  << experiment[vec[i + 1].first][vec[i + 1].second].getMZ() << std::endl;
222  }
223  }
224 #endif
225  if (vec.empty())
226  {
227 #ifdef DEBUG_OPS
228  std::cout << "According to the convex hulls no mass traces found for this feature->estimate!"
229  << features[f].getRT() << " " << features[f].getMZ() << " " << features[f].getCharge() << std::endl;
230 #endif
231  // we estimate the convex hull
232  typename MSExperiment<InputPeakType>::ConstIterator spec_iter = experiment.RTBegin(features[f].getRT());
233  if (spec_iter == experiment.end())
234  --spec_iter;
235 
236  double dist1 = fabs(spec_iter->getRT() - features[f].getRT());
237  double dist2 = std::numeric_limits<double>::max();
238  double dist3 = std::numeric_limits<double>::max();
239  if ((spec_iter + 1) != experiment.end())
240  {
241  dist2 = fabs((spec_iter + 1)->getRT() - features[f].getRT());
242  }
243  if (spec_iter != experiment.begin())
244  {
245  dist3 = fabs((spec_iter - 1)->getRT() - features[f].getRT());
246  }
247  if (dist3 <= dist1 && dist3 <= dist2)
248  {
249  --spec_iter;
250  }
251  else if (dist2 <= dist3 && dist2 <= dist1)
252  {
253  ++spec_iter;
254  }
255  std::pair<Size, Size> start;
256  std::pair<Size, Size> end;
257  start.first = distance(experiment.begin(), spec_iter);
258  end.first = start.first;
259 
260  typename MSSpectrum<InputPeakType>::ConstIterator mz_iter = spec_iter->MZBegin(features[f].getMZ());
261  if (spec_iter->begin() == spec_iter->end())
262  {
263  indices.push_back(vec);
264  continue;
265  }
266  if (mz_iter == spec_iter->end() || (mz_iter->getMZ() > features[f].getMZ() && mz_iter != spec_iter->begin()))
267  --mz_iter;
268  while (mz_iter != spec_iter->begin())
269  {
270  if (fabs((mz_iter - 1)->getMZ() - features[f].getMZ()) < 0.5)
271  --mz_iter;
272  else
273  break;
274  }
275  start.second = distance(spec_iter->begin(), mz_iter);
276  typename MSSpectrum<InputPeakType>::ConstIterator mz_end = mz_iter;
277 #ifdef DEBUG_OPS
278  std::cout << features[f].getMZ() << " Start: " << experiment[start.first].getRT() << " " << experiment[start.first][start.second].getMZ();
279 #endif
280  Int charge = features[f].getCharge();
281  if (charge == 0)
282  charge = 1;
283  while (mz_end + 1 != spec_iter->end())
284  {
285  if (fabs((mz_end + 1)->getMZ() - features[f].getMZ()) < 3.0 / (double)charge)
286  ++mz_end;
287  else
288  break;
289  }
290  end.second = distance(spec_iter->begin(), mz_end);
291 #ifdef DEBUG_OPS
292  std::cout << "\tEnd: " << experiment[end.first].getRT() << " " << experiment[end.first][end.second].getMZ() << std::endl;
293 #endif
294  vec.push_back(start);
295  vec.push_back(end);
296  }
297 
298  indices.push_back(vec);
299  }
300  // eliminate nearby peaks
301  if (param_.getValue("exclude_overlapping_peaks") == "true")
302  checkMassRanges_(indices, experiment);
303  }
304 
305  template <typename InputPeakType>
307  const std::vector<std::vector<std::pair<Size, Size> > >& mass_ranges,
308  const MSExperiment<InputPeakType>& experiment,
309  const std::set<Int>& charges_set,
310  std::vector<std::vector<std::pair<Size, double> > >& xics)
311  {
312  xics.clear();
313  xics.resize(experiment.size());
314  // for each feature
315  for (Size f = 0; f < mass_ranges.size(); ++f)
316  {
317  // is charge valid
318  if (charges_set.count(features[f].getCharge()) < 1)
319  {
320  continue;
321  }
322  // go through all scans where the feature occurs
323  for (Size s = 0; s < mass_ranges[f].size(); s += 2)
324  {
325  // sum intensity over all raw datapoints belonging to the feature in the current scan
326  double weight = 0.;
327  for (Size j = mass_ranges[f][s].second; j <= mass_ranges[f][s + 1].second; ++j)
328  {
329  weight += experiment[mass_ranges[f][s].first][j].getIntensity();
330  }
331  // enter xic in the vector for scan
332  xics[mass_ranges[f][s].first].push_back(std::make_pair(f, weight));
333  }
334  }
335 
336  for (Size s = 0; s < xics.size(); ++s)
337  {
338  sort(xics[s].begin(), xics[s].end(), PairComparatorSecondElement<std::pair<Size, double> >());
339  }
340  }
341 
342  template <typename InputPeakType>
344  const MSExperiment<InputPeakType>& experiment,
346  std::set<Int>& charges_set,
347  bool feature_based)
348  {
349 
350  const double window = param_.getValue("selection_window");
351  const double excl_window = param_.getValue("min_peak_distance");
352 
353  // get the mass ranges for each features for each scan it occurs in
354  std::vector<std::vector<std::pair<Size, Size> > > indices;
355  getMassRanges(features, experiment, indices);
356  double rt_dist = 0.;
357  if (experiment.size() > 1)
358  {
359  rt_dist = experiment[1].getRT() - experiment[0].getRT();
360  }
361 
362  // feature based selection (e.g. with LC-MALDI)
363  if (feature_based)
364  {
365  // create ILP
366  PSLPFormulation ilp_wrapper;
367 
368  std::vector<IndexTriple> variable_indices;
369  std::vector<int> solution_indices;
370  ilp_wrapper.createAndSolveILPForKnownLCMSMapFeatureBased(features, experiment, variable_indices,
371  indices, charges_set,
372  param_.getValue("ms2_spectra_per_rt_bin"),
373  solution_indices);
374 
375  sort(variable_indices.begin(), variable_indices.end(), PSLPFormulation::IndexLess());
376 #ifdef DEBUG_OPS
377  std::cout << "best_solution " << std::endl;
378 #endif
379  // print best solution
380  // create inclusion list
381  for (Size i = 0; i < solution_indices.size(); ++i)
382  {
383  Size feature_index = variable_indices[solution_indices[i]].feature;
384  Size feature_scan_idx = variable_indices[solution_indices[i]].scan;
385  typename MSExperiment<InputPeakType>::ConstIterator scan = experiment.begin() + feature_scan_idx;
387  Precursor p;
388  std::vector<Precursor> pcs;
389  p.setIntensity(features[feature_index].getIntensity());
390  p.setMZ(features[feature_index].getMZ());
391  p.setCharge(features[feature_index].getCharge());
392  pcs.push_back(p);
393  ms2_spec.setPrecursors(pcs);
394  ms2_spec.setRT(scan->getRT() + rt_dist / 2.0);
395  ms2_spec.setMSLevel(2);
396  // link ms2 spectrum with features overlapping its precursor
397  // Warning: this depends on the current order of features in the map
398  // Attention: make sure to name ALL features that overlap, not only one!
399  ms2_spec.setMetaValue("parent_feature_ids", ListUtils::create<Int>(String(feature_index)));
400  ms2.addSpectrum(ms2_spec);
401  std::cout << " MS2 spectra generated at: " << scan->getRT() << " x " << p.getMZ() << "\n";
402 
403  }
404 #ifdef DEBUG_OPS
405  std::cout << solution_indices.size() << " out of " << features.size()
406  << " precursors are in best solution.\n";
407 #endif
408  }
409  else // scan based selection (take the x highest signals for each spectrum)
410  {
411 #ifdef DEBUG_OPS
412  std::cout << "scan based precursor selection" << std::endl;
413 #endif
414  // if the highest signals for each scan shall be selected we don't need an ILP formulation
415 
416  //cache the values for each feature
417  std::vector<DoubleList> feature_elution_bounds;
418  std::vector<DoubleList> elution_profile_intensities;
419  std::vector<DoubleList> isotope_intensities;
420 
421  bool meta_values_present = false;
422 
423  if (!features.empty() &&
424  features[0].metaValueExists("elution_profile_bounds") &&
425  features[0].metaValueExists("elution_profile_intensities") &&
426  features[0].metaValueExists("isotope_intensities"))
427  {
428  for (Size feat = 0; feat < features.size(); ++feat)
429  {
430  feature_elution_bounds.push_back(features[feat].getMetaValue("elution_profile_bounds"));
431  elution_profile_intensities.push_back(features[feat].getMetaValue("elution_profile_intensities"));
432  isotope_intensities.push_back(features[feat].getMetaValue("isotope_intensities"));
433  }
434  meta_values_present = true;
435  }
436 
437  //for each feature cache for which scans it has to be considered
438  std::vector<std::vector<Size> > scan_features(experiment.size());
439 
440  for (Size feat = 0; feat < features.size(); ++feat)
441  {
442  if (charges_set.count(features[feat].getCharge()))
443  {
444  Size lower_rt = features[feat].getConvexHull().getBoundingBox().minX();
445  Size upper_rt = features[feat].getConvexHull().getBoundingBox().maxX();
447  for (it = experiment.RTBegin(lower_rt); it != experiment.RTEnd(upper_rt); ++it)
448  {
449  scan_features[it - experiment.begin()].push_back(feat);
450  }
451  }
452  }
453 
454  bool dynamic_exclusion = param_.getValue("Exclusion:use_dynamic_exclusion") == "true" ? true : false;
455  typedef std::map<std::pair<double, double>, Size, PairComparatorSecondElement<std::pair<double, double> > > ExclusionListType;
456  ExclusionListType exclusion_list;
457  Size exclusion_specs = (Size)(floor((double)param_.getValue("Exclusion:exclusion_time") / (double) rt_dist));
458  if (!dynamic_exclusion)
459  {
460  //if the dynamic exclusion if not active we use the exclusion list to guarantee no two peaks within min_peak_distance are selected for single scan
461  exclusion_specs = 0;
462  }
463 
464  //cache bounding boxes of features and mass traces (mass trace bb are also widened for effective discovery of enclosing peaks in intervals)
465  std::map<Size, typename OpenMS::DBoundingBox<2> > bounding_boxes_f;
466  std::map<std::pair<Size, Size>, typename OpenMS::DBoundingBox<2> > bounding_boxes;
467  for (Size feature_num = 0; feature_num < features.size(); ++feature_num)
468  {
469  if (charges_set.count(features[feature_num].getCharge()))
470  {
471  bounding_boxes_f.insert(std::make_pair(feature_num, features[feature_num].getConvexHull().getBoundingBox()));
472  const std::vector<ConvexHull2D> mass_traces = features[feature_num].getConvexHulls();
473  for (Size mass_trace_num = 0; mass_trace_num < mass_traces.size(); ++mass_trace_num)
474  {
475  typename OpenMS::DBoundingBox<2> tmp_bbox = mass_traces[mass_trace_num].getBoundingBox();
476  tmp_bbox.setMinY(tmp_bbox.minY() - window);
477  tmp_bbox.setMaxY(tmp_bbox.maxY() + window);
478  bounding_boxes.insert(std::make_pair(std::make_pair(feature_num, mass_trace_num), tmp_bbox));
479  }
480  }
481  }
482 
483  Size max_spec = (Int)param_.getValue("ms2_spectra_per_rt_bin");
484  // get best x signals for each scan
485  for (Size i = 0; i < experiment.size(); ++i)
486  {
487 #ifdef DEBUG_OPS
488  std::cout << "scan " << experiment[i].getRT() << ":";
489 #endif
490 
491  updateExclusionList_(exclusion_list);
492  MSSpectrum<InputPeakType> scan = experiment[i];
493  scan.sortByIntensity(true);
494  Size selected_peaks = 0, j = 0;
495 
496  while (selected_peaks < max_spec && j < scan.size())
497  {
498  double peak_mz = scan[j].getMZ();
499  double peak_rt = scan.getRT();
500 
501  ExclusionListType::iterator it_low = exclusion_list.lower_bound(std::make_pair(peak_mz, peak_mz));
502  if (it_low != exclusion_list.end() && it_low->first.first <= peak_mz)
503  {
504  ++j;
505  continue;
506  }
507  ++selected_peaks;
508 
509  //find all features (mass traces that are in the window around peak_mz)
511  std::vector<Precursor> pcs;
512  std::set<std::pair<Size, Size> > selected_mt;
513  IntList parent_feature_ids;
514 
515  double local_mz = peak_mz;
516  //std::cerr<<"MZ pos: "<<local_mz<<std::endl;
517  for (Size scan_feat_id = 0; scan_feat_id < scan_features[i].size(); ++scan_feat_id)
518  {
519  Size feature_num = scan_features[i][scan_feat_id];
520  if (bounding_boxes_f[feature_num].encloses(peak_rt, local_mz))
521  {
522  //find a mass trace enclosing the point
523  double feature_intensity = 0;
524  for (Size mass_trace_num = 0; mass_trace_num < features[feature_num].getConvexHulls().size(); ++mass_trace_num)
525  {
526  if (bounding_boxes[std::make_pair(feature_num, mass_trace_num)].encloses(DPosition<2>(peak_rt, local_mz)))
527  {
528  double elu_factor = 1.0, iso_factor = 1.0;
529  //get the intensity factor for the position in the elution profile
530  if (meta_values_present)
531  {
532  DoubleList xxx = elution_profile_intensities[feature_num];
533  DoubleList yyy = feature_elution_bounds[feature_num];
534 // std::cout << "PEAKRT: " << peak_rt << std::endl;
535 // std::cout << "Max: " << yyy[3] << " vs. " << bounding_boxes_f[feature_num].maxX() << std::endl;
536 // std::cout << "Min: " << yyy[1] << " vs. " << bounding_boxes_f[feature_num].minX() << std::endl;
537  OPENMS_PRECONDITION(i - yyy[0] < xxx.size(), "Tried to access invalid index for elution factor");
538  elu_factor = xxx[i - yyy[0]]; // segfault here: "i-yyy[0]" yields invalid index
539  iso_factor = isotope_intensities[feature_num][mass_trace_num];
540  }
541  feature_intensity += features[feature_num].getIntensity() * iso_factor * elu_factor;
542  }
543  }
544  Precursor p;
545  p.setIntensity(feature_intensity);
546  p.setMZ(features[feature_num].getMZ());
547  p.setCharge(features[feature_num].getCharge());
548  pcs.push_back(p);
549  parent_feature_ids.push_back((Int)feature_num);
550  }
551  }
552 
553  if (!pcs.empty())
554  {
555  //std::cerr<<"scan "<<i<<" added spectrum for features: "<<parent_feature_ids<<std::endl;
556  ms2_spec.setPrecursors(pcs);
557  ms2_spec.setMSLevel(2);
558  ms2_spec.setRT(experiment[i].getRT() + rt_dist / 2.0); //(selected_peaks+1)*rt_dist/(max_spec+1) );
559  ms2_spec.setMetaValue("parent_feature_ids", parent_feature_ids);
560  ms2.addSpectrum(ms2_spec);
561  }
562 
563  //add m/z window to exclusion list
564  exclusion_list.insert(std::make_pair(std::make_pair(peak_mz - excl_window, peak_mz + excl_window), exclusion_specs + 1));
565 
566  ++j;
567  }
568  }
569  }
570  }
571 
572  template <typename InputPeakType>
573  void OfflinePrecursorIonSelection::checkMassRanges_(std::vector<std::vector<std::pair<Size, Size> > >& mass_ranges,
574  const MSExperiment<InputPeakType>& experiment)
575  {
576  std::vector<std::vector<std::pair<Size, Size> > > checked_mass_ranges;
577  double min_peak_distance = param_.getValue("min_peak_distance");
578  checked_mass_ranges.reserve(mass_ranges.size());
579  for (Size f = 0; f < mass_ranges.size(); ++f)
580  {
581  std::vector<std::pair<Size, Size> > checked_mass_ranges_f;
582  for (Size s_idx = 0; s_idx < mass_ranges[f].size(); s_idx += 2)
583  {
584  Size s = mass_ranges[f][s_idx].first;
585  bool overlapping_features = false;
587  // check if other features overlap with this feature in the current scan
589  const InputPeakType& peak_left_border = experiment[s][mass_ranges[f][s_idx].second];
590  const InputPeakType& peak_right_border = experiment[s][mass_ranges[f][s_idx + 1].second];
591  for (Size fmr = 0; fmr < mass_ranges.size(); ++fmr)
592  {
593  if (fmr == f)
594  continue;
595  for (Size mr = 0; mr < mass_ranges[fmr].size(); mr += 2)
596  {
597  if (mass_ranges[fmr][mr].first == s) // same spectrum
598  {
599  const InputPeakType& tmp_peak_left = experiment[s][mass_ranges[fmr][mr].second];
600  const InputPeakType& tmp_peak_right = experiment[s][mass_ranges[fmr][mr + 1].second];
601 #ifdef DEBUG_OPS
602  std::cout << tmp_peak_left.getMZ() << " < "
603  << peak_left_border.getMZ() - min_peak_distance << " && "
604  << tmp_peak_right.getMZ() << " < "
605  << peak_left_border.getMZ() - min_peak_distance << " ? "
606  << (tmp_peak_left.getMZ() < peak_left_border.getMZ() - min_peak_distance &&
607  tmp_peak_right.getMZ() < peak_left_border.getMZ() - min_peak_distance)
608  << " || "
609  << tmp_peak_left.getMZ() << " > "
610  << peak_right_border.getMZ() + min_peak_distance << " && "
611  << tmp_peak_right.getMZ() << " > "
612  << peak_right_border.getMZ() + min_peak_distance << " ? "
613  << (tmp_peak_left.getMZ() > peak_right_border.getMZ() + min_peak_distance &&
614  tmp_peak_right.getMZ() > peak_right_border.getMZ() + min_peak_distance)
615  << std::endl;
616 #endif
617  // all other features have to be either completely left or
618  // right of the current feature
619  if (!((tmp_peak_left.getMZ() < peak_left_border.getMZ() - min_peak_distance &&
620  tmp_peak_right.getMZ() < peak_left_border.getMZ() - min_peak_distance) ||
621  (tmp_peak_left.getMZ() > peak_right_border.getMZ() + min_peak_distance &&
622  tmp_peak_right.getMZ() > peak_right_border.getMZ() + min_peak_distance)))
623  {
624 #ifdef DEBUG_OPS
625  std::cout << "found overlapping peak" << std::endl;
626 #endif
627  overlapping_features = true;
628  break;
629  }
630  }
631  }
632  }
633  if (!overlapping_features)
634  {
635 #ifdef DEBUG_OPS
636  std::cout << "feature in spec ok" << mass_ranges[f][s_idx].second << " in spec "
637  << mass_ranges[f][s_idx].first << std::endl;
638 #endif
639  checked_mass_ranges_f.insert(checked_mass_ranges_f.end(),
640  mass_ranges[f].begin() + s_idx,
641  mass_ranges[f].begin() + s_idx + 2);
642  }
643  }
644  checked_mass_ranges.push_back(checked_mass_ranges_f);
645  }
646  mass_ranges.swap(checked_mass_ranges);
647  }
648 
649  template <typename T>
650  void OfflinePrecursorIonSelection::updateExclusionList_(std::vector<std::pair<T, Size> >& exclusion_list)
651  {
652  for (Size i = 0; i < exclusion_list.size(); ++i)
653  {
654  if (exclusion_list[i].second > 0)
655  --exclusion_list[i].second;
656  }
657  sort(exclusion_list.begin(), exclusion_list.end(), PairComparatorSecondElementMore<std::pair<T, Size> >());
658  typename std::vector<std::pair<T, Size> >::iterator iter = exclusion_list.begin();
659  while (iter != exclusion_list.end() && iter->second != 0)
660  ++iter;
661  exclusion_list.erase(iter, exclusion_list.end());
662  }
663 
664  inline void OfflinePrecursorIonSelection::updateExclusionList_(std::map<std::pair<double, double>, Size, PairComparatorSecondElement<std::pair<double, double> > >& exclusion_list)
665  {
666  std::map<std::pair<double, double>, Size, PairComparatorSecondElement<std::pair<double, double> > >::iterator it;
667 
668  it = exclusion_list.begin();
669 
670  while (it != exclusion_list.end())
671  {
672  if ((it->second--) == 1)
673  {
674  exclusion_list.erase(it++);
675  }
676  else
677  {
678  ++it;
679  }
680  }
681  }
682 
683 }
684 
685 #endif // OPENMS_ANALYSIS_ID_OFFLINEPRECURSORIONSELECTION_H
CoordinateType maxY() const
Accessor for max_ coordinate maximum.
Definition: DIntervalBase.h:259
void setMetaValue(const String &name, const DataValue &value)
sets the DataValue corresponding to a name
A more convenient string class.
Definition: String.h:57
Precursor meta information.
Definition: Precursor.h:56
void setMinY(CoordinateType const c)
Mutator for max_ coordinate of the smaller point.
Definition: DIntervalBase.h:272
Size size() const
Definition: MSExperiment.h:117
std::vector< double > DoubleList
Vector of double precision real types.
Definition: ListUtils.h:66
Param param_
Container for current parameters.
Definition: DefaultParamHandler.h:150
SOLVER
Definition: LPWrapper.h:129
Definition: PSLPFormulation.h:141
A container for features.
Definition: FeatureMap.h:93
#define OPENMS_PRECONDITION(condition, message)
Precondition macro.
Definition: openms/include/OpenMS/CONCEPT/Macros.h:107
void setMaxY(CoordinateType const c)
Mutator for max_ coordinate of the larger point.
Definition: DIntervalBase.h:286
LPWrapper::SOLVER getLPSolver()
Definition: OfflinePrecursorIonSelection.h:107
CoordinateType getMZ() const
Non-mutable access to m/z.
Definition: Peak1D.h:114
ContainerType::const_iterator ConstIterator
Non-mutable iterator.
Definition: MSSpectrum.h:123
Iterator begin()
Definition: MSExperiment.h:147
bool enclosesBoundingBox(const Feature &f, typename MSExperiment< InputPeakType >::CoordinateType rt, typename MSExperiment< InputPeakType >::CoordinateType mz)
Definition: OfflinePrecursorIonSelection.h:139
Struct that holds the indices of the precursors in the feature map and the ilp formulation.
Definition: PSLPFormulation.h:69
std::vector< Int > IntList
Vector of signed integers.
Definition: ListUtils.h:59
void createAndSolveILPForKnownLCMSMapFeatureBased(const FeatureMap &features, const MSExperiment< InputPeakType > &experiment, std::vector< IndexTriple > &variable_indices, std::vector< std::vector< std::pair< Size, Size > > > &mass_ranges, std::set< Int > &charges_set, UInt ms2_spectra_per_rt_bin, std::vector< int > &solution_indices)
Encode ILP formulation for a given LC-MS map, but unknown protein sample.
Definition: PSLPFormulation.h:292
void setMZ(CoordinateType mz)
Mutable access to m/z.
Definition: Peak1D.h:120
void setIntensity(IntensityType intensity)
Mutable access to the data point intensity (height)
Definition: Peak1D.h:111
ConstIterator RTEnd(CoordinateType rt) const
Fast search for spectrum range end (returns the past-the-end iterator)
Definition: MSExperiment.h:388
Class for comparison of std::pair using second ONLY e.g. for use with std::sort.
Definition: ComparatorUtils.h:368
void makePrecursorSelectionForKnownLCMSMap(const FeatureMap &features, const MSExperiment< InputPeakType > &experiment, MSExperiment< InputPeakType > &ms2, std::set< Int > &charges_set, bool feature_based)
Makes the precursor selection for a given feature map, either feature or scan based.
Definition: OfflinePrecursorIonSelection.h:343
Iterator end()
Definition: MSExperiment.h:157
ConstIterator RTBegin(CoordinateType rt) const
Fast search for spectrum range begin.
Definition: MSExperiment.h:374
CoordinateType minY() const
Accessor for max_ coordinate minimum.
Definition: DIntervalBase.h:247
void getMassRanges(const FeatureMap &features, const MSExperiment< InputPeakType > &experiment, std::vector< std::vector< std::pair< Size, Size > > > &indices)
Calculates the mass ranges for each feature and stores them as indices of the raw data...
Definition: OfflinePrecursorIonSelection.h:155
const DataValue & getValue(const String &key) const
Returns a value of a parameter.
void calculateXICs_(const FeatureMap &features, const std::vector< std::vector< std::pair< Size, Size > > > &mass_ranges, const MSExperiment< InputPeakType > &experiment, const std::set< Int > &charges_set, std::vector< std::vector< std::pair< Size, double > > > &xics)
Calculate the sum of intensities of relevant features for each scan separately.
Definition: OfflinePrecursorIonSelection.h:306
void updateExclusionList_(std::vector< std::pair< T, Size > > &exclusion_list)
Definition: OfflinePrecursorIonSelection.h:650
const std::vector< ConvexHull2D > & getConvexHulls() const
Non-mutable access to the convex hulls.
void sortByIntensity(bool reverse=false)
Lexicographically sorts the peaks by their intensity.
Definition: MSSpectrum.h:342
LPWrapper::SOLVER solver_
Definition: OfflinePrecursorIonSelection.h:135
double getRT() const
Definition: MSSpectrum.h:243
void setMSLevel(UInt ms_level)
Sets the MS level.
Definition: MSSpectrum.h:265
An LC-MS feature.
Definition: Feature.h:70
void setLPSolver(LPWrapper::SOLVER solver)
Definition: OfflinePrecursorIonSelection.h:101
void addSpectrum(const MSSpectrum< PeakT > &spectrum)
adds a spectra to the list
Definition: MSExperiment.h:758
void setRT(double rt)
Sets the absolute retention time (is seconds)
Definition: MSSpectrum.h:249
Invalid UInt exception.
Definition: Exception.h:304
void checkMassRanges_(std::vector< std::vector< std::pair< Size, Size > > > &mass_ranges, const MSExperiment< InputPeakType > &experiment)
Eliminates overlapping peaks.
Definition: OfflinePrecursorIonSelection.h:573
In-Memory representation of a mass spectrometry experiment.
Definition: MSExperiment.h:69
std::vector< SpectrumType >::const_iterator ConstIterator
Non-mutable iterator.
Definition: MSExperiment.h:103
Implements different algorithms for precursor ion selection.
Definition: OfflinePrecursorIonSelection.h:62
void setPrecursors(const std::vector< Precursor > &precursors)
sets the precursors
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:121
Class for comparison of std::pair using second ONLY e.g. for use with std::sort.
Definition: ComparatorUtils.h:340
bool empty() const
Definition: MSExperiment.h:127
void setCharge(Int charge)
Mutable access to the charge.
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:92
PSLPFormulation::IndexTriple IndexTriple
Definition: OfflinePrecursorIonSelection.h:66
int Int
Signed integer type.
Definition: Types.h:96
Implements ILP formulation of precursor selection problems.
Definition: PSLPFormulation.h:54

OpenMS / TOPP release 2.0.0 Documentation generated on Wed Mar 30 2016 16:18:40 using doxygen 1.8.5