Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages

OptimizePeakDeconvolution.h (Maintainer: Alexandra Zerck)

Go to the documentation of this file.
00001 // -*- Mode: C++; tab-width: 2; -*-
00002 // vi: set ts=2:
00003 //
00004 // --------------------------------------------------------------------------
00005 //                   OpenMS Mass Spectrometry Framework
00006 // --------------------------------------------------------------------------
00007 //  Copyright (C) 2003-2008 -- Oliver Kohlbacher, Knut Reinert
00008 //
00009 //  This library is free software; you can redistribute it and/or
00010 //  modify it under the terms of the GNU Lesser General Public
00011 //  License as published by the Free Software Foundation; either
00012 //  version 2.1 of the License, or (at your option) any later version.
00013 //
00014 //  This library is distributed in the hope that it will be useful,
00015 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
00016 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017 //  Lesser General Public License for more details.
00018 //
00019 //  You should have received a copy of the GNU Lesser General Public
00020 //  License along with this library; if not, write to the Free Software
00021 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 //
00023 // --------------------------------------------------------------------------
00024 // $Maintainer: Alexandra Zerck $
00025 // --------------------------------------------------------------------------
00026 //
00027 
00028 #ifndef OPENMS_TRANSFORMATIONS_RAW2PEAK_OPTIMIZEPEAKDECONVOLUTION_H
00029 #define OPENMS_TRANSFORMATIONS_RAW2PEAK_OPTIMIZEPEAKDECONVOLUTION_H
00030 
00031 #include <OpenMS/TRANSFORMATIONS/RAW2PEAK/PeakShape.h>
00032 #include <OpenMS/DATASTRUCTURES/Param.h>
00033 #include <gsl/gsl_vector.h>
00034 #include <gsl/gsl_multifit_nlin.h>
00035 #include <gsl/gsl_blas.h>
00036 #include <OpenMS/TRANSFORMATIONS/RAW2PEAK/OptimizePick.h>
00037 #include <OpenMS/DATASTRUCTURES/DefaultParamHandler.h>
00038 
00039 //#define DEBUG_DECONV
00040 #include <iostream>
00041 #ifdef DEBUG_DECONV
00042 #include <iostream>
00043 #include <fstream>
00044 #endif
00045 #include <vector>
00046 
00047 namespace OpenMS
00048 {
00049 
00050   namespace OptimizationFunctions
00051   {
00052     extern std::vector<PeakShape> peaks_DC_;
00053     extern std::vector<double> positions_DC_;
00054     extern std::vector<double> signal_DC_;
00055 
00062     struct PenaltyFactorsIntensity : public PenaltyFactors
00063     {
00064       PenaltyFactorsIntensity():PenaltyFactors(),height(0){}
00065       PenaltyFactorsIntensity(const PenaltyFactorsIntensity& p) : PenaltyFactors(p), height(p.height) {}
00066       inline PenaltyFactorsIntensity& operator=(const PenaltyFactorsIntensity& p)
00067       {
00068         height=p.height;
00069         pos=p.pos;
00070         lWidth=p.lWidth;
00071         rWidth=p.rWidth;
00072         
00073         return *this;
00074       }
00075       ~PenaltyFactorsIntensity(){}
00076 
00077       double height;
00078 
00079 
00080     };
00081     
00082 
00083     
00084   }//namespace OptimizationFunctions
00085 
00099   class OptimizePeakDeconvolution : public DefaultParamHandler
00100   {
00101   public:
00105     typedef std::vector<DRawDataPoint<1> > RawDataVector;
00106     typedef RawDataVector::iterator RawDataPointIterator;
00108 
00109 
00113 
00114     OptimizePeakDeconvolution( );
00115 
00117     OptimizePeakDeconvolution(const OptimizePeakDeconvolution& opt)
00118       :DefaultParamHandler(opt),
00119        penalties_(opt.penalties_),
00120        charge_(opt.charge_){}
00121 
00123     virtual ~OptimizePeakDeconvolution(){}
00125 
00129     inline OptimizePeakDeconvolution& operator=(const OptimizePeakDeconvolution& opt)
00130     {
00131       DefaultParamHandler::operator=(opt);
00132       penalties_=opt.penalties_;
00133       charge_=opt.charge_;
00134 
00135       return *this;
00136     }
00138 
00139 
00143 
00144     inline const OptimizationFunctions::PenaltyFactorsIntensity& getPenalties() const { return penalties_; }
00146     inline void setPenalties(const OptimizationFunctions::PenaltyFactorsIntensity& penalties)
00147     {
00148       penalties_ = penalties;
00149       param_.setValue("penalties:left_width",penalties_.lWidth);
00150       param_.setValue("penalties:right_width",penalties_.rWidth);
00151       param_.setValue("penalties:height",penalties_.height);
00152       param_.setValue("penalties:position",penalties_.pos);
00153     }
00154     
00156     inline const int getCharge() const { return charge_; }
00158     inline void setCharge(const int charge) { charge_ = charge; }
00160 
00161 
00163     bool optimize(std::vector<PeakShape>& peaks,int failure);
00164 
00165   protected:
00166     // Penalty factors for some paramter in the optimization
00167     OptimizationFunctions::PenaltyFactorsIntensity penalties_;
00168 
00170     int charge_;
00171 
00173     static const double dist_;
00174 
00176     int getNumberOfPeaks_(int charge, std::vector<PeakShape>& temp_shapes);
00177 
00178     // After each iteration the fwhm of all peaks is checked whether it isn't too large
00179     bool checkFWHM_(std::vector<PeakShape>& peaks,gsl_multifit_fdfsolver *& fit);
00180 
00181     void updateMembers_();
00182   };// class
00183   
00184 }// namespace OpenMS
00185 
00186 
00187 #endif

Generated Tue Apr 1 15:36:36 2008 -- using doxygen 1.5.4 OpenMS / TOPP 1.1