OfflineValidation/interface/FitWithRooFit.h

0001 #ifndef Alignment_OfflineValidation_FitWithRooFit_h
0002 #define Alignment_OfflineValidation_FitWithRooFit_h
0003
0004 #include "RooGlobalFunc.h"
0005 #include "TCanvas.h"
0006 #include "TTree.h"
0007 #include "TH1D.h"
0008 #include "TRandom.h"
0009 #include "RooRealVar.h"
0010 #include "RooDataSet.h"
0011 #include "RooGaussian.h"
0012 #include "RooVoigtian.h"
0013 #include "RooExponential.h"
0014 #include "RooPlot.h"
0015 #include "RooDataHist.h"
0016 #include "RooAddPdf.h"
0017 #include "RooChebychev.h"
0018 #include "RooFitResult.h"
0019 #include "RooGenericPdf.h"
0020 #include "RooGaussModel.h"
0021 #include "RooAddModel.h"
0022 #include "RooPolynomial.h"
0023 #include "RooCBShape.h"
0024 #include "RooMinimizer.h"
0025 #include "RooBreitWigner.h"
0026 #include "RooFFTConvPdf.h"
0027
0028 /**
0029  * This class allows to use RooFit to perform a fit on a TH1 histogram. <br>
0030  * The currently implemented functions are: <br>
0031  * - signal: <br>
0032  * -- gaussian <br>
0033  * -- double gaussian <br>
0034  * -- voigtian <br>
0035  * - background: <br>
0036  * -- exponential <br>
0037  * It is possible to select any combination of signal and background. <br>
0038  * The fit() method receives the TH1, two strings specifying the signal and background type
0039  * and the min and max x of the histogram over which to perform the fit. <br>
0040  * The variables of the fit must be initialized separately. For example when doing a gaussian+exponential
0041  * fit the initMean, initSigma, initFSig and initExpCoeff methods must be used to create and initialize
0042  * the corresponding variables. <br>
0043  * The methods names after the variables return the fit result.
0044  */
0045
0046 class FitWithRooFit {
0047 public:
0048   std::unique_ptr<RooDataHist> importTH1(TH1* histo, double xMin, double xMax);
0049
0050   void fit(TH1* histo,
0051            const TString signalType,
0052            const TString backgroundType,
0053            double xMin = 0.,
0054            double xMax = 0.,
0055            bool sumW2Error = false);
0056
0057   void initMean(double value, double min, double max, const TString& name = "mean", const TString& title = "mean");
0058   void initMean2(double value, double min, double max, const TString& name = "mean2", const TString& title = "mean2");
0059   void initMean3(double value, double min, double max, const TString& name = "mean3", const TString& title = "mean3");
0060   void initSigma(double value, double min, double max, const TString& name = "sigma", const TString& title = "sigma");
0061   void initSigma2(double value, double min, double max, const TString& name = "sigma2", const TString& title = "sigma2");
0062   void initSigma3(double value, double min, double max, const TString& name = "sigma3", const TString& title = "sigma3");
0063   void initGamma(double value, double min, double max, const TString& name = "gamma", const TString& title = "gamma");
0064   void initGaussFrac(
0065       double value, double min, double max, const TString& name = "GaussFrac", const TString& title = "GaussFrac");
0066   void initGaussFrac2(
0067       double value, double min, double max, const TString& name = "GaussFrac2", const TString& title = "GaussFrac2");
0068   void initExpCoeffA0(
0069       double value, double min, double max, const TString& name = "expCoeffa0", const TString& title = "expCoeffa0");
0070   void initExpCoeffA1(
0071       double value, double min, double max, const TString& name = "expCoeffa1", const TString& title = "expCoeffa1");
0072   void initExpCoeffA2(
0073       double value, double min, double max, const TString& name = "expCoeffa2", const TString& title = "expCoeffa2");
0074   void initFsig(
0075       double value, double min, double max, const TString& name = "fsig", const TString& title = "signal fraction");
0076   void initA0(double value, double min, double max, const TString& name = "a0", const TString& title = "a0");
0077   void initA1(double value, double min, double max, const TString& name = "a1", const TString& title = "a1");
0078   void initA2(double value, double min, double max, const TString& name = "a2", const TString& title = "a2");
0079   void initA3(double value, double min, double max, const TString& name = "a3", const TString& title = "a3");
0080   void initA4(double value, double min, double max, const TString& name = "a4", const TString& title = "a4");
0081   void initA5(double value, double min, double max, const TString& name = "a5", const TString& title = "a5");
0082   void initA6(double value, double min, double max, const TString& name = "a6", const TString& title = "a6");
0083   void initAlpha(double value, double min, double max, const TString& name = "alpha", const TString& title = "alpha");
0084   void initN(double value, double min, double max, const TString& name = "n", const TString& title = "n");
0085   void initFGCB(double value, double min, double max, const TString& name = "fGCB", const TString& title = "fGCB");
0086
0087   inline RooRealVar* mean() { return mean_.get(); }
0088   inline RooRealVar* mean2() { return mean2_.get(); }
0089   inline RooRealVar* mean3() { return mean3_.get(); }
0090   inline RooRealVar* sigma() { return sigma_.get(); }
0091   inline RooRealVar* sigma2() { return sigma2_.get(); }
0092   inline RooRealVar* sigma3() { return sigma3_.get(); }
0093   inline RooRealVar* gamma() { return gamma_.get(); }
0094   inline RooRealVar* gaussFrac() { return gaussFrac_.get(); }
0095   inline RooRealVar* gaussFrac2() { return gaussFrac2_.get(); }
0096   inline RooRealVar* expCoeffa0() { return expCoeffa0_.get(); }
0097   inline RooRealVar* expCoeffa1() { return expCoeffa1_.get(); }
0098   inline RooRealVar* expCoeffa2() { return expCoeffa2_.get(); }
0099   inline RooRealVar* fsig() { return fsig_.get(); }
0100   inline RooRealVar* a0() { return a0_.get(); }
0101   inline RooRealVar* a1() { return a1_.get(); }
0102   inline RooRealVar* a2() { return a2_.get(); }
0103   inline RooRealVar* a3() { return a3_.get(); }
0104   inline RooRealVar* a4() { return a4_.get(); }
0105   inline RooRealVar* a5() { return a5_.get(); }
0106   inline RooRealVar* a6() { return a6_.get(); }
0107   inline RooRealVar* alpha() { return alpha_.get(); }
0108   inline RooRealVar* n() { return n_.get(); }
0109   inline RooRealVar* fGCB() { return fGCB_.get(); }
0110
0111   void reinitializeParameters();
0112
0113   std::unique_ptr<RooAbsPdf> buildSignalModel(RooRealVar* x, const TString& signalType);
0114   std::unique_ptr<RooAbsPdf> buildBackgroundModel(RooRealVar* x, const TString& backgroundType);
0115   std::unique_ptr<RooAbsPdf> buildModel(RooRealVar* x, const TString& signalType, const TString& backgroundType);
0116
0117   bool useChi2_ = false;
0118
0119 protected:
0120   // Declare all variables
0121   std::unique_ptr<RooRealVar> mean_;
0122   std::unique_ptr<RooRealVar> mean2_;
0123   std::unique_ptr<RooRealVar> mean3_;
0124   std::unique_ptr<RooRealVar> sigma_;
0125   std::unique_ptr<RooRealVar> sigma2_;
0126   std::unique_ptr<RooRealVar> sigma3_;
0127   std::unique_ptr<RooRealVar> gamma_;
0128   std::unique_ptr<RooRealVar> gaussFrac_;
0129   std::unique_ptr<RooRealVar> gaussFrac2_;
0130   std::unique_ptr<RooRealVar> expCoeffa0_;
0131   std::unique_ptr<RooRealVar> expCoeffa1_;
0132   std::unique_ptr<RooRealVar> expCoeffa2_;
0133   std::unique_ptr<RooRealVar> fsig_;
0134   std::unique_ptr<RooRealVar> a0_;
0135   std::unique_ptr<RooRealVar> a1_;
0136   std::unique_ptr<RooRealVar> a2_;
0137   std::unique_ptr<RooRealVar> a3_;
0138   std::unique_ptr<RooRealVar> a4_;
0139   std::unique_ptr<RooRealVar> a5_;
0140   std::unique_ptr<RooRealVar> a6_;
0141   std::unique_ptr<RooRealVar> alpha_;
0142   std::unique_ptr<RooRealVar> n_;
0143   std::unique_ptr<RooRealVar> fGCB_;
0144
0145   // Initial values
0146   double initVal_mean = 0.0;
0147   double initVal_mean2 = 0.0;
0148   double initVal_mean3 = 0.0;
0149   double initVal_sigma = 0.0;
0150   double initVal_sigma2 = 0.0;
0151   double initVal_sigma3 = 0.0;
0152   double initVal_gamma = 0.0;
0153   double initVal_gaussFrac = 0.0;
0154   double initVal_gaussFrac2 = 0.0;
0155   double initVal_expCoeffa0 = 0.0;
0156   double initVal_expCoeffa1 = 0.0;
0157   double initVal_expCoeffa2 = 0.0;
0158   double initVal_fsig = 0.0;
0159   double initVal_a0 = 0.0;
0160   double initVal_a1 = 0.0;
0161   double initVal_a2 = 0.0;
0162   double initVal_a3 = 0.0;
0163   double initVal_a4 = 0.0;
0164   double initVal_a5 = 0.0;
0165   double initVal_a6 = 0.0;
0166   double initVal_alpha = 0.0;
0167   double initVal_n = 0.0;
0168   double initVal_fGCB = 0.0;
0169 };
0170
0171 #endif