Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooChiSquarePdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * Kyle Cranmer
7 * *
8 *****************************************************************************/
9
10/** \class RooChiSquarePdf
11 \ingroup Roofit
12
13The PDF of the Chi Square distribution for n degrees of freedom.
14Oddly, this is hard to find in ROOT (except via relation to GammaDist).
15Here we also implement the analytic integral.
16**/
17
18#include "RooChiSquarePdf.h"
19#include "RooFit.h"
20#include "RooAbsReal.h"
21#include "RooRealVar.h"
22#include "RooBatchCompute.h"
23
24#include "TMath.h"
25
26#include <cmath>
27using namespace std;
28
30
31////////////////////////////////////////////////////////////////////////////////
32
34{
35}
36
37////////////////////////////////////////////////////////////////////////////////
38
39RooChiSquarePdf::RooChiSquarePdf(const char* name, const char* title,
40 RooAbsReal& x, RooAbsReal& ndof):
41 RooAbsPdf(name, title),
42 _x("x", "Dependent", this, x),
43 _ndof("ndof","ndof", this, ndof)
44{
45}
46
47////////////////////////////////////////////////////////////////////////////////
48
50 RooAbsPdf(other, name),
51 _x("x", this, other._x),
52 _ndof("ndof",this,other._ndof)
53{
54}
55
56////////////////////////////////////////////////////////////////////////////////
57
59{
60 if(_x <= 0) return 0;
61
62 return pow(_x,(_ndof/2.)-1.) * exp(-_x/2.) / TMath::Gamma(_ndof/2.) / pow(2.,_ndof/2.);
63}
64
65////////////////////////////////////////////////////////////////////////////////
66/// Compute multiple values of ChiSquare distribution.
68 return RooBatchCompute::dispatch->computeChiSquare(this, evalData, _x->getValues(evalData, normSet), _ndof->getValues(evalData, normSet));
69}
70
71////////////////////////////////////////////////////////////////////////////////
72/// No analytical calculation available (yet) of integrals over subranges
73
74Int_t RooChiSquarePdf::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* rangeName) const
75{
76 if (rangeName && strlen(rangeName)) {
77 return 0 ;
78 }
79
80 if (matchArgs(allVars, analVars, _x)) return 1;
81 return 0;
82}
83
84////////////////////////////////////////////////////////////////////////////////
85
86Double_t RooChiSquarePdf::analyticalIntegral(Int_t code, const char* rangeName) const
87{
88 assert(1 == code); (void)code;
89 Double_t xmin = _x.min(rangeName); Double_t xmax = _x.max(rangeName);
90
91 // TMath::Prob needs ndof to be an integer, or it returns 0.
92 // return TMath::Prob(xmin, _ndof) - TMath::Prob(xmax,_ndof);
93
94 // cumulative is known based on lower incomplete gamma function, or regularized gamma function
95 // Wikipedia defines lower incomplete gamma function without the normalization 1/Gamma(ndof),
96 // but it is included in the ROOT implementation.
97 Double_t pmin = TMath::Gamma(_ndof/2,xmin/2);
98 Double_t pmax = TMath::Gamma(_ndof/2,xmax/2);
99
100 // only use this if range is appropriate
101 return pmax-pmin;
102}
#define ClassImp(name)
Definition Rtypes.h:364
char name[80]
Definition TGX11.cxx:110
float xmin
float xmax
double pow(double, double)
double exp(double)
typedef void((*Func_t)())
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:61
virtual RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet=nullptr) const
by this change, please consult the release notes for ROOT 6.24 for guidance on how to make this trans...
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:29
virtual RooSpan< double > computeChiSquare(const RooAbsReal *, RunContext &, RooSpan< const double > x, RooSpan< const double > ndof)=0
The PDF of the Chi Square distribution for n degrees of freedom.
Double_t evaluate() const
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
RooRealProxy _ndof
RooSpan< double > evaluateSpan(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Compute multiple values of ChiSquare distribution.
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
No analytical calculation available (yet) of integrals over subranges.
RooRealProxy _x
A simple container to hold a batch of data values.
Definition RooSpan.h:34
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
Double_t x[n]
Definition legend1.C:17
R__EXTERN RooBatchComputeInterface * dispatch
This dispatch pointer points to an implementation of the compute library, provided one has been loade...
Double_t Gamma(Double_t z)
Computation of gamma(z) for all z.
Definition TMath.cxx:348
This struct enables passing computation data around between elements of a computation graph.
Definition RunContext.h:31