00001 ////////////////////////////////////////////////////////////////////////// 00002 // 00003 // 'BASIC FUNCTIONALITY' RooFit tutorial macro #103 00004 // 00005 // Interpreted functions and p.d.f.s 00006 // 00007 // 00008 // 00009 // 07/2008 - Wouter Verkerke 00010 // 00011 ///////////////////////////////////////////////////////////////////////// 00012 00013 #ifndef __CINT__ 00014 #include "RooGlobalFunc.h" 00015 #endif 00016 #include "RooRealVar.h" 00017 #include "RooDataSet.h" 00018 #include "RooGaussian.h" 00019 #include "TCanvas.h" 00020 #include "RooPlot.h" 00021 #include "RooFitResult.h" 00022 #include "RooGenericPdf.h" 00023 00024 using namespace RooFit ; 00025 00026 00027 class TestBasic103 : public RooFitTestUnit 00028 { 00029 public: 00030 TestBasic103(TFile* refFile, Bool_t writeRef, Int_t verbose) : RooFitTestUnit("Interpreted expression p.d.f.",refFile,writeRef,verbose) {} ; 00031 Bool_t testCode() { 00032 00033 ///////////////////////////////////////////////////////// 00034 // G e n e r i c i n t e r p r e t e d p . d . f . // 00035 ///////////////////////////////////////////////////////// 00036 00037 // Declare observable x 00038 RooRealVar x("x","x",-20,20) ; 00039 00040 // C o n s t r u c t g e n e r i c p d f f r o m i n t e r p r e t e d e x p r e s s i o n 00041 // ------------------------------------------------------------------------------------------------- 00042 00043 // To construct a proper p.d.f, the formula expression is explicitly normalized internally by dividing 00044 // it by a numeric integral of the expresssion over x in the range [-20,20] 00045 // 00046 RooRealVar alpha("alpha","alpha",5,0.1,10) ; 00047 RooGenericPdf genpdf("genpdf","genpdf","(1+0.1*abs(x)+sin(sqrt(abs(x*alpha+0.1))))",RooArgSet(x,alpha)) ; 00048 00049 00050 // S a m p l e , f i t a n d p l o t g e n e r i c p d f 00051 // --------------------------------------------------------------- 00052 00053 // Generate a toy dataset from the interpreted p.d.f 00054 RooDataSet* data = genpdf.generate(x,10000) ; 00055 00056 // Fit the interpreted p.d.f to the generated data 00057 genpdf.fitTo(*data) ; 00058 00059 // Make a plot of the data and the p.d.f overlaid 00060 RooPlot* xframe = x.frame(Title("Interpreted expression pdf")) ; 00061 data->plotOn(xframe) ; 00062 genpdf.plotOn(xframe) ; 00063 00064 00065 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00066 // S t a n d a r d p . d . f a d j u s t w i t h i n t e r p r e t e d h e l p e r f u n c t i o n // 00067 // // 00068 // Make a gauss(x,sqrt(mean2),sigma) from a standard RooGaussian // 00069 // // 00070 ///////////////////////////////////////////////////////////////////////////////////////////////////////////////// 00071 00072 00073 // C o n s t r u c t s t a n d a r d p d f w i t h f o r m u l a r e p l a c i n g p a r a m e t e r 00074 // ------------------------------------------------------------------------------------------------------------ 00075 00076 // Construct parameter mean2 and sigma 00077 RooRealVar mean2("mean2","mean^2",10,0,200) ; 00078 RooRealVar sigma("sigma","sigma",3,0.1,10) ; 00079 00080 // Construct interpreted function mean = sqrt(mean^2) 00081 RooFormulaVar mean("mean","mean","sqrt(mean2)",mean2) ; 00082 00083 // Construct a gaussian g2(x,sqrt(mean2),sigma) ; 00084 RooGaussian g2("g2","h2",x,mean,sigma) ; 00085 00086 00087 // G e n e r a t e t o y d a t a 00088 // --------------------------------- 00089 00090 // Construct a separate gaussian g1(x,10,3) to generate a toy Gaussian dataset with mean 10 and width 3 00091 RooGaussian g1("g1","g1",x,RooConst(10),RooConst(3)) ; 00092 RooDataSet* data2 = g1.generate(x,1000) ; 00093 00094 00095 // F i t a n d p l o t t a i l o r e d s t a n d a r d p d f 00096 // ------------------------------------------------------------------- 00097 00098 // Fit g2 to data from g1 00099 RooFitResult* r = g2.fitTo(*data2,Save()) ; 00100 00101 // Plot data on frame and overlay projection of g2 00102 RooPlot* xframe2 = x.frame(Title("Tailored Gaussian pdf")) ; 00103 data2->plotOn(xframe2) ; 00104 g2.plotOn(xframe2) ; 00105 00106 regPlot(xframe,"rf103_plot1") ; 00107 regPlot(xframe2,"rf103_plot2") ; 00108 regResult(r,"rf103_fit1") ; 00109 00110 delete data ; 00111 delete data2 ; 00112 00113 return kTRUE ; 00114 } 00115 } ;