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 { return this -> Fit(fHistogram, fFitFunction); }; 00152 00153 /** 00154 * Performs the fit. 00155 * @param hist The histogram (TH1D). 00156 * @param func The fit function. 00157 * @return An error code. */ 00158 int Fit(TH1D * hist, TF1 * func); 00159 00160 /** 00161 * Draw the fit in the current pad. */ 00162 void DrawFit(const char * options = "", bool flaglegend = false); 00163 00164 /** 00165 * Calculate the p-value using fast-MCMC. 00166 * @param par A set of parameter values 00167 * @param pvalue The pvalue 00168 * @param nIterations number of pseudo experiments generated by the Markov chain 00169 * @return An error code */ 00170 int CalculatePValueFast(std::vector<double> par, double &pvalue, int nIterations = 100000); 00171 00172 /** 00173 * Calculate the p-value using approximate chi^2 distribution of scaled likelihood. 00174 * Approximation is valid for bin contents >5, see eq. (32.12), 00175 * PDG: Statistics, Monte Carlo, Group Theory. Physics Letters B 667, 316-339(2008). 00176 * @param par The set of parameter values used in the model, usually the best fit parameters 00177 * @param pvalue The pvalue 00178 * @return An error code */ 00179 int CalculatePValueLikelihood(std::vector<double> par, double &pvalue); 00180 00181 /** 00182 * Calculate the p-value using approximate chi^2 distribution of squared difference 00183 * for conventional weights. 00184 * Approximation is valid for bin contents >5 and not as as good for little data as 00185 * CalculatePValueLikelihood, see eq. (32.13), 00186 * PDG: Statistics, Monte Carlo, Group Theory. Physics Letters B 667, 316-339(2008). 00187 * @param par The set of parameter values used in the model, usually the best fit parameters 00188 * @param pvalue The pvalue 00189 * @param weight use the variance from the expected #counts (true) or the measured counts (false) 00190 * @return An error code */ 00191 int CalculatePValueLeastSquares(std::vector<double> par, double &pvalue, bool weightExpect=true); 00192 00193 /** 00194 * Calculate the p-value using Kolmogorov-Smirnov test. Note that the 00195 * reference distribution is known only asymptotically. Some explanation is given in 00196 * http://root.cern.ch/root/htmldoc/TMath.html 00197 * @param par The set of parameter values used in the model, usually the best fit parameters 00198 * @param pvalue The pvalue 00199 * @return An error code */ 00200 int CalculatePValueKolmogorov(std::vector<double> par, double &pvalue); 00201 00202 00203 double CDF(const std::vector<double>& parameters, int index, bool lower=false); 00204 00205 /* @} */ 00206 00207 private: 00208 00209 /** 00210 * The histogram containing the data. 00211 */ 00212 TH1D * fHistogram; 00213 00214 /** 00215 * The fit function */ 00216 TF1 * fFitFunction; 00217 00218 /** 00219 * Flag for using the ROOT TH1D::Integral method (true), or linear 00220 * interpolation (false) */ 00221 bool fFlagIntegration; 00222 00223 /** 00224 * Pointer to the error band (for legend) */ 00225 TGraph * fErrorBand; 00226 00227 /** 00228 * Pointer to a graph for displaying the fit function */ 00229 TGraph * fGraphFitFunction; 00230 00231 /** 00232 * The histogram containing the expected data. 00233 */ 00234 TH1D * fHistogramExpected; 00235 }; 00236 00237 // --------------------------------------------------------- 00238 00239 #endif