35 #ifndef OPENMS_TRANSFORMATIONS_FEATUREFINDER_EMGSCORING_H 
   36 #define OPENMS_TRANSFORMATIONS_FEATUREFINDER_EMGSCORING_H 
   39 #include <boost/math/special_functions/fpclassify.hpp>  
   80     template<
typename SpectrumType, 
class TransitionT>
 
   84       std::vector<double> fit_scores;
 
   86       bool smooth_data = 
false;
 
   96         fit_scores.push_back(fscore);
 
  100       avg_score /= transition_group.
size();
 
  109       if (current_section.size() < 2)
 
  116       typedef Peak1D LocalPeakType;
 
  120       std::vector<LocalPeakType> data_to_fit;
 
  121       prepareFit_(current_section, data_to_fit, smooth_data);
 
  132     template<
class LocalPeakType>
 
  154       if (boost::math::isnan(quality)) quality = -1.0;
 
  160     template<
class LocalPeakType>
 
  166       for (ConvexHull2D::PointArrayType::const_iterator it = current_section.begin(); it != current_section.end(); it++)
 
  170         p.setIntensity(it->getY());
 
  171         filter_spec.push_back(p);
 
  176       std::vector<DoubleReal> distances;
 
  177       for (
Size j = 1; j < filter_spec.size(); ++j)
 
  179         distances.push_back(filter_spec[j].getMZ() - filter_spec[j - 1].getMZ());
 
  181       DoubleReal dist_average = std::accumulate(distances.begin(), distances.end(), 0.0) / (
DoubleReal) distances.size();
 
  186       new_peak.
setMZ(filter_spec.back().getMZ() + dist_average);
 
  187       filter_spec.push_back(new_peak);
 
  188       new_peak.
setMZ(filter_spec.back().getMZ() + dist_average);
 
  189       filter_spec.push_back(new_peak);
 
  190       new_peak.
setMZ(filter_spec.back().getMZ() + dist_average);
 
  191       filter_spec.push_back(new_peak);
 
  194       new_peak.
setMZ(filter_spec.front().getMZ() - dist_average);
 
  195       filter_spec.insert(filter_spec.begin(), new_peak);
 
  196       new_peak.
setMZ(filter_spec.front().getMZ() - dist_average);
 
  197       filter_spec.insert(filter_spec.begin(), new_peak);
 
  198       new_peak.
setMZ(filter_spec.front().getMZ() - dist_average);
 
  199       filter_spec.insert(filter_spec.begin(), new_peak);
 
  208         filter_param.setValue(
"gaussian_width", 4 * dist_average);
 
  210         filter.
filter(filter_spec);
 
  214       for (
Size j = 0; j != filter_spec.size(); ++j)
 
  217         p.setPosition(filter_spec[j].getMZ());
 
  218         p.setIntensity(filter_spec[j].getIntensity());
 
  219         data_to_fit.push_back(p);
 
const Param & getDefaults() const 
Non-mutable access to the default parameters. 
A more convenient string class. 
Definition: String.h:56
Size size() const 
Definition: MRMTransitionGroup.h:107
Scoring of an elution peak using an exponentially modified gaussian distribution model. 
Definition: EmgScoring.h:60
const std::vector< SpectrumType > & getChromatograms() const 
Definition: MRMTransitionGroup.h:148
#define OPENMS_PRECONDITION(condition, message)
Precondition macro. 
Definition: Macros.h:107
std::vector< PointType > PointArrayType
Definition: ConvexHull2D.h:79
EmgFitter1D fitter_emg1D_
Definition: EmgScoring.h:223
Exponentially modified gaussian distribution fitter (1-dim.) using Levenberg-Marquardt algorithm (GSL...
Definition: EmgFitter1D.h:47
QualityType fit1d(const RawDataArrayType &range, InterpolationModel *&model)
return interpolation model 
void setMZ(CoordinateType mz)
Mutable access to m/z. 
Definition: Peak1D.h:114
void setIntensity(IntensityType intensity)
Mutable access to the data point intensity (height) 
Definition: Peak1D.h:105
void setParameters(const Param ¶m)
Sets the parameters. 
The representation of a transition group that has information about the individual chromatograms as w...
Definition: MRMTransitionGroup.h:56
Abstract class for 1D-models that are approximated using linear interpolation. 
Definition: InterpolationModel.h:55
const std::vector< ConvexHull2D > & getConvexHulls() const 
Non-mutable access to the convex hulls. 
A 1-dimensional raw data point or peak. 
Definition: Peak1D.h:55
void setFitterParam(Param param)
Definition: EmgScoring.h:69
An LC-MS feature. 
Definition: Feature.h:66
~EmgScoring()
Definition: EmgScoring.h:67
double elutionModelFit(ConvexHull2D::PointArrayType current_section, bool smooth_data)
Definition: EmgScoring.h:106
Management and storage of parameters / INI files. 
Definition: Param.h:69
This class represents a Gaussian lowpass-filter which works on uniform as well as on non-uniform prof...
Definition: GaussFilter.h:73
Feature & getFeature(String key)
get a specified feature 
double calcElutionFitScore(MRMFeature &mrmfeature, MRMTransitionGroup< SpectrumType, TransitionT > &transition_group)
calculate the elution profile fit score 
Definition: EmgScoring.h:81
size_t Size
Size type e.g. used as variable which can hold result of size() 
Definition: Types.h:144
void prepareFit_(const ConvexHull2D::PointArrayType ¤t_section, std::vector< LocalPeakType > &data_to_fit, bool smooth_data)
Definition: EmgScoring.h:161
EmgScoring()
Definition: EmgScoring.h:65
A multi-chromatogram MRM feature. 
Definition: MRMFeature.h:50
Param getDefaults()
Definition: EmgScoring.h:74
double fitRT_(std::vector< LocalPeakType > &rt_input_data, InterpolationModel *&model)
Definition: EmgScoring.h:133
void filter(MSSpectrum< PeakType > &spectrum)
Smoothes an MSSpectrum containing profile data. 
Definition: GaussFilter.h:92
const Param & getParameters() const 
Non-mutable access to the parameters. 
double DoubleReal
Double-precision real type. 
Definition: Types.h:118