00001 #ifndef __BCMODELHISTOGRAMFITTER__H 00002 #define __BCMODELHISTOGRAMFITTER__H 00003 00004 /*! 00005 * \class BCHistogramFitter 00006 * \brief A class for fitting histograms with functions 00007 * \author Daniel Kollar 00008 * \author Kevin Kröninger 00009 * \version 1.0 00010 * \date 11.2008 00011 * \detail This class allows fitting of a TH1D histogram using 00012 * a TF1 function. 00013 */ 00014 00015 /* 00016 * Copyright (C) 2008, Daniel Kollar and Kevin Kroeninger. 00017 * All rights reserved. 00018 * 00019 * For the licensing terms see doc/COPYING. 00020 */ 00021 00022 // --------------------------------------------------------- 00023 00024 #include <vector> 00025 00026 #include <BAT/BCModel.h> 00027 00028 // ROOT classes 00029 class TH1D; 00030 class TF1; 00031 00032 // --------------------------------------------------------- 00033 00034 class BCHistogramFitter : public BCModel 00035 { 00036 public: 00037 00038 /** \name Constructors and destructors */ 00039 /* @{ */ 00040 00041 /** 00042 * The default constructor. */ 00043 BCHistogramFitter(); 00044 00045 /** 00046 * A constructor. 00047 * @param hist The histogram (TH1D). 00048 * @param func The fit function. */ 00049 BCHistogramFitter(TH1D * hist, TF1 * func); 00050 00051 /** 00052 * The default destructor. */ 00053 ~BCHistogramFitter(); 00054 00055 /* @} */ 00056 00057 /** \name Member functions (get) */ 00058 /* @{ */ 00059 00060 /** 00061 * @return The data histogram */ 00062 TH1D * GetHistogram() 00063 { return fHistogram; }; 00064 00065 /** 00066 * @return The histogram of expected counts*/ 00067 TH1D * GetHistogramExpected() 00068 { return fHistogramExpected; }; 00069 00070 00071 /** 00072 * @return The fit function */ 00073 TF1 * GetFitFunction() 00074 { return fFitFunction; }; 00075 00076 /** 00077 * @return pointer to the error band */ 00078 TGraph * GetErrorBand() 00079 { return fErrorBand; }; 00080 00081 /** 00082 * @return pointer to a graph for the fit function */ 00083 TGraph * GetGraphFitFunction() 00084 { return fGraphFitFunction; }; 00085 00086 /* @} */ 00087 00088 /** \name Member functions (set) */ 00089 /* @{ */ 00090 00091 /** 00092 * @param hist The histogram containing the data 00093 * @ return An error code (1:pass, 0:fail). 00094 */ 00095 int SetHistogram(TH1D * hist); 00096 00097 /** 00098 * @param hist The histogram with the expected counts (typically non-integer values!) 00099 * @ return An error code (1:pass, 0:fail). 00100 */ 00101 int SetHistogramExpected(const std::vector <double>& parameters); 00102 00103 /** 00104 * @param func The fit function 00105 * @ return An error code (1:pass, 0:fail). 00106 */ 00107 int SetFitFunction(TF1 * func); 00108 00109 /** 00110 * Sets the flag for integration. \n 00111 * true: use ROOT's TH1D::Integrate() \n 00112 * false: use linear interpolation */ 00113 void SetFlagIntegration(bool flag) 00114 { fFlagIntegration = flag; }; 00115 00116 /* @} */ 00117 /** \name Member functions (miscellaneous methods) */ 00118 /* @{ */ 00119 00120 /** 00121 * The log of the prior probability. Overloaded from BCModel. 00122 * @param parameters A vector of doubles containing the parameter values. */ 00123 // virtual double LogAPrioriProbability(std::vector <double> parameters); 00124 00125 /** 00126 * The log of the conditional probability. Overloaded from BCModel. 00127 * @param parameters A vector of doubles containing the parameter values. */ 00128 virtual double LogLikelihood(std::vector <double> parameters); 00129 00130 /** 00131 * Plots the histogram 00132 * @param options Options for plotting. 00133 * @param filename Name of the file which the histogram is printed into. 00134 * The following options are available:\n 00135 * F : plots the fit function on top of the data 00136 * E0 : plots the fit function and the 68% prob. uncertainty band of the fit function on top of the data 00137 * E1 : plots the expectation from the fit function and the uncertainty bin-by-bin as error bars. */ 00138 // void PrintHistogram(const char * options = "", const char * filename = ""); 00139 00140 /** 00141 * Returns the y-value of the 1-dimensional fit function at an x and 00142 * for a set of parameters. 00143 * @param x A vector with the x-value. 00144 * @param parameters A set of parameters. */ 00145 double FitFunction(std::vector <double> x, std::vector <double> parameters); 00146 00147 /** 00148 * Performs the fit. 00149 * @return An error code. */ 00150 int Fit(); 00151 00152 /** 00153 * Performs the fit. 00154 * @param hist The histogram (TH1D). 00155 * @param func The fit function. 00156 * @return An error code. */ 00157 int Fit(TH1D * hist, TF1 * func); 00158 00159 /** 00160 * Draw the fit in the current pad. */ 00161 void DrawFit(const char * options = "", bool flaglegend = false); 00162 00163 /** 00164 * Calculate the p-value using fast-MCMC. 00165 * @param par A set of parameter values 00166 * @param pvalue The pvalue 00167 * @param nIterations number of pseudo experiments generated by the Markov chain 00168 * @return An error code */ 00169 int CalculatePValueFast(std::vector<double> par, double &pvalue, int nIterations = 100000); 00170 00171 /** 00172 * Calculate the p-value using approximate chi^2 distribution of scaled likelihood. 00173 * Approximation is valid for bin contents >5, see eq. (32.12), 00174 * PDG: Statistics, Monte Carlo, Group Theory. Physics Letters B 667, 316-339(2008). 00175 * @param par The set of parameter values used in the model, usually the best fit parameters 00176 * @param pvalue The pvalue 00177 * @return An error code */ 00178 int CalculatePValueLikelihood(std::vector<double> par, double &pvalue); 00179 00180 /** 00181 * Calculate the p-value using approximate chi^2 distribution of squared difference 00182 * for conventional weights. 00183 * Approximation is valid for bin contents >5 and not as as good for little data as 00184 * CalculatePValueLikelihood, see eq. (32.13), 00185 * PDG: Statistics, Monte Carlo, Group Theory. Physics Letters B 667, 316-339(2008). 00186 * @param par The set of parameter values used in the model, usually the best fit parameters 00187 * @param pvalue The pvalue 00188 * @param weight use the variance from the expected #counts (true) or the measured counts (false) 00189 * @return An error code */ 00190 int CalculatePValueLeastSquares(std::vector<double> par, double &pvalue, bool weightExpect=true); 00191 00192 /** 00193 * Calculate the p-value using Kolmogorov-Smirnov test. Note that the 00194 * reference distribution is known only asymptotically. Some explanation is given in 00195 * http://root.cern.ch/root/htmldoc/TMath.html 00196 * @param par The set of parameter values used in the model, usually the best fit parameters 00197 * @param pvalue The pvalue 00198 * @return An error code */ 00199 int CalculatePValueKolmogorov(std::vector<double> par, double &pvalue); 00200 00201 00202 double CDF(const std::vector<double>& parameters, int index, bool lower=false); 00203 00204 /* @} */ 00205 00206 private: 00207 00208 /** 00209 * The histogram containing the data. 00210 */ 00211 TH1D * fHistogram; 00212 00213 /** 00214 * The fit function */ 00215 TF1 * fFitFunction; 00216 00217 /** 00218 * Flag for using the ROOT TH1D::Integral method (true), or linear 00219 * interpolation (false) */ 00220 bool fFlagIntegration; 00221 00222 /** 00223 * Pointer to the error band (for legend) */ 00224 TGraph * fErrorBand; 00225 00226 /** 00227 * Pointer to a graph for displaying the fit function */ 00228 TGraph * fGraphFitFunction; 00229 00230 /** 00231 * The histogram containing the expected data. 00232 */ 00233 TH1D * fHistogramExpected; 00234 }; 00235 00236 // --------------------------------------------------------- 00237 00238 #endif