Logo ROOT  
Reference Guide
RooAbsPdf.h
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * File: $Id: RooAbsPdf.h,v 1.90 2007/07/21 21:32:52 wouter Exp $
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#ifndef ROO_ABS_PDF
17#define ROO_ABS_PDF
18
19#include "RooAbsReal.h"
20#include "RooObjCacheManager.h"
21#include "RooCmdArg.h"
22
23class RooDataSet;
24class RooDataHist ;
25class RooArgSet ;
26class RooAbsGenContext ;
27class RooFitResult ;
28class RooExtendPdf ;
29class RooCategory ;
30class TPaveText;
31class TH1F;
32class TH2F;
33class TList ;
34class RooLinkedList ;
35class RooMinimizer ;
36class RooNumGenConfig ;
37class RooRealIntegral ;
38namespace RooBatchCompute {
39struct RunContext;
40}
41
42class RooAbsPdf : public RooAbsReal {
43public:
44
45 // Constructors, assignment etc
46 RooAbsPdf() ;
47 RooAbsPdf(const char *name, const char *title=0) ;
48 RooAbsPdf(const char *name, const char *title, Double_t minVal, Double_t maxVal) ;
49 // RooAbsPdf(const RooAbsPdf& other, const char* name=0);
50 virtual ~RooAbsPdf();
51
52 // Toy MC generation
53
54 ////////////////////////////////////////////////////////////////////////////////
55 /// See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&)
56 /// \param[in] nEvents How many events to generate
57 RooDataSet *generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg& arg1,
58 const RooCmdArg& arg2=RooCmdArg::none(), const RooCmdArg& arg3=RooCmdArg::none(),
59 const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none()) {
60 return generate(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5) ;
61 }
62 RooDataSet *generate(const RooArgSet &whatVars,
63 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
64 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
65 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) ;
66 RooDataSet *generate(const RooArgSet &whatVars, Double_t nEvents = 0, Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE,
67 const char* binnedTag="", Bool_t expectedData=kFALSE, Bool_t extended = kFALSE) const;
68 RooDataSet *generate(const RooArgSet &whatVars, const RooDataSet &prototype, Int_t nEvents= 0,
69 Bool_t verbose=kFALSE, Bool_t randProtoOrder=kFALSE, Bool_t resampleProto=kFALSE) const;
70
71
72 class GenSpec {
73 public:
74 virtual ~GenSpec() ;
76 private:
77 GenSpec(RooAbsGenContext* context, const RooArgSet& whatVars, RooDataSet* protoData, Int_t nGen, Bool_t extended,
78 Bool_t randProto, Bool_t resampleProto, TString dsetName, Bool_t init=kFALSE) ;
79 GenSpec(const GenSpec& other) ;
80
81 friend class RooAbsPdf ;
91 ClassDef(GenSpec,0) // Generation specification
92 } ;
93
94 ///Prepare GenSpec configuration object for efficient generation of multiple datasets from identical specification.
95 GenSpec* prepareMultiGen(const RooArgSet &whatVars,
96 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
97 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
98 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) ;
99 ///Generate according to GenSpec obtained from prepareMultiGen().
100 RooDataSet* generate(GenSpec&) const ;
101
102
103 ////////////////////////////////////////////////////////////////////////////////
104 /// As RooAbsPdf::generateBinned(const RooArgSet&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&)
105 /// \param[in] nEvents How many events to generate
106 virtual RooDataHist *generateBinned(const RooArgSet &whatVars, Double_t nEvents, const RooCmdArg& arg1,
107 const RooCmdArg& arg2=RooCmdArg::none(), const RooCmdArg& arg3=RooCmdArg::none(),
108 const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none()) const {
109 return generateBinned(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5);
110 }
111 virtual RooDataHist *generateBinned(const RooArgSet &whatVars,
112 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
113 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
114 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) const;
115 virtual RooDataHist *generateBinned(const RooArgSet &whatVars, Double_t nEvents, Bool_t expectedData=kFALSE, Bool_t extended=kFALSE) const;
116
117 virtual RooDataSet* generateSimGlobal(const RooArgSet& whatVars, Int_t nEvents) ;
118
119 ///Helper calling plotOn(RooPlot*, RooLinkedList&) const
120 virtual RooPlot* plotOn(RooPlot* frame,
121 const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
122 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
123 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
124 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none(),
125 const RooCmdArg& arg9=RooCmdArg::none(), const RooCmdArg& arg10=RooCmdArg::none()
126 ) const {
127 return RooAbsReal::plotOn(frame,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10) ;
128 }
129 virtual RooPlot* plotOn(RooPlot* frame, RooLinkedList& cmdList) const ;
130
131 /// Add a box with parameter values (and errors) to the specified frame
132 virtual RooPlot* paramOn(RooPlot* frame,
133 const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
134 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
135 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
136 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
137
138 virtual RooPlot* paramOn(RooPlot* frame, const RooAbsData* data, const char *label= "", Int_t sigDigits = 2,
139 Option_t *options = "NELU", Double_t xmin=0.65,
140 Double_t xmax = 0.9, Double_t ymax = 0.9) ;
141
142 // Built-in generator support
143 virtual Int_t getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const;
144 virtual void initGenerator(Int_t code) ;
145 virtual void generateEvent(Int_t code);
146 virtual Bool_t isDirectGenSafe(const RooAbsArg& arg) const ;
147
148 // Configuration of MC generators used for this pdf
149 const RooNumGenConfig* getGeneratorConfig() const ;
153 void setGeneratorConfig() ;
154 void setGeneratorConfig(const RooNumGenConfig& config) ;
155
156 // -log(L) fits to binned and unbinned data
157 virtual RooFitResult* fitTo(RooAbsData& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
158 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
159 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
160 virtual RooFitResult* fitTo(RooAbsData& data, const RooLinkedList& cmdList) ;
161
162 /// Configuration struct for RooAbsPdf::minimizeNLL with all the default
163 //values that also should be taked as the default values for
164 //RooAbsPdf::fitTo.
166 double recoverFromNaN = 10.;
167 std::string fitOpt = "";
168 int optConst = 2;
169 int verbose = 0;
170 int doSave = 0;
171 int doTimer = 0;
172 int printLevel = 1;
173 int strat = 1;
174 int initHesse = 0;
175 int hesse = 1;
176 int minos = 0;
177 int numee = 10;
178 int doEEWall = 1;
179 int doWarn = 1;
180 int doSumW2 = -1;
181 int doAsymptotic = -1;
182 const RooArgSet* minosSet = nullptr;
183 std::string minType = "Minuit";
184 std::string minAlg = "minuit";
185 };
186 std::unique_ptr<RooFitResult> minimizeNLL(RooAbsReal & nll, RooAbsData const& data, MinimizerConfig const& cfg);
187
188 virtual RooAbsReal* createNLL(RooAbsData& data, const RooLinkedList& cmdList) ;
189 virtual RooAbsReal* createNLL(RooAbsData& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
190 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
191 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
192
193 // Chi^2 fits to histograms
196 virtual RooFitResult* chi2FitTo(RooDataHist& data, const RooLinkedList& cmdList) ;
197 virtual RooAbsReal* createChi2(RooDataHist& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
198 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
199 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
200
201 // Chi^2 fits to X-Y datasets
202 virtual RooAbsReal* createChi2(RooDataSet& data, const RooLinkedList& cmdList) ;
203
204
205
206
207
208 // Constraint management
209 virtual RooArgSet* getConstraints(const RooArgSet& /*observables*/, RooArgSet& /*constrainedParams*/, Bool_t /*stripDisconnected*/) const {
210 // Interface to retrieve constraint terms on this pdf. Default implementation returns null
211 return 0 ;
212 }
213 virtual RooArgSet* getAllConstraints(const RooArgSet& observables, RooArgSet& constrainedParams, Bool_t stripDisconnected=kTRUE) const ;
214
215 // Project p.d.f into lower dimensional p.d.f
216 virtual RooAbsPdf* createProjection(const RooArgSet& iset) ;
217
218 // Create cumulative density function from p.d.f
219 RooAbsReal* createCdf(const RooArgSet& iset, const RooArgSet& nset=RooArgSet()) ;
220 RooAbsReal* createCdf(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2=RooCmdArg::none(),
221 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
222 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
223 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
224 RooAbsReal* createScanCdf(const RooArgSet& iset, const RooArgSet& nset, Int_t numScanBins, Int_t intOrder) ;
225
226 // Function evaluation support
227 virtual Double_t getValV(const RooArgSet* set=0) const ;
228 virtual Double_t getLogVal(const RooArgSet* set=0) const ;
229
231 RooSpan<const double> getLogValBatch(std::size_t begin, std::size_t batchSize,
232 const RooArgSet* normSet = nullptr) const;
233 RooSpan<const double> getLogProbabilities(RooBatchCompute::RunContext& evalData, const RooArgSet* normSet = nullptr) const;
234
235 /// \copydoc getNorm(const RooArgSet*) const
236 Double_t getNorm(const RooArgSet& nset) const {
237 return getNorm(&nset) ;
238 }
239 virtual Double_t getNorm(const RooArgSet* set=0) const ;
240
241 virtual void resetErrorCounters(Int_t resetValue=10) ;
242 void setTraceCounter(Int_t value, Bool_t allNodes=kFALSE) ;
243private:
244 Bool_t traceEvalPdf(Double_t value) const;
245
246public:
247
248 Double_t analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=0) const ;
249
250 /// Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
251 /// Always returns false, unless a PDF overrides this function.
252 virtual Bool_t selfNormalized() const {
253 return kFALSE ;
254 }
255
256 // Support for extended maximum likelihood, switched off by default
258 /// Returns ability of PDF to provide extended likelihood terms. Possible
259 /// answers are in the enumerator RooAbsPdf::ExtendMode.
260 /// This default implementation always returns CanNotBeExtended.
261 virtual ExtendMode extendMode() const { return CanNotBeExtended; }
262 /// If true, PDF can provide extended likelihood term.
263 inline Bool_t canBeExtended() const {
264 return (extendMode() != CanNotBeExtended) ;
265 }
266 /// If true PDF must provide extended likelihood term.
267 inline Bool_t mustBeExtended() const {
268 return (extendMode() == MustBeExtended) ;
269 }
270 /// Return expected number of events to be used in calculation of extended
271 /// likelihood.
272 virtual Double_t expectedEvents(const RooArgSet* nset) const ;
273 /// Return expected number of events to be used in calculation of extended
274 /// likelihood. This function should not be overridden, as it just redirects
275 /// to the actual virtual function but takes a RooArgSet reference instead of
276 /// pointer (\see expectedEvents(const RooArgSet*) const).
277 double expectedEvents(const RooArgSet& nset) const {
278 return expectedEvents(&nset) ;
279 }
280
281 // Printing interface (human readable)
282 virtual void printValue(std::ostream& os) const ;
283 virtual void printMultiline(std::ostream& os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const ;
284
285 static void verboseEval(Int_t stat) ;
286 static int verboseEval() ;
287
288 virtual Double_t extendedTerm(Double_t observedEvents, const RooArgSet* nset=0) const ;
289
290 void setNormRange(const char* rangeName) ;
291 const char* normRange() const {
292 return _normRange.Length()>0 ? _normRange.Data() : 0 ;
293 }
294 void setNormRangeOverride(const char* rangeName) ;
295
296 const RooAbsReal* getNormIntegral(const RooArgSet& nset) const { return getNormObj(0,&nset,0) ; }
297
298protected:
299
300public:
301 virtual const RooAbsReal* getNormObj(const RooArgSet* set, const RooArgSet* iset, const TNamed* rangeName=0) const ;
302
303 virtual RooAbsGenContext* binnedGenContext(const RooArgSet &vars, Bool_t verbose= kFALSE) const ;
304
305 virtual RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=0,
306 const RooArgSet* auxProto=0, Bool_t verbose= kFALSE) const ;
307
308 virtual RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=0, const RooArgSet* auxProto=0,
309 Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE, const char* binnedTag="") const ;
310
311private:
312
313 RooDataSet *generate(RooAbsGenContext& context, const RooArgSet& whatVars, const RooDataSet* prototype,
314 Double_t nEvents, Bool_t verbose, Bool_t randProtoOrder, Bool_t resampleProto, Bool_t skipInit=kFALSE,
315 Bool_t extended=kFALSE) const ;
316
317 // Implementation version
318 virtual RooPlot* paramOn(RooPlot* frame, const RooArgSet& params, Bool_t showConstants=kFALSE,
319 const char *label= "", Int_t sigDigits = 2, Option_t *options = "NELU", Double_t xmin=0.65,
320 Double_t xmax= 0.99,Double_t ymax=0.95, const RooCmdArg* formatCmd=0) ;
321
322 void logBatchComputationErrors(RooSpan<const double>& outputs, std::size_t begin) const;
323
324protected:
325 virtual RooPlot *plotOn(RooPlot *frame, PlotOpt o) const;
326
327 friend class RooEffGenContext ;
328 friend class RooAddGenContext ;
329 friend class RooProdGenContext ;
330 friend class RooSimGenContext ;
332 friend class RooConvGenContext ;
333 friend class RooSimultaneous ;
334 friend class RooAddGenContextOrig ;
335 friend class RooProdPdf ;
336 friend class RooMCStudy ;
337
338 Int_t* randomizeProtoOrder(Int_t nProto,Int_t nGen,Bool_t resample=kFALSE) const ;
339
340 friend class RooExtendPdf ;
341 // This also forces the definition of a copy ctor in derived classes
342 RooAbsPdf(const RooAbsPdf& other, const char* name = 0);
343
344 friend class RooRealIntegral ;
346
347 virtual Bool_t syncNormalization(const RooArgSet* dset, Bool_t adjustProxies=kTRUE) const ;
348
349 friend class RooAbsAnaConvPdf ;
351 mutable RooAbsReal* _norm ; //! Normalization integral (owned by _normMgr)
352 mutable RooArgSet* _normSet ; //! Normalization set with for above integral
353
355 public:
356 CacheElem(RooAbsReal& norm) : _norm(&norm) {} ;
357 virtual ~CacheElem() ;
360 } ;
361 mutable RooObjCacheManager _normMgr ; //! The cache manager
362
363 friend class CacheElem ; // Cache needs to be able to clear _norm pointer
364
366 // Hook function intercepting redirectServer calls. Discard current normalization
367 // object if any server is redirected
368
369 // Object is own by _normCacheManager that will delete object as soon as cache
370 // is sterilized by server redirect
371 _norm = 0 ;
372 return kFALSE ;
373 } ;
374
375
376 mutable Int_t _errorCount ; // Number of errors remaining to print
377 mutable Int_t _traceCount ; // Number of traces remaining to print
378 mutable Int_t _negCount ; // Number of negative probablities remaining to print
379
380 Bool_t _selectComp ; // Component selection flag for RooAbsPdf::plotCompOn
381
382 RooNumGenConfig* _specGeneratorConfig ; //! MC generator configuration specific for this object
383
384 TString _normRange ; // Normalization range
386
387private:
388 int calcAsymptoticCorrectedCovariance(RooMinimizer& minimizer, RooAbsData const& data);
389 int calcSumW2CorrectedCovariance(RooMinimizer& minimizer, RooAbsReal const& nll) const;
390
391 ClassDef(RooAbsPdf,5) // Abstract PDF with normalization support
392};
393
394
395#endif
int Int_t
Definition: RtypesCore.h:45
const Bool_t kFALSE
Definition: RtypesCore.h:101
bool Bool_t
Definition: RtypesCore.h:63
double Double_t
Definition: RtypesCore.h:59
const Bool_t kTRUE
Definition: RtypesCore.h:100
const char Option_t
Definition: RtypesCore.h:66
#define ClassDef(name, id)
Definition: Rtypes.h:325
static void indent(ostringstream &buf, int indent_level)
char name[80]
Definition: TGX11.cxx:110
float xmin
Definition: THbookFile.cxx:95
float xmax
Definition: THbookFile.cxx:95
float ymax
Definition: THbookFile.cxx:95
RooAbsAnaConvPdf is the base class for PDFs that represent a physics model that can be analytically c...
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
friend class RooArgSet
Definition: RooAbsArg.h:600
RooAbsCacheElement is the abstract base class for objects to be stored in RooAbsCache cache manager o...
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:79
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
Normalization set with for above integral.
Definition: RooAbsPdf.h:354
virtual ~CacheElem()
Destructor of normalization cache element.
Definition: RooAbsPdf.cxx:3269
virtual RooArgList containedArgs(Action)
Definition: RooAbsPdf.h:358
CacheElem(RooAbsReal &norm)
Definition: RooAbsPdf.h:356
RooAbsReal * _norm
Definition: RooAbsPdf.h:359
RooArgSet _whatVars
Definition: RooAbsPdf.h:83
RooAbsGenContext * _genContext
Definition: RooAbsPdf.h:82
Bool_t _resampleProto
Definition: RooAbsPdf.h:88
GenSpec(const GenSpec &other)
TString _dsetName
Definition: RooAbsPdf.h:89
RooDataSet * _protoData
Definition: RooAbsPdf.h:84
GenSpec * prepareMultiGen(const RooArgSet &whatVars, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none())
Prepare GenSpec configuration object for efficient generation of multiple datasets from identical spe...
Definition: RooAbsPdf.cxx:2081
RooSpan< const double > getLogValBatch(std::size_t begin, std::size_t batchSize, const RooArgSet *normSet=nullptr) const
RooObjCacheManager _normMgr
Definition: RooAbsPdf.h:361
Double_t getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition: RooAbsPdf.h:236
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooAbsPdf.cxx:2330
virtual ~RooAbsPdf()
Destructor.
Definition: RooAbsPdf.cxx:265
double expectedEvents(const RooArgSet &nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.h:277
void logBatchComputationErrors(RooSpan< const double > &outputs, std::size_t begin) const
Scan through outputs and fix+log all nans and negative values.
Definition: RooAbsPdf.cxx:718
RooSpan< const double > getLogProbabilities(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet=nullptr) const
Compute the log-likelihoods for all events in the requested batch.
Definition: RooAbsPdf.cxx:741
RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const
Compute batch of values for given input data, and normalise by integrating over the observables in no...
Definition: RooAbsPdf.cxx:355
void setGeneratorConfig()
Remove the specialized numeric MC generator configuration associated with this object.
Definition: RooAbsPdf.cxx:3503
virtual void resetErrorCounters(Int_t resetValue=10)
Reset error counter to given value, limiting the number of future error messages for this pdf to 'res...
Definition: RooAbsPdf.cxx:623
virtual RooFitResult * chi2FitTo(RooDataHist &data, const RooLinkedList &cmdList)
Calls RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList) and returns fit result.
Definition: RooAbsPdf.cxx:1712
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
Definition: RooAbsPdf.cxx:3257
Bool_t traceEvalPdf(Double_t value) const
Check that passed value is positive and not 'not-a-number'.
Definition: RooAbsPdf.cxx:433
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
Definition: RooAbsPdf.cxx:988
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE) const
Interface function to create a generator context from a p.d.f.
Definition: RooAbsPdf.cxx:1909
RooAbsReal * createScanCdf(const RooArgSet &iset, const RooArgSet &nset, Int_t numScanBins, Int_t intOrder)
Definition: RooAbsPdf.cxx:3399
TString _normRange
MC generator configuration specific for this object.
Definition: RooAbsPdf.h:384
void setNormRange(const char *rangeName)
Definition: RooAbsPdf.cxx:3534
const RooAbsReal * getNormIntegral(const RooArgSet &nset) const
Definition: RooAbsPdf.h:296
virtual RooDataHist * generateBinned(const RooArgSet &whatVars, Double_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none()) const
As RooAbsPdf::generateBinned(const RooArgSet&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,...
Definition: RooAbsPdf.h:106
virtual RooFitResult * fitTo(RooAbsData &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Fit PDF to given dataset.
Definition: RooAbsPdf.cxx:1456
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooAbsPdf.cxx:2308
RooNumGenConfig * specialGeneratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
Definition: RooAbsPdf.cxx:3448
Bool_t canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:263
virtual Bool_t selfNormalized() const
Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
Definition: RooAbsPdf.h:252
virtual void printMultiline(std::ostream &os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const
Print multi line detailed information of this RooAbsPdf.
Definition: RooAbsPdf.cxx:1882
Double_t analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const
Analytical integral with normalization (see RooAbsReal::analyticalIntegralWN() for further informatio...
Definition: RooAbsPdf.cxx:413
virtual RooPlot * paramOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Add a box with parameter values (and errors) to the specified frame.
Definition: RooAbsPdf.cxx:3085
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,...
Definition: RooAbsPdf.h:57
virtual RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE, const char *binnedTag="") const
Definition: RooAbsPdf.cxx:1918
Int_t _traceCount
Definition: RooAbsPdf.h:377
RooAbsReal * _norm
Definition: RooAbsPdf.h:351
int calcAsymptoticCorrectedCovariance(RooMinimizer &minimizer, RooAbsData const &data)
Use the asymptotically correct approach to estimate errors in the presence of weights.
Definition: RooAbsPdf.cxx:1184
virtual Double_t extendedTerm(Double_t observedEvents, const RooArgSet *nset=0) const
Return the extended likelihood term ( ) of this PDF for the given number of observed events.
Definition: RooAbsPdf.cxx:784
virtual void printValue(std::ostream &os) const
Print value of p.d.f, also print normalization integral that was last used, if any.
Definition: RooAbsPdf.cxx:1863
virtual RooArgSet * getAllConstraints(const RooArgSet &observables, RooArgSet &constrainedParams, Bool_t stripDisconnected=kTRUE) const
This helper function finds and collects all constraints terms of all component p.d....
Definition: RooAbsPdf.cxx:3416
virtual Bool_t syncNormalization(const RooArgSet *dset, Bool_t adjustProxies=kTRUE) const
Verify that the normalization integral cached with this PDF is valid for given set of normalization o...
Definition: RooAbsPdf.cxx:527
virtual Double_t getValV(const RooArgSet *set=0) const
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:283
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none(), const RooCmdArg &arg9=RooCmdArg::none(), const RooCmdArg &arg10=RooCmdArg::none()) const
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition: RooAbsPdf.h:120
Int_t _errorCount
Definition: RooAbsPdf.h:373
Bool_t _selectComp
Definition: RooAbsPdf.h:380
virtual RooArgSet * getConstraints(const RooArgSet &, RooArgSet &, Bool_t) const
Definition: RooAbsPdf.h:209
@ CanBeExtended
Definition: RooAbsPdf.h:257
@ MustBeExtended
Definition: RooAbsPdf.h:257
@ CanNotBeExtended
Definition: RooAbsPdf.h:257
void setTraceCounter(Int_t value, Bool_t allNodes=kFALSE)
Reset trace counter to given value, limiting the number of future trace messages for this pdf to 'val...
Definition: RooAbsPdf.cxx:635
const char * normRange() const
Definition: RooAbsPdf.h:291
RooAbsReal * createCdf(const RooArgSet &iset, const RooArgSet &nset=RooArgSet())
Create a cumulative distribution function of this p.d.f in terms of the observables listed in iset.
Definition: RooAbsPdf.cxx:3322
virtual const RooAbsReal * getNormObj(const RooArgSet *set, const RooArgSet *iset, const TNamed *rangeName=0) const
Return pointer to RooAbsReal object that implements calculation of integral over observables iset in ...
Definition: RooAbsPdf.cxx:492
virtual Double_t getLogVal(const RooArgSet *set=0) const
Return the log of the current value with given normalization An error message is printed if the argum...
Definition: RooAbsPdf.cxx:661
Int_t * randomizeProtoOrder(Int_t nProto, Int_t nGen, Bool_t resample=kFALSE) const
Return lookup table with randomized order for nProto prototype events.
Definition: RooAbsPdf.cxx:2271
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.cxx:3237
Int_t _negCount
Definition: RooAbsPdf.h:378
std::unique_ptr< RooFitResult > minimizeNLL(RooAbsReal &nll, RooAbsData const &data, MinimizerConfig const &cfg)
Minimizes a given NLL variable by finding the optimal parameters with the RooMinimzer.
Definition: RooAbsPdf.cxx:1479
virtual RooAbsReal * createChi2(RooDataHist &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Create a from a histogram and this function.
Definition: RooAbsPdf.cxx:1758
Bool_t mustBeExtended() const
If true PDF must provide extended likelihood term.
Definition: RooAbsPdf.h:267
void setNormRangeOverride(const char *rangeName)
Definition: RooAbsPdf.cxx:3551
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of 'global observables' – for RooStats tools.
Definition: RooAbsPdf.cxx:2587
virtual Bool_t isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
Definition: RooAbsPdf.cxx:2343
const RooNumGenConfig * getGeneratorConfig() const
Return the numeric MC generator configuration used for this object.
Definition: RooAbsPdf.cxx:3476
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
Definition: RooAbsPdf.cxx:2318
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition: RooAbsPdf.h:261
RooAbsPdf()
Default constructor.
Definition: RooAbsPdf.cxx:207
virtual RooAbsGenContext * binnedGenContext(const RooArgSet &vars, Bool_t verbose=kFALSE) const
Return a binned generator context.
Definition: RooAbsPdf.cxx:1899
static RooNumGenConfig * defaultGeneratorConfig()
Returns the default numeric MC generator configuration for all RooAbsReals.
Definition: RooAbsPdf.cxx:3438
RooNumGenConfig * _specGeneratorConfig
Definition: RooAbsPdf.h:382
RooArgSet * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:352
friend class RooAddGenContextOrig
Definition: RooAbsPdf.h:334
Double_t _rawValue
Definition: RooAbsPdf.h:350
static TString _normRangeOverride
Definition: RooAbsPdf.h:385
static Int_t _verboseEval
Definition: RooAbsPdf.h:345
virtual Bool_t redirectServersHook(const RooAbsCollection &, Bool_t, Bool_t, Bool_t)
Definition: RooAbsPdf.h:365
int calcSumW2CorrectedCovariance(RooMinimizer &minimizer, RooAbsReal const &nll) const
Apply correction to errors and covariance matrix.
Definition: RooAbsPdf.cxx:1264
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Return a p.d.f that represent a projection of this p.d.f integrated over given observables.
Definition: RooAbsPdf.cxx:3287
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:61
virtual RooAbsReal * createChi2(RooDataHist &data, const RooLinkedList &cmdList)
Create a variable from a histogram and this function.
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg(), const RooCmdArg &arg10=RooCmdArg()) const
Plot (project) PDF on specified frame.
virtual Double_t maxVal(Int_t code) const
Return maximum value for set of observables identified by code assigned in getMaxVal.
virtual RooFitResult * chi2FitTo(RooDataHist &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Perform a fit to given histogram.
RooAddGenContext is an efficient implementation of the generator context specific for RooAddPdf PDFs.
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:35
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:27
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:27
static const RooCmdArg & none()
Return reference to null argument.
Definition: RooCmdArg.cxx:52
RooConvGenContext is an efficient implementation of the generator context specific for RooAbsAnaConvP...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:44
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:36
RooEffGenContext is a specialized generator context for p.d.fs represented by class RooEffProd,...
RooExtendPdf is a wrapper around an existing PDF that adds a parameteric extended likelihood term to ...
Definition: RooExtendPdf.h:22
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
Definition: RooFitResult.h:40
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
Definition: RooLinkedList.h:38
RooMCStudy is a helper class to facilitate Monte Carlo studies such as 'goodness-of-fit' studies,...
Definition: RooMCStudy.h:32
RooMinimizer is a wrapper class around ROOT::Fit:Fitter that provides a seamless interface between th...
Definition: RooMinimizer.h:55
RooNumGenConfig holds the configuration parameters of the various numeric integrators used by RooReal...
Class RooObjCacheManager is an implementation of class RooCacheManager<RooAbsCacheElement> and specia...
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition: RooPlot.h:44
RooProdGenContext is an efficient implementation of the generator context specific for RooProdPdf PDF...
RooProdPdf is an efficient implementation of a product of PDFs of the form.
Definition: RooProdPdf.h:37
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
RooSimGenContext is an efficient implementation of the generator context specific for RooSimultaneous...
RooSimSplitGenContext is an efficient implementation of the generator context specific for RooSimulta...
RooSimultaneous facilitates simultaneous fitting of multiple PDFs to subsets of a given dataset.
A simple container to hold a batch of data values.
Definition: RooSpan.h:34
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:575
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:251
A doubly linked list.
Definition: TList.h:44
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:29
A Pave (see TPave) with text, lines or/and boxes inside.
Definition: TPaveText.h:21
Basic string class.
Definition: TString.h:136
Ssiz_t Length() const
Definition: TString.h:410
const char * Data() const
Definition: TString.h:369
RooCmdArg NumEvents(Int_t numEvents)
Namespace for dispatching RooFit computations to various backends.
Configuration struct for RooAbsPdf::minimizeNLL with all the default.
Definition: RooAbsPdf.h:165
const RooArgSet * minosSet
Definition: RooAbsPdf.h:182
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:31