OpenMS  2.8.0
TraceFitter.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-2021.
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: Timo Sachsenberg$
32 // $Authors: Stephan Aiche, Marc Sturm $
33 // --------------------------------------------------------------------------
34 
35 #pragma once
36 
39 #include <OpenMS/KERNEL/Peak1D.h>
40 
41 // forward decl
42 namespace Eigen
43 {
44  template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
45  class Matrix;
46  using MatrixXd = Matrix<double, -1, -1, 0, -1, -1>;
47  using VectorXd = Matrix<double, -1, 1, 0, -1, 1>;
48 }
49 
50 namespace OpenMS
51 {
52 
62  class OPENMS_DLLAPI TraceFitter :
63  public DefaultParamHandler
64  {
65 
66 public:
68  //TODO: This is copy and paste from LevMarqFitter1d.h. Make a generic wrapper for LM optimization
70  {
71 public:
72  int inputs() const;
73  int values() const;
74 
75  GenericFunctor(int dimensions, int num_data_points);
76 
77  virtual ~GenericFunctor();
78 
79  virtual int operator()(const Eigen::VectorXd& x, Eigen::VectorXd& fvec) = 0;
80  // compute Jacobian matrix for the different parameters
81  virtual int df(const Eigen::VectorXd& x, Eigen::MatrixXd& J) = 0;
82 
83 protected:
84  const int m_inputs, m_values;
85  };
86 
89 
91  TraceFitter(const TraceFitter& source);
92 
95 
97  ~TraceFitter() override;
98 
103 
107  virtual double getLowerRTBound() const = 0;
108 
112  virtual double getUpperRTBound() const = 0;
113 
117  virtual double getHeight() const = 0;
118 
122  virtual double getCenter() const = 0;
123 
127  virtual double getFWHM() const = 0;
128 
132  virtual double getValue(double rt) const = 0;
133 
141 
148  virtual bool checkMinimalRTSpan(const std::pair<double, double>& rt_bounds, const double min_rt_span) = 0;
149 
155  virtual bool checkMaximalRTSpan(const double max_rt_span) = 0;
156 
160  virtual double getArea() = 0;
161 
171  virtual String getGnuplotFormula(const FeatureFinderAlgorithmPickedHelperStructs::MassTrace& trace, const char function_name, const double baseline, const double rt_shift) = 0;
172 
173 protected:
174  struct ModelData
175  {
177  bool weighted;
178  };
179 
180  void updateMembers_() override;
181 
187  virtual void getOptimizedParameters_(const Eigen::VectorXd&) = 0;
191  void optimize_(Eigen::VectorXd& x_init, GenericFunctor& functor);
192 
196  bool weighted_;
197 
198  };
199 
200 }
201 
A base class for all classes handling default parameters.
Definition: DefaultParamHandler.h:93
A more convenient string class.
Definition: String.h:60
Definition: TraceFitter.h:70
virtual int df(const Eigen::VectorXd &x, Eigen::MatrixXd &J)=0
GenericFunctor(int dimensions, int num_data_points)
virtual int operator()(const Eigen::VectorXd &x, Eigen::VectorXd &fvec)=0
const int m_inputs
Definition: TraceFitter.h:84
Abstract fitter for RT profile fitting.
Definition: TraceFitter.h:64
TraceFitter()
default constructor
virtual double getLowerRTBound() const =0
virtual void fit(FeatureFinderAlgorithmPickedHelperStructs::MassTraces &traces)=0
FeatureFinderAlgorithmPickedHelperStructs::MassTraces * traces_ptr
Definition: TraceFitter.h:176
bool weighted_
Whether to weight mass traces by theoretical intensity during the optimization.
Definition: TraceFitter.h:196
virtual double getHeight() const =0
void optimize_(Eigen::VectorXd &x_init, GenericFunctor &functor)
virtual String getGnuplotFormula(const FeatureFinderAlgorithmPickedHelperStructs::MassTrace &trace, const char function_name, const double baseline, const double rt_shift)=0
bool weighted
Definition: TraceFitter.h:177
virtual double getValue(double rt) const =0
virtual double getCenter() const =0
TraceFitter & operator=(const TraceFitter &source)
assignment operator
double computeTheoretical(const FeatureFinderAlgorithmPickedHelperStructs::MassTrace &trace, Size k) const
virtual double getUpperRTBound() const =0
virtual double getFWHM() const =0
TraceFitter(const TraceFitter &source)
copy constructor
void updateMembers_() override
This method is used to update extra member variables at the end of the setParameters() method.
virtual bool checkMaximalRTSpan(const double max_rt_span)=0
virtual double getArea()=0
SignedSize max_iterations_
Maximum number of iterations.
Definition: TraceFitter.h:194
~TraceFitter() override
destructor
virtual void getOptimizedParameters_(const Eigen::VectorXd &)=0
virtual bool checkMinimalRTSpan(const std::pair< double, double > &rt_bounds, const double min_rt_span)=0
Definition: TraceFitter.h:175
ptrdiff_t SignedSize
Signed Size type e.g. used as pointer difference.
Definition: Types.h:134
size_t Size
Size type e.g. used as variable which can hold result of size()
Definition: Types.h:127
Definition: IsobaricIsotopeCorrector.h:41
Matrix< double, -1, -1, 0, -1, -1 > MatrixXd
Definition: IsobaricIsotopeCorrector.h:44
Matrix< double, -1, 1, 0, -1, 1 > VectorXd
Definition: IsobaricIsotopeCorrector.h:45
Definition: IsobaricIsotopeCorrector.h:43
const double k
Definition: Constants.h:153
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:47
Helper struct for mass traces used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:80
Helper struct for a collection of mass traces used in FeatureFinderAlgorithmPicked.
Definition: FeatureFinderAlgorithmPickedHelperStructs.h:111