Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooExtendPdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
5 * Authors: *
6 * WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu *
7 * DK, David Kirkby, UC Irvine, dkirkby@uci.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * and Stanford University. All rights reserved. *
11 * *
12 * Redistribution and use in source and binary forms, *
13 * with or without modification, are permitted according to the terms *
14 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
15 *****************************************************************************/
16
17/** \class RooExtendPdf
18RooExtendPdf is a wrapper around an existing PDF that adds a
19parameteric extended likelihood term to the PDF, optionally divided by a
20fractional term from a partial normalization of the PDF:
21\f[
22 n_\mathrm{Expected} = N \quad \text{or} \quad n_\mathrm{Expected} = N / \mathrm{frac}
23\f]
24where \f$ N \f$ is supplied as a RooAbsReal to RooExtendPdf.
25The fractional term is defined as
26\f[
27 \mathrm{frac} = \frac{\int_\mathrm{cutRegion[x]} \mathrm{pdf}(x,y) \; \mathrm{d}x \mathrm{d}y}{
28 \int_\mathrm{normRegion[x]} \mathrm{pdf}(x,y) \; \mathrm{d}x \mathrm{d}y}
29\f]
30
31where \f$ x \f$ is the set of dependents involved in the selection region and \f$ y \f$
32is the set of remaining dependents.
33
34\f$ \mathrm{cutRegion}[x] \f$ is a limited integration range that is contained in
35the nominal integration range \f$ \mathrm{normRegion}[x] \f$.
36*/
37
38#include "Riostream.h"
39
40#include "RooExtendPdf.h"
41#include "RooArgList.h"
42#include "RooRealVar.h"
43#include "RooFormulaVar.h"
44#include "RooNameReg.h"
45#include "RooConstVar.h"
46#include "RooProduct.h"
47#include "RooRatio.h"
48#include "RooMsgService.h"
49
50
51
52using std::endl;
53
55
56RooExtendPdf::RooExtendPdf(const char *name, const char *title, RooAbsPdf& pdf,
57 RooAbsReal& norm, const char* rangeName)
58 : RooExtendPdf{name, title, pdf, RooAbsReal::Ref{norm}, rangeName} {}
59
60/// Constructor. The ExtendPdf behaves identical to the supplied input pdf,
61/// but adds an extended likelihood term. expectedEvents() will return
62/// `norm` if `rangeName` remains empty. If `rangeName` is not empty,
63/// `norm` will refer to this range, and expectedEvents will return the
64/// total number of events over the full range of the observables.
65/// \param[in] name Name of the pdf
66/// \param[in] title Title of the pdf (for plotting)
67/// \param[in] pdf The pdf to be extended
68/// \param[in] norm Expected number of events
69/// \param[in] rangeName If given, the number of events denoted by `norm` is interpreted as
70/// the number of events in this range only
71RooExtendPdf::RooExtendPdf(const char *name, const char *title, RooAbsPdf& pdf,
72 RooAbsReal::Ref norm, const char* rangeName) :
73 RooAbsPdf(name,title),
74 _pdf("pdf", "PDF", this, pdf),
75 _n("n","Normalization",this,norm),
76 _rangeName(RooNameReg::ptr(rangeName))
77{
78
79 // Copy various setting from pdf
80 setUnit(_pdf->getUnit()) ;
82}
83
84
85
86RooExtendPdf::RooExtendPdf(const RooExtendPdf& other, const char* name) :
87 RooAbsPdf(other,name),
88 _pdf("pdf",this,other._pdf),
89 _n("n",this,other._n),
91{
92 // Copy constructor
93}
94
95
96/// Return the number of expected events over the full range of all variables.
97/// `norm`, the variable set as normalisation constant in the constructor,
98/// will yield the number of events in the range set in the constructor. That is, the function returns
99/// \f[
100/// N = \mathrm{norm} \; \cdot \; \frac{\int_{(x_F,y_F)} \mathrm{pdf}(x,y) }{\int_{(x_C,y_F)} \mathrm{pdf}(x,y)}
101/// \f]
102/// Where \f$ x \f$ is the set of dependents with a restricted range (defined by `rangeName` in the constructor),
103/// and \f$ y \f$ are the other dependents. \f$ x_C \f$ is the integration
104/// of \f$ x \f$ over the restricted range, and \f$ x_F \f$ is the integration of
105/// \f$ x \f$ over the full range. `norm` is the number of events given as parameter to the constructor.
106///
107/// If the nested PDF can be extended, \f$ N \f$ is further scaled by its expected number of events.
109{
110 const RooAbsPdf& pdf = *_pdf;
111
112 if (_rangeName && (!nset || nset->empty())) {
113 coutW(InputArguments) << "RooExtendPdf::expectedEvents(" << GetName() << ") WARNING: RooExtendPdf needs non-null normalization set to calculate fraction in range "
114 << _rangeName << ". Results may be nonsensical" << endl ;
115 }
116
117 double nExp = _n ;
118
119 // Optionally multiply with fractional normalization
120 if (_rangeName) {
121
122 double fracInt = pdf.getNormObj(nset,nset,_rangeName)->getVal();
123
124
125 if ( fracInt == 0. || _n == 0.) {
126 coutW(Eval) << "RooExtendPdf(" << GetName() << ") WARNING: nExpected = " << _n << " / "
127 << fracInt << " for nset = " << (nset?*nset:RooArgSet()) << endl ;
128 }
129
130 nExp /= fracInt ;
131 }
132
133 // Multiply with original Nexpected, if defined
134 if (pdf.canBeExtended()) nExp *= pdf.expectedEvents(nset) ;
135
136 return nExp ;
137}
138
139
140std::unique_ptr<RooAbsReal> RooExtendPdf::createExpectedEventsFunc(const RooArgSet *nset) const
141{
142 const RooAbsPdf& pdf = *_pdf;
143
144 RooArgList prodList;
145 prodList.add(*_n);
146
147 // Optionally multiply with fractional normalization
148 std::unique_ptr<RooAbsReal> rangeFactor;
149 if (_rangeName) {
150 std::unique_ptr<RooAbsReal> fracInteg{pdf.createIntegral(*nset, *nset, RooNameReg::str(_rangeName))};
151 // Create one over integral term
152 auto rangeFactorName = std::string("one_over_") + fracInteg->GetName();
153 rangeFactor = std::make_unique<RooRatio>(rangeFactorName.c_str(), rangeFactorName.c_str(), RooFit::RooConst(1.0), *fracInteg);
154 rangeFactor->addOwnedComponents(std::move(fracInteg));
155 prodList.add(*rangeFactor);
156 }
157
158 // Multiply with original Nexpected, if defined
159 std::unique_ptr<RooAbsReal> pdfExpectedEvents;
160 if (pdf.canBeExtended()) {
161 pdfExpectedEvents = pdf.createExpectedEventsFunc(nset);
162 prodList.add(*pdfExpectedEvents);
163 }
164
165 auto name = std::string(GetName()) + "_expectedEvents";
166 auto out = std::make_unique<RooProduct>(name.c_str(), name.c_str(), prodList);
167 if(rangeFactor) {
168 out->addOwnedComponents(std::move(rangeFactor));
169 }
170 if(pdfExpectedEvents) {
171 out->addOwnedComponents(std::move(pdfExpectedEvents));
172 }
173 return out;
174}
175
176
178{
179 // Use the result of the underlying pdf.
180 ctx.addResult(this, ctx.getResult(_pdf));
181}
std::string _rangeName
Name of range in which to calculate test statistic.
#define coutW(a)
#define ClassImp(name)
Definition Rtypes.h:377
char name[80]
Definition TGX11.cxx:110
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Abstract interface for all probability density functions.
Definition RooAbsPdf.h:40
virtual double expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition RooAbsPdf.h:219
virtual std::unique_ptr< RooAbsReal > createExpectedEventsFunc(const RooArgSet *nset) const
Returns an object that represents the expected number of events for a given normalization set,...
virtual const RooAbsReal * getNormObj(const RooArgSet *set, const RooArgSet *iset, const TNamed *rangeName=nullptr) const
Return pointer to RooAbsReal object that implements calculation of integral over observables iset in ...
A RooAbsReal::Ref can be constructed from a RooAbsReal& or a double that will be implicitly converted...
Definition RooAbsReal.h:68
Abstract base class for objects that represent a real value and implements functionality common to al...
Definition RooAbsReal.h:59
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, const RooCmdArg &arg2={}, const RooCmdArg &arg3={}, const RooCmdArg &arg4={}, const RooCmdArg &arg5={}, const RooCmdArg &arg6={}, const RooCmdArg &arg7={}, const RooCmdArg &arg8={}) const
Create an object that represents the integral of the function over one or more observables listed in ...
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:103
void setUnit(const char *unit)
Definition RooAbsReal.h:147
const char * getPlotLabel() const
Get the label associated with the variable.
const Text_t * getUnit() const
Definition RooAbsReal.h:143
void setPlotLabel(const char *label)
Set the label associated with this variable.
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
RooExtendPdf is a wrapper around an existing PDF that adds a parameteric extended likelihood term to ...
std::unique_ptr< RooAbsReal > createExpectedEventsFunc(const RooArgSet *nset) const override
Returns an object that represents the expected number of events for a given normalization set,...
RooTemplateProxy< RooAbsPdf > _pdf
Input p.d.f.
double expectedEvents(const RooArgSet *nset) const override
Return the number of expected events over the full range of all variables.
const TNamed * _rangeName
Name of subset range.
RooExtendPdf()=default
RooTemplateProxy< RooAbsReal > _n
Number of expected events.
void translate(RooFit::Detail::CodeSquashContext &ctx) const override
This function defines a translation for each RooAbsReal based object that can be used to express the ...
A class to maintain the context for squashing of RooFit models into code.
void addResult(RooAbsArg const *key, std::string const &value)
A function to save an expression that includes/depends on the result of the input node.
std::string const & getResult(RooAbsArg const &arg)
Gets the result for the given node using the node name.
Registry for const char* names.
Definition RooNameReg.h:26
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition RooNameReg.h:39
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
RooConstVar & RooConst(double val)