Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages
EGHTraceFitter.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: Stephan Aiche $
33 // --------------------------------------------------------------------------
34 
35 #ifndef OPENMS_TRANSFORMATIONS_FEATUREFINDER_EGHTRACEFITTER_H
36 #define OPENMS_TRANSFORMATIONS_FEATUREFINDER_EGHTRACEFITTER_H
37 
38 #include <sstream>
41 
42 namespace OpenMS
43 {
44 
57  template <class PeakType>
59  public TraceFitter<PeakType>
60  {
61 public:
63  {
64  //setName("EGHTraceFitter");
65  }
66 
67  EGHTraceFitter(const EGHTraceFitter & other) :
68  TraceFitter<PeakType>(other)
69  {
70  this->height_ = other.height_;
71  this->apex_rt_ = other.apex_rt_;
72  this->sigma_square_ = other.sigma_square_;
73  this->tau_ = other.tau_;
74 
75  this->sigma_5_bound_ = other.sigma_5_bound_;
76  this->fwhm_bound_ = other.fwhm_bound_;
77 
79  }
80 
82  {
84 
85  this->height_ = source.height_;
86  this->apex_rt_ = source.apex_rt_;
87  this->sigma_square_ = source.sigma_square_;
88  this->tau_ = source.tau_;
89 
90  this->sigma_5_bound_ = source.sigma_5_bound_;
91  this->fwhm_bound_ = source.fwhm_bound_;
92 
94 
95  return *this;
96  }
97 
98  virtual ~EGHTraceFitter()
99  {
100  }
101 
102  // override important methods
104  {
105  setInitialParameters_(traces);
106 
107  double x_init[NUM_PARAMS_] = {height_, apex_rt_, sigma_square_, tau_};
108 
109  Size num_params = NUM_PARAMS_;
110 
111  TraceFitter<PeakType>::optimize_(traces, num_params, x_init,
115  }
116 
118  {
119  return sigma_5_bound_.first;
120  }
121 
123  {
124  return tau_;
125  }
126 
128  {
129  return sigma_5_bound_.second;
130  }
131 
133  {
134  return height_;
135  }
136 
138  {
139  return sigma_square_;
140  }
141 
143  {
144  return apex_rt_;
145  }
146 
147  bool checkMaximalRTSpan(const DoubleReal max_rt_span)
148  {
149  return (sigma_5_bound_.second - sigma_5_bound_.first) > max_rt_span * region_rt_span_;
150  }
151 
152  virtual bool checkMinimalRTSpan(const std::pair<DoubleReal, DoubleReal> & rt_bounds, const DoubleReal min_rt_span)
153  {
154  return (rt_bounds.second - rt_bounds.first) < min_rt_span * (sigma_5_bound_.second - sigma_5_bound_.first);
155  }
156 
158  {
159  double rt = trace.peaks[k].first;
160  double t_diff, t_diff2, denominator = 0.0;
161  double fegh = 0.0;
162 
163  t_diff = rt - apex_rt_;
164  t_diff2 = t_diff * t_diff; // -> (t - t_R)^2
165 
166  denominator = 2 * sigma_square_ + tau_ * t_diff; // -> 2\sigma_{g}^{2} + \tau \left(t - t_R\right)
167 
168  if (denominator > 0.0)
169  {
170  fegh = trace.theoretical_int * height_ * exp(-t_diff2 / denominator);
171  }
172 
173  return fegh;
174  }
175 
177  {
178  return height_ * (fwhm_bound_.second - fwhm_bound_.first);
179  }
180 
182  {
183 
184  std::pair<DoubleReal, DoubleReal> bounds = getAlphaBoundaries_(0.5);
185  return bounds.second - bounds.first;
186  }
187 
188  virtual String getGnuplotFormula(FeatureFinderAlgorithmPickedHelperStructs::MassTrace<PeakType> const & trace, const char function_name, const DoubleReal baseline, const DoubleReal rt_shift)
189  {
190  std::stringstream s;
191  s << String(function_name) << "(x)= " << baseline << " + ";
192  s << "("; // the overall bracket
193  s << "((" << 2 * sigma_square_ << " + " << tau_ << " * (x - " << (rt_shift + apex_rt_) << " )) > 0) ? "; // condition
194  s << (trace.theoretical_int * height_) << " * exp(-1 * (x - " << (rt_shift + apex_rt_) << ")**2 " <<
195  "/" <<
196  " ( " << 2 * sigma_square_ << " + " << tau_ << " * (x - " << (rt_shift + apex_rt_) << " )))";
197  s << " : 0)";
198  return String(s.str());
199  }
200 
201 protected:
204 
207 
208  std::pair<DoubleReal, DoubleReal> sigma_5_bound_;
209  std::pair<DoubleReal, DoubleReal> fwhm_bound_;
210 
212 
213  static const Size NUM_PARAMS_ = 4;
214 
220  std::pair<DoubleReal, DoubleReal> getAlphaBoundaries_(const DoubleReal alpha) const
221  {
222  std::pair<DoubleReal, DoubleReal> bounds;
223  DoubleReal L = log(alpha);
224  DoubleReal s = sqrt(
225  ((L * tau_) * (L * tau_) / 4) - 2 * L * sigma_square_
226  );
227 
228  DoubleReal s1, s2;
229  s1 = (-1 * (L * tau_) / 2) + s;
230  s2 = (-1 * (L * tau_) / 2) - s;
231 
232  // the smaller one (should be < 0) = lower bound
233  bounds.first = apex_rt_ + std::min(s1, s2);
234  // bigger one (should be > 0) = upper bound
235  bounds.second = apex_rt_ + std::max(s1, s2);
236 
237  return bounds;
238  }
239 
240  void getOptimizedParameters_(gsl_multifit_fdfsolver * fdfsolver)
241  {
242  height_ = gsl_vector_get(fdfsolver->x, 0);
243  apex_rt_ = gsl_vector_get(fdfsolver->x, 1);
244  sigma_square_ = gsl_vector_get(fdfsolver->x, 2);
245  tau_ = gsl_vector_get(fdfsolver->x, 3);
246 
247  // we set alpha to 0.04 which is conceptually equal to
248  // 2.5 sigma for lower and upper bound
250  // this is needed for the intensity contribution -> this is the 1.25 sigma region
251  fwhm_bound_ = getAlphaBoundaries_(0.45783);
252  }
253 
254  static Int residual_(const gsl_vector * param, void * data, gsl_vector * f)
255  {
257 
258  double H = gsl_vector_get(param, 0);
259  double tR = gsl_vector_get(param, 1);
260  double sigma_square = gsl_vector_get(param, 2);
261  double tau = gsl_vector_get(param, 3);
262 
263  double t_diff, t_diff2, denominator = 0.0;
264 
265  double fegh = 0.0;
266 
267  UInt count = 0;
268  for (Size t = 0; t < traces->size(); ++t)
269  {
271  for (Size i = 0; i < trace.peaks.size(); ++i)
272  {
273  DoubleReal rt = trace.peaks[i].first;
274 
275  t_diff = rt - tR;
276  t_diff2 = t_diff * t_diff; // -> (t - t_R)^2
277 
278  denominator = 2 * sigma_square + tau * t_diff; // -> 2\sigma_{g}^{2} + \tau \left(t - t_R\right)
279 
280  if (denominator > 0.0)
281  {
282  fegh = traces->baseline + trace.theoretical_int * H * exp(-t_diff2 / denominator);
283  }
284  else
285  {
286  fegh = 0.0;
287  }
288 
289  gsl_vector_set(f, count, (fegh - trace.peaks[i].second->getIntensity()));
290  ++count;
291  }
292  }
293  return GSL_SUCCESS;
294  }
295 
296  static Int jacobian_(const gsl_vector * param, void * data, gsl_matrix * J)
297  {
299 
300  double H = gsl_vector_get(param, 0);
301  double tR = gsl_vector_get(param, 1);
302  double sigma_square = gsl_vector_get(param, 2);
303  double tau = gsl_vector_get(param, 3);
304 
305  double derivative_H, derivative_tR, derivative_sigma_square, derivative_tau = 0.0;
306  double t_diff, t_diff2, exp1, denominator = 0.0;
307 
308  UInt count = 0;
309  for (Size t = 0; t < traces->size(); ++t)
310  {
312  for (Size i = 0; i < trace.peaks.size(); ++i)
313  {
314  DoubleReal rt = trace.peaks[i].first;
315 
316  t_diff = rt - tR;
317  t_diff2 = t_diff * t_diff; // -> (t - t_R)^2
318 
319  denominator = 2 * sigma_square + tau * t_diff; // -> 2\sigma_{g}^{2} + \tau \left(t - t_R\right)
320 
321  if (denominator > 0)
322  {
323  exp1 = exp(-t_diff2 / denominator);
324 
325  // \partial H f_{egh}(t) = \exp\left( \frac{-\left(t-t_R \right)}{2\sigma_{g}^{2} + \tau \left(t - t_R\right)} \right)
326  derivative_H = trace.theoretical_int * exp1;
327 
328  // \partial t_R f_{egh}(t) &=& H \exp \left( \frac{-\left(t-t_R \right)}{2\sigma_{g}^{2} + \tau \left(t - t_R\right)} \right) \left( \frac{\left( 4 \sigma_{g}^{2} + \tau \left(t-t_R \right) \right) \left(t-t_R \right)}{\left( 2\sigma_{g}^{2} + \tau \left(t - t_R\right) \right)^2} \right)
329  derivative_tR = trace.theoretical_int * H * exp1 * (((4 * sigma_square + tau * t_diff) * t_diff) / (denominator * denominator));
330 
331  // \partial \sigma_{g}^{2} f_{egh}(t) &=& H \exp \left( \frac{-\left(t-t_R \right)^2}{2\sigma_{g}^{2} + \tau \left(t - t_R\right)} \right) \left( \frac{ 2 \left(t - t_R\right)^2}{\left( 2\sigma_{g}^{2} + \tau \left(t - t_R\right) \right)^2} \right)
332  derivative_sigma_square = trace.theoretical_int * H * exp1 * ((2 * t_diff2) / (denominator * denominator));
333 
334  // \partial \tau f_{egh}(t) &=& H \exp \left( \frac{-\left(t-t_R \right)^2}{2\sigma_{g}^{2} + \tau \left(t - t_R\right)} \right) \left( \frac{ \left(t - t_R\right)^3}{\left( 2\sigma_{g}^{2} + \tau \left(t - t_R\right) \right)^2} \right)
335  derivative_tau = trace.theoretical_int * H * exp1 * ((t_diff * t_diff2) / (denominator * denominator));
336  }
337  else
338  {
339  derivative_H = 0.0;
340  derivative_tR = 0.0;
341  derivative_sigma_square = 0.0;
342  derivative_tau = 0.0;
343  }
344 
345  // set the jacobian matrix
346  gsl_matrix_set(J, count, 0, derivative_H);
347  gsl_matrix_set(J, count, 1, derivative_tR);
348  gsl_matrix_set(J, count, 2, derivative_sigma_square);
349  gsl_matrix_set(J, count, 3, derivative_tau);
350 
351  ++count;
352  }
353  }
354  return GSL_SUCCESS;
355  }
356 
357  static Int evaluate_(const gsl_vector * param, void * data, gsl_vector * f, gsl_matrix * J)
358  {
359  residual_(param, data, f);
360  jacobian_(param, data, J);
361  return GSL_SUCCESS;
362  }
363 
365  {
366  LOG_DEBUG << "EGHTraceFitter->setInitialParameters(..)" << std::endl;
367  LOG_DEBUG << "Traces length: " << traces.size() << std::endl;
368  LOG_DEBUG << "Max trace: " << traces.max_trace << std::endl;
369 
370  // initial values for externals
371  height_ = traces[traces.max_trace].max_peak->getIntensity() - traces.baseline;
372  LOG_DEBUG << "height: " << height_ << std::endl;
373  apex_rt_ = traces[traces.max_trace].max_rt;
374  LOG_DEBUG << "apex_rt: " << apex_rt_ << std::endl;
375  region_rt_span_ = traces[traces.max_trace].peaks.back().first - traces[traces.max_trace].peaks[0].first;
376  LOG_DEBUG << "region_rt_span_: " << region_rt_span_ << std::endl;
377 
378  const PeakType * max_peak = traces[traces.max_trace].peaks.begin()->second;
379  Size max_pos = 0;
380 
381  for (Size i = 1; i < traces[traces.max_trace].peaks.size(); ++i)
382  {
383  if (traces[traces.max_trace].peaks[i].second->getIntensity() > max_peak->getIntensity())
384  {
385  max_peak = traces[traces.max_trace].peaks[i].second;
386  max_pos = i;
387  }
388  }
389 
390  Size i = max_pos;
391  LOG_DEBUG << "max_pos: " << max_pos << std::endl;
392  if (traces[traces.max_trace].peaks.size() < 3)
393  {
394  // TODO: abort the whole thing here??
395  // because below we REQUIRE at least three peaks!!!
396  }
397 
398  Size filter_max_pos = traces[traces.max_trace].peaks.size() - 2;
399 
400  // compute a smoothed value for the maxima
401  // if the maximum is close to the borders, we need to think of something...
402  DoubleReal smoothed_height;
403  if ((max_pos < 2) || (max_pos + 2 >= traces[traces.max_trace].peaks.size()))
404  {
405  // ... too close to border... no smoothing
406  smoothed_height = traces[traces.max_trace].peaks[max_pos].second->getIntensity();
407  // TODO: does this trace even make sense?! why wasn't it extended it further? or should we have skipped it beforehand?
408  }
409  else
410  {
411  smoothed_height = (traces[traces.max_trace].peaks[max_pos - 2].second->getIntensity()
412  + traces[traces.max_trace].peaks[max_pos - 1].second->getIntensity()
413  + traces[traces.max_trace].peaks[max_pos].second->getIntensity()
414  + traces[traces.max_trace].peaks[max_pos + 1].second->getIntensity()
415  + traces[traces.max_trace].peaks[max_pos + 2].second->getIntensity()) / 5.0;
416  }
417 
418  // use moving average filter to avoid bad initial values
419  // moving average of size 5
420  // TODO: optimize windows size
421  while (i > 2 && i < filter_max_pos)
422  {
423  // compute smoothed
424  DoubleReal smoothed = (traces[traces.max_trace].peaks[i - 2].second->getIntensity()
425  + traces[traces.max_trace].peaks[i - 1].second->getIntensity()
426  + traces[traces.max_trace].peaks[i].second->getIntensity()
427  + traces[traces.max_trace].peaks[i + 1].second->getIntensity()
428  + traces[traces.max_trace].peaks[i + 2].second->getIntensity()) / 5.0;
429 
430  if (smoothed / smoothed_height < 0.5) break;
431  else --i;
432  }
433  LOG_DEBUG << "Left alpha at " << i << " with " << traces[traces.max_trace].peaks[i].first << std::endl;
434  double A = apex_rt_ - traces[traces.max_trace].peaks[i].first;
435 
436  i = max_pos;
437  while (i < filter_max_pos && i > 2)
438  {
439  DoubleReal smoothed = (traces[traces.max_trace].peaks[i - 2].second->getIntensity()
440  + traces[traces.max_trace].peaks[i - 1].second->getIntensity()
441  + traces[traces.max_trace].peaks[i].second->getIntensity()
442  + traces[traces.max_trace].peaks[i + 1].second->getIntensity()
443  + traces[traces.max_trace].peaks[i + 2].second->getIntensity()) / 5.0;
444 
445  if (smoothed / smoothed_height < 0.5) break;
446  else ++i;
447  }
448  LOG_DEBUG << "Right alpha at " << i << " with " << traces[traces.max_trace].peaks[i].first << std::endl;
449  double B = traces[traces.max_trace].peaks[i].first - apex_rt_;
450 
451  //LOG_DEBUG << "A: " << A << std::endl;
452  //LOG_DEBUG << "B: " << B << std::endl;
453 
454  // compute estimates for tau / sigma_square based on A/B
455  double log_alpha = log(0.5);
456 
457  tau_ = (-1 / log_alpha) * (B - A);
458  LOG_DEBUG << "tau: " << tau_ << std::endl;
459  sigma_square_ = (-1 / (2 * log_alpha)) * (B * A);
460  LOG_DEBUG << "sigma_square: " << sigma_square_ << std::endl;
461  }
462 
463  virtual void updateMembers_()
464  {
466  }
467 
468  void printState_(SignedSize iter, gsl_multifit_fdfsolver * s)
469  {
470  LOG_DEBUG << "iter: " << iter << " "
471  << "height: " << gsl_vector_get(s->x, 0) << " "
472  << "apex_rt: " << gsl_vector_get(s->x, 1) << " "
473  << "sigma_square: " << gsl_vector_get(s->x, 2) << " "
474  << "tau: " << gsl_vector_get(s->x, 3) << " "
475  << "|f(x)| = " << gsl_blas_dnrm2(s->f) << std::endl;
476  }
477 
478  };
479 
480 } // namespace OpenMS
481 
482 #endif // #ifndef OPENMS_TRANSFORMATIONS_FEATUREFINDER_FEATUREFINDERALGORITHMPICKEDTRACEFITTERGAUSS_H
const double k
DoubleReal computeTheoretical(const FeatureFinderAlgorithmPickedHelperStructs::MassTrace< PeakType > &trace, Size k)
Definition: EGHTraceFitter.h:157
A more convenient string class.
Definition: String.h:56
EGHTraceFitter()
Definition: EGHTraceFitter.h:62
DoubleReal tau_
Definition: EGHTraceFitter.h:206
void getOptimizedParameters_(gsl_multifit_fdfsolver *fdfsolver)
Definition: EGHTraceFitter.h:240
A 2-dimensional raw data point or peak.
Definition: Peak2D.h:55
DoubleReal getTau() const
Definition: EGHTraceFitter.h:122
Size max_trace
Maximum intensity trace.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:269
void optimize_(FeatureFinderAlgorithmPickedHelperStructs::MassTraces< PeakType > &traces, const Size num_params, double x_init[], Int(*residual)(const gsl_vector *x, void *params, gsl_vector *f), Int(*jacobian)(const gsl_vector *x, void *params, gsl_matrix *J), Int(*evaluate)(const gsl_vector *x, void *params, gsl_vector *f, gsl_matrix *J))
Definition: TraceFitter.h:199
virtual void updateMembers_()
This method is used to update extra member variables at the end of the setParameters() method...
Definition: TraceFitter.h:182
void fit(FeatureFinderAlgorithmPickedHelperStructs::MassTraces< PeakType > &traces)
Definition: EGHTraceFitter.h:103
IntensityType getIntensity() const
Definition: Peak2D.h:161
Helper struct for a collection of mass traces used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:153
Abstract fitter for RT profile fitting.
Definition: TraceFitter.h:62
static const Size NUM_PARAMS_
Definition: EGHTraceFitter.h:213
DoubleReal getFWHM() const
Definition: EGHTraceFitter.h:181
ptrdiff_t SignedSize
Signed Size type e.g. used as pointer difference.
Definition: Types.h:151
virtual DoubleReal getFeatureIntensityContribution()
Definition: EGHTraceFitter.h:176
std::pair< DoubleReal, DoubleReal > getAlphaBoundaries_(const DoubleReal alpha) const
Return an ordered pair of the positions where the EGH reaches a height of alpha * height of the EGH...
Definition: EGHTraceFitter.h:220
DoubleReal getSigmaSquare() const
Definition: EGHTraceFitter.h:137
#define LOG_DEBUG
Macro for general debugging information.
Definition: LogStream.h:459
DoubleReal height_
Definition: EGHTraceFitter.h:203
static Int jacobian_(const gsl_vector *param, void *data, gsl_matrix *J)
Definition: EGHTraceFitter.h:296
virtual TraceFitter & operator=(const TraceFitter &source)
assignment operator
Definition: TraceFitter.h:86
virtual ~EGHTraceFitter()
Definition: EGHTraceFitter.h:98
DoubleReal theoretical_int
Theoretical intensity value (scaled to [0,1])
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:91
std::pair< DoubleReal, DoubleReal > sigma_5_bound_
Definition: EGHTraceFitter.h:208
EGHTraceFitter & operator=(const EGHTraceFitter &source)
Definition: EGHTraceFitter.h:81
DoubleReal getUpperRTBound() const
Definition: EGHTraceFitter.h:127
std::pair< DoubleReal, DoubleReal > fwhm_bound_
Definition: EGHTraceFitter.h:209
DoubleReal getCenter() const
Definition: EGHTraceFitter.h:142
virtual void updateMembers_()
This method is used to update extra member variables at the end of the setParameters() method...
Definition: EGHTraceFitter.h:463
DoubleReal baseline
Estimated baseline in the region of the feature (used for the fit)
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:271
virtual bool checkMinimalRTSpan(const std::pair< DoubleReal, DoubleReal > &rt_bounds, const DoubleReal min_rt_span)
Definition: EGHTraceFitter.h:152
DoubleReal region_rt_span_
Definition: EGHTraceFitter.h:211
DoubleReal apex_rt_
Definition: EGHTraceFitter.h:202
void printState_(SignedSize iter, gsl_multifit_fdfsolver *s)
Definition: EGHTraceFitter.h:468
DoubleReal getHeight() const
Definition: EGHTraceFitter.h:132
DoubleReal sigma_square_
Definition: EGHTraceFitter.h:205
static Int residual_(const gsl_vector *param, void *data, gsl_vector *f)
Definition: EGHTraceFitter.h:254
bool checkMaximalRTSpan(const DoubleReal max_rt_span)
Definition: EGHTraceFitter.h:147
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
EGHTraceFitter(const EGHTraceFitter &other)
Definition: EGHTraceFitter.h:67
virtual String getGnuplotFormula(FeatureFinderAlgorithmPickedHelperStructs::MassTrace< PeakType > const &trace, const char function_name, const DoubleReal baseline, const DoubleReal rt_shift)
Definition: EGHTraceFitter.h:188
A RT Profile fitter using an Exponential Gaussian Hybrid background model.
Definition: EGHTraceFitter.h:58
void setInitialParameters_(FeatureFinderAlgorithmPickedHelperStructs::MassTraces< PeakType > &traces)
Definition: EGHTraceFitter.h:364
DoubleReal getLowerRTBound() const
Definition: EGHTraceFitter.h:117
int Int
Signed integer type.
Definition: Types.h:100
static Int evaluate_(const gsl_vector *param, void *data, gsl_vector *f, gsl_matrix *J)
Definition: EGHTraceFitter.h:357
Helper struct for mass traces used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:83

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