OpenMS
OptimizePick.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-2023.
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: Eva Lange $
33 // --------------------------------------------------------------------------
34 
35 #pragma once
36 
38 #include <OpenMS/KERNEL/Peak1D.h>
39 
40 // forward decl
41 namespace Eigen
42 {
43  template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
44  class Matrix;
45  using MatrixXd = Matrix<double, -1, -1, 0, -1, -1>;
46  using VectorXd = Matrix<double, -1, 1, 0, -1, 1>;
47 }
48 
49 #include <vector>
50 
51 namespace OpenMS
52 {
53  namespace OptimizationFunctions
54  {
56  typedef std::vector<Peak1D> RawDataVector;
58  typedef RawDataVector::iterator PeakIterator;
59 
68  struct OPENMS_DLLAPI PenaltyFactors
69  {
71  pos(0), lWidth(0), rWidth(0) {}
73  pos(p.pos), lWidth(p.lWidth), rWidth(p.rWidth) {}
75  {
76  pos = p.pos;
77  lWidth = p.lWidth;
78  rWidth = p.rWidth;
79 
80  return *this;
81  }
82 
84 
86  double pos;
88  double lWidth;
90  double rWidth;
91  };
92  }
93 
94 
101  class OPENMS_DLLAPI OptimizePick
102  {
103 public:
104 
105  struct Data
106  {
108  std::vector<double> positions;
109  std::vector<double> signal;
111  std::vector<PeakShape> peaks;
112 
114 
115  };
116 
118  {
119  public:
120  int inputs() const { return m_inputs; }
121  int values() const { return m_values; }
122 
123  OptPeakFunctor(unsigned dimensions, unsigned num_data_points, const OptimizePick::Data * data)
124  : m_inputs(dimensions), m_values(num_data_points), m_data(data) {}
125 
127  // compute Jacobian matrix for the different parameters
128  int df(const Eigen::VectorXd &x, Eigen::MatrixXd &J);
129 
130  private:
131  const int m_inputs, m_values;
132  const Data * m_data;
133  };
134 
136  typedef std::vector<Peak1D> RawDataVector;
138  typedef RawDataVector::iterator PeakIterator;
139 
140 
143  max_iteration_(400)
144  {}
145 
148  const int max_iteration_);
149 
152 
154  inline const struct OptimizationFunctions::PenaltyFactors & getPenalties() const { return penalties_; }
156  inline struct OptimizationFunctions::PenaltyFactors & getPenalties() { return penalties_; }
158  inline void setPenalties(const struct OptimizationFunctions::PenaltyFactors & penalties) { penalties_ = penalties; }
159 
161  inline UInt getNumberIterations() const { return max_iteration_; }
163  inline unsigned int & getNumberIterations() { return max_iteration_; }
165  inline void setNumberIterations(const int max_iteration) { max_iteration_ = max_iteration; }
166 
168  void optimize(std::vector<PeakShape> & peaks, Data & data);
169 
170 
171 protected:
173  struct OptimizationFunctions::PenaltyFactors penalties_;
174 
176  unsigned int max_iteration_;
177  };
178 }
179 
Definition: OptimizePick.h:118
OptPeakFunctor(unsigned dimensions, unsigned num_data_points, const OptimizePick::Data *data)
Definition: OptimizePick.h:123
int values() const
Definition: OptimizePick.h:121
int df(const Eigen::VectorXd &x, Eigen::MatrixXd &J)
const Data * m_data
Definition: OptimizePick.h:132
const int m_inputs
Definition: OptimizePick.h:131
int operator()(const Eigen::VectorXd &x, Eigen::VectorXd &fvec)
int inputs() const
Definition: OptimizePick.h:120
This class provides the non-linear optimization of the peak parameters.
Definition: OptimizePick.h:102
void setNumberIterations(const int max_iteration)
Mutable access to the number of iterations.
Definition: OptimizePick.h:165
void setPenalties(const struct OptimizationFunctions::PenaltyFactors &penalties)
Mutable access to the penalty factors.
Definition: OptimizePick.h:158
std::vector< Peak1D > RawDataVector
Profile data vector type.
Definition: OptimizePick.h:136
OptimizationFunctions::PenaltyFactors penalties
Definition: OptimizePick.h:113
OptimizePick()
Constructor.
Definition: OptimizePick.h:142
RawDataVector::iterator PeakIterator
Profile data iterator type.
Definition: OptimizePick.h:138
unsigned int max_iteration_
Maximum number of iterations during optimization.
Definition: OptimizePick.h:176
std::vector< double > positions
Positions and intensity values of the profile data.
Definition: OptimizePick.h:108
OptimizePick(const struct OptimizationFunctions::PenaltyFactors &penalties_, const int max_iteration_)
Constructor to set the penalty factors, the number of optimization iterations as well as the threshol...
unsigned int & getNumberIterations()
Mutable access to the number of iterations.
Definition: OptimizePick.h:163
struct OptimizationFunctions::PenaltyFactors & getPenalties()
Mutable access to the penalty factors.
Definition: OptimizePick.h:156
std::vector< PeakShape > peaks
This container contains the peak shapes to be optimized.
Definition: OptimizePick.h:111
const struct OptimizationFunctions::PenaltyFactors & getPenalties() const
Non-mutable access to the penalty factors.
Definition: OptimizePick.h:154
~OptimizePick()
Destructor.
void optimize(std::vector< PeakShape > &peaks, Data &data)
Start the optimization of the peak shapes peaks. The original peak shapes will be substituted by the ...
UInt getNumberIterations() const
Non-mutable access to the number of iterations.
Definition: OptimizePick.h:161
std::vector< double > signal
Definition: OptimizePick.h:109
Definition: OptimizePick.h:106
unsigned int UInt
Unsigned integer type.
Definition: Types.h:94
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
std::vector< Peak1D > RawDataVector
Profile data vector type.
Definition: OptimizePick.h:56
RawDataVector::iterator PeakIterator
Profile data iterator type.
Definition: OptimizePick.h:58
Main OpenMS namespace.
Definition: FeatureDeconvolution.h:48
Class for the penalty factors used during the optimization.
Definition: OptimizePick.h:69
double rWidth
Penalty factor for the peak shape's right width parameter.
Definition: OptimizePick.h:90
double pos
Penalty factor for the peak shape's position.
Definition: OptimizePick.h:86
PenaltyFactors()
Definition: OptimizePick.h:70
double lWidth
Penalty factor for the peak shape's left width parameter.
Definition: OptimizePick.h:88
PenaltyFactors & operator=(const PenaltyFactors &p)
Definition: OptimizePick.h:74
PenaltyFactors(const PenaltyFactors &p)
Definition: OptimizePick.h:72
~PenaltyFactors()
Definition: OptimizePick.h:83