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 ;
38
39namespace RooBatchCompute {
40struct RunContext;
41}
42
43class RooAbsPdf : public RooAbsReal {
44public:
45
46 // Constructors, assignment etc
47 RooAbsPdf() ;
48 RooAbsPdf(const char *name, const char *title=0) ;
49 RooAbsPdf(const char *name, const char *title, Double_t minVal, Double_t maxVal) ;
50 // RooAbsPdf(const RooAbsPdf& other, const char* name=0);
51 virtual ~RooAbsPdf();
52
53 // Toy MC generation
54
55 ////////////////////////////////////////////////////////////////////////////////
56 /// See RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&)
57 /// \param[in] nEvents How many events to generate
58 RooDataSet *generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg& arg1,
59 const RooCmdArg& arg2=RooCmdArg::none(), const RooCmdArg& arg3=RooCmdArg::none(),
60 const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none()) {
61 return generate(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5) ;
62 }
63 RooDataSet *generate(const RooArgSet &whatVars,
64 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
65 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
66 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) ;
67 RooDataSet *generate(const RooArgSet &whatVars, Double_t nEvents = 0, Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE,
68 const char* binnedTag="", Bool_t expectedData=kFALSE, Bool_t extended = kFALSE) const;
69 RooDataSet *generate(const RooArgSet &whatVars, const RooDataSet &prototype, Int_t nEvents= 0,
70 Bool_t verbose=kFALSE, Bool_t randProtoOrder=kFALSE, Bool_t resampleProto=kFALSE) const;
71
72
73 class GenSpec {
74 public:
75 virtual ~GenSpec() ;
77 private:
78 GenSpec(RooAbsGenContext* context, const RooArgSet& whatVars, RooDataSet* protoData, Int_t nGen, Bool_t extended,
79 Bool_t randProto, Bool_t resampleProto, TString dsetName, Bool_t init=kFALSE) ;
80 GenSpec(const GenSpec& other) ;
81
82 friend class RooAbsPdf ;
92 ClassDef(GenSpec,0) // Generation specification
93 } ;
94
95 ///Prepare GenSpec configuration object for efficient generation of multiple datasets from identical specification.
96 GenSpec* prepareMultiGen(const RooArgSet &whatVars,
97 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
98 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
99 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) ;
100 ///Generate according to GenSpec obtained from prepareMultiGen().
101 RooDataSet* generate(GenSpec&) const ;
102
103
104 ////////////////////////////////////////////////////////////////////////////////
105 /// As RooAbsPdf::generateBinned(const RooArgSet&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&)
106 /// \param[in] nEvents How many events to generate
107 virtual RooDataHist *generateBinned(const RooArgSet &whatVars, Double_t nEvents, const RooCmdArg& arg1,
108 const RooCmdArg& arg2=RooCmdArg::none(), const RooCmdArg& arg3=RooCmdArg::none(),
109 const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none()) const {
110 return generateBinned(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5);
111 }
112 virtual RooDataHist *generateBinned(const RooArgSet &whatVars,
113 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
114 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
115 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) const;
116 virtual RooDataHist *generateBinned(const RooArgSet &whatVars, Double_t nEvents, Bool_t expectedData=kFALSE, Bool_t extended=kFALSE) const;
117
118 virtual RooDataSet* generateSimGlobal(const RooArgSet& whatVars, Int_t nEvents) ;
119
120 ///Helper calling plotOn(RooPlot*, RooLinkedList&) const
121 virtual RooPlot* plotOn(RooPlot* frame,
122 const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
123 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
124 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
125 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none(),
126 const RooCmdArg& arg9=RooCmdArg::none(), const RooCmdArg& arg10=RooCmdArg::none()
127 ) const {
128 return RooAbsReal::plotOn(frame,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10) ;
129 }
130 virtual RooPlot* plotOn(RooPlot* frame, RooLinkedList& cmdList) const ;
131
132 /// Add a box with parameter values (and errors) to the specified frame
133 virtual RooPlot* paramOn(RooPlot* frame,
134 const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
135 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
136 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
137 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
138
139 virtual RooPlot* paramOn(RooPlot* frame, const RooAbsData* data, const char *label= "", Int_t sigDigits = 2,
140 Option_t *options = "NELU", Double_t xmin=0.65,
141 Double_t xmax = 0.9, Double_t ymax = 0.9) ;
142
143 // Built-in generator support
144 virtual Int_t getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, Bool_t staticInitOK=kTRUE) const;
145 virtual void initGenerator(Int_t code) ;
146 virtual void generateEvent(Int_t code);
147 virtual Bool_t isDirectGenSafe(const RooAbsArg& arg) const ;
148
149 // Configuration of MC generators used for this pdf
150 const RooNumGenConfig* getGeneratorConfig() const ;
154 void setGeneratorConfig() ;
155 void setGeneratorConfig(const RooNumGenConfig& config) ;
156
157 // -log(L) fits to binned and unbinned data
158 virtual RooFitResult* fitTo(RooAbsData& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
159 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
160 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
161 virtual RooFitResult* fitTo(RooAbsData& data, const RooLinkedList& cmdList) ;
162
163 /// Configuration struct for RooAbsPdf::minimizeNLL with all the default
164 //values that also should be taked as the default values for
165 //RooAbsPdf::fitTo.
167 double recoverFromNaN = 10.;
168 std::string fitOpt = "";
169 int optConst = 2;
170 int verbose = 0;
171 int doSave = 0;
172 int doTimer = 0;
173 int printLevel = 1;
174 int strat = 1;
175 int initHesse = 0;
176 int hesse = 1;
177 int minos = 0;
178 int numee = 10;
179 int doEEWall = 1;
180 int doWarn = 1;
181 int doSumW2 = -1;
182 int doAsymptotic = -1;
183 const RooArgSet* minosSet = nullptr;
184 std::string minType = "Minuit";
185 std::string minAlg = "minuit";
186 };
187 std::unique_ptr<RooFitResult> minimizeNLL(RooAbsReal & nll, RooAbsData const& data, MinimizerConfig const& cfg);
188
189 virtual RooAbsReal* createNLL(RooAbsData& data, const RooLinkedList& cmdList) ;
190 virtual RooAbsReal* createNLL(RooAbsData& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
191 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
192 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
193
194 // Chi^2 fits to histograms
197 virtual RooFitResult* chi2FitTo(RooDataHist& data, const RooLinkedList& cmdList) ;
198 virtual RooAbsReal* createChi2(RooDataHist& data, const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
199 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
200 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
201
202 // Chi^2 fits to X-Y datasets
203 virtual RooAbsReal* createChi2(RooDataSet& data, const RooLinkedList& cmdList) ;
204
205
206
207
208
209 // Constraint management
210 virtual RooArgSet* getConstraints(const RooArgSet& /*observables*/, RooArgSet& /*constrainedParams*/, Bool_t /*stripDisconnected*/) const {
211 // Interface to retrieve constraint terms on this pdf. Default implementation returns null
212 return 0 ;
213 }
214 virtual RooArgSet* getAllConstraints(const RooArgSet& observables, RooArgSet& constrainedParams, Bool_t stripDisconnected=kTRUE) const ;
215
216 // Project p.d.f into lower dimensional p.d.f
217 virtual RooAbsPdf* createProjection(const RooArgSet& iset) ;
218
219 // Create cumulative density function from p.d.f
220 RooAbsReal* createCdf(const RooArgSet& iset, const RooArgSet& nset=RooArgSet()) ;
221 RooAbsReal* createCdf(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2=RooCmdArg::none(),
222 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
223 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
224 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
225 RooAbsReal* createScanCdf(const RooArgSet& iset, const RooArgSet& nset, Int_t numScanBins, Int_t intOrder) ;
226
227 // Function evaluation support
228 virtual Double_t getValV(const RooArgSet* set=0) const ;
229 virtual Double_t getLogVal(const RooArgSet* set=0) const ;
230
233 RooSpan<const double> getLogValBatch(std::size_t begin, std::size_t batchSize,
234 const RooArgSet* normSet = nullptr) const;
235 RooSpan<const double> getLogProbabilities(RooBatchCompute::RunContext& evalData, const RooArgSet* normSet = nullptr) const;
236 void getLogProbabilities(RooSpan<const double> pdfValues, double * output) const;
237
238 void computeBatch(cudaStream_t*, double* output, size_t size, RooBatchCompute::DataMap&) const;
239
240 /// \copydoc getNorm(const RooArgSet*) const
241 Double_t getNorm(const RooArgSet& nset) const {
242 return getNorm(&nset) ;
243 }
244 virtual Double_t getNorm(const RooArgSet* set=0) const ;
245 inline RooAbsReal* getIntegral(RooArgSet const& set) const {
246 syncNormalization(&set,true) ;
247 getVal(set);
248 assert(_norm != nullptr);
249 return _norm;
250 }
252 return _norm;
253 }
254
255
256 virtual void resetErrorCounters(Int_t resetValue=10) ;
257 void setTraceCounter(Int_t value, Bool_t allNodes=kFALSE) ;
258
259 Double_t analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=0) const ;
260
261 /// Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
262 /// Always returns false, unless a PDF overrides this function.
263 virtual Bool_t selfNormalized() const {
264 return kFALSE ;
265 }
266
267 // Support for extended maximum likelihood, switched off by default
269 /// Returns ability of PDF to provide extended likelihood terms. Possible
270 /// answers are in the enumerator RooAbsPdf::ExtendMode.
271 /// This default implementation always returns CanNotBeExtended.
272 virtual ExtendMode extendMode() const { return CanNotBeExtended; }
273 /// If true, PDF can provide extended likelihood term.
274 inline Bool_t canBeExtended() const {
275 return (extendMode() != CanNotBeExtended) ;
276 }
277 /// If true PDF must provide extended likelihood term.
278 inline Bool_t mustBeExtended() const {
279 return (extendMode() == MustBeExtended) ;
280 }
281 /// Return expected number of events to be used in calculation of extended
282 /// likelihood.
283 virtual Double_t expectedEvents(const RooArgSet* nset) const ;
284 /// Return expected number of events to be used in calculation of extended
285 /// likelihood. This function should not be overridden, as it just redirects
286 /// to the actual virtual function but takes a RooArgSet reference instead of
287 /// pointer (\see expectedEvents(const RooArgSet*) const).
288 double expectedEvents(const RooArgSet& nset) const {
289 return expectedEvents(&nset) ;
290 }
291
292 // Printing interface (human readable)
293 virtual void printValue(std::ostream& os) const ;
294 virtual void printMultiline(std::ostream& os, Int_t contents, Bool_t verbose=kFALSE, TString indent="") const ;
295
296 static void verboseEval(Int_t stat) ;
297 static int verboseEval() ;
298
299 double extendedTerm(double observedEvents, const RooArgSet* nset=0) const;
300 double extendedTerm(RooAbsData const& data, bool weightSquared) const;
301
302 void setNormRange(const char* rangeName) ;
303 const char* normRange() const {
304 return _normRange.Length()>0 ? _normRange.Data() : 0 ;
305 }
306 void setNormRangeOverride(const char* rangeName) ;
307
308 const RooAbsReal* getNormIntegral(const RooArgSet& nset) const { return getNormObj(0,&nset,0) ; }
309
310 virtual const RooAbsReal* getNormObj(const RooArgSet* set, const RooArgSet* iset, const TNamed* rangeName=0) const ;
311
312 virtual RooAbsGenContext* binnedGenContext(const RooArgSet &vars, Bool_t verbose= kFALSE) const ;
313
314 virtual RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=0,
315 const RooArgSet* auxProto=0, Bool_t verbose= kFALSE) const ;
316
317 virtual RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=0, const RooArgSet* auxProto=0,
318 Bool_t verbose=kFALSE, Bool_t autoBinned=kTRUE, const char* binnedTag="") const ;
319
320private:
321
322 RooDataSet *generate(RooAbsGenContext& context, const RooArgSet& whatVars, const RooDataSet* prototype,
323 Double_t nEvents, Bool_t verbose, Bool_t randProtoOrder, Bool_t resampleProto, Bool_t skipInit=kFALSE,
324 Bool_t extended=kFALSE) const ;
325
326 // Implementation version
327 virtual RooPlot* paramOn(RooPlot* frame, const RooArgSet& params, Bool_t showConstants=kFALSE,
328 const char *label= "", Int_t sigDigits = 2, Option_t *options = "NELU", Double_t xmin=0.65,
329 Double_t xmax= 0.99,Double_t ymax=0.95, const RooCmdArg* formatCmd=0) ;
330
331 void logBatchComputationErrors(RooSpan<const double>& outputs, std::size_t begin) const;
332 Bool_t traceEvalPdf(Double_t value) const;
333
334
335protected:
336 double normalizeWithNaNPacking(double rawVal, double normVal) const;
337
338 virtual RooPlot *plotOn(RooPlot *frame, PlotOpt o) const;
339
340 friend class RooEffGenContext ;
341 friend class RooAddGenContext ;
342 friend class RooProdGenContext ;
343 friend class RooSimGenContext ;
345 friend class RooConvGenContext ;
346 friend class RooSimultaneous ;
347 friend class RooAddGenContextOrig ;
348 friend class RooProdPdf ;
349 friend class RooMCStudy ;
350
351 Int_t* randomizeProtoOrder(Int_t nProto,Int_t nGen,Bool_t resample=kFALSE) const ;
352
353 friend class RooExtendPdf ;
354 // This also forces the definition of a copy ctor in derived classes
355 RooAbsPdf(const RooAbsPdf& other, const char* name = 0);
356
357 friend class RooRealIntegral ;
359
360 virtual Bool_t syncNormalization(const RooArgSet* dset, Bool_t adjustProxies=kTRUE) const ;
361
362 friend class RooAbsAnaConvPdf ;
364 mutable RooAbsReal* _norm = nullptr; //! Normalization integral (owned by _normMgr)
365 mutable RooArgSet const* _normSet = nullptr; //! Normalization set with for above integral
366 inline const RooArgSet* getNormSet() { return _normSet; }
367
369 public:
370 CacheElem(RooAbsReal& norm) : _norm(&norm) {} ;
371 virtual ~CacheElem() ;
374 } ;
375 mutable RooObjCacheManager _normMgr ; //! The cache manager
376
377 friend class CacheElem ; // Cache needs to be able to clear _norm pointer
378
380 // Hook function intercepting redirectServer calls. Discard current normalization
381 // object if any server is redirected
382
383 // Object is own by _normCacheManager that will delete object as soon as cache
384 // is sterilized by server redirect
385 _norm = nullptr ;
386
387 // Similar to the situation with the normalization integral above: if a
388 // server is redirected, the cached normalization set might not point to
389 // the right observables anymore. We need to reset it.
390 _normSet = nullptr ;
391 return false ;
392 } ;
393
394
395 mutable Int_t _errorCount ; // Number of errors remaining to print
396 mutable Int_t _traceCount ; // Number of traces remaining to print
397 mutable Int_t _negCount ; // Number of negative probablities remaining to print
398
399 Bool_t _selectComp ; // Component selection flag for RooAbsPdf::plotCompOn
400
401 RooNumGenConfig* _specGeneratorConfig ; //! MC generator configuration specific for this object
402
403 TString _normRange ; // Normalization range
405
406private:
407 int calcAsymptoticCorrectedCovariance(RooMinimizer& minimizer, RooAbsData const& data);
408 int calcSumW2CorrectedCovariance(RooMinimizer& minimizer, RooAbsReal const& nll) const;
409
410 ClassDef(RooAbsPdf,5) // Abstract PDF with normalization support
411};
412
413
414
415
416#endif
size_t size(const MatrixT &matrix)
retrieve the size of a square matrix
const Bool_t kFALSE
Definition: RtypesCore.h:101
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:74
friend class RooArgSet
Definition: RooAbsArg.h:624
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:80
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
virtual ~CacheElem()
Destructor of normalization cache element.
Definition: RooAbsPdf.cxx:3351
virtual RooArgList containedArgs(Action)
Definition: RooAbsPdf.h:372
CacheElem(RooAbsReal &norm)
Definition: RooAbsPdf.h:370
RooAbsReal * _norm
Definition: RooAbsPdf.h:373
RooArgSet _whatVars
Definition: RooAbsPdf.h:84
RooAbsGenContext * _genContext
Definition: RooAbsPdf.h:83
Bool_t _resampleProto
Definition: RooAbsPdf.h:89
GenSpec(const GenSpec &other)
TString _dsetName
Definition: RooAbsPdf.h:90
RooDataSet * _protoData
Definition: RooAbsPdf.h:85
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:2163
RooSpan< const double > getLogValBatch(std::size_t begin, std::size_t batchSize, const RooArgSet *normSet=nullptr) const
RooObjCacheManager _normMgr
Definition: RooAbsPdf.h:375
Double_t getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition: RooAbsPdf.h:241
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooAbsPdf.cxx:2412
virtual ~RooAbsPdf()
Destructor.
Definition: RooAbsPdf.cxx:284
double expectedEvents(const RooArgSet &nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.h:288
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:739
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:762
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:379
void setGeneratorConfig()
Remove the specialized numeric MC generator configuration associated with this object.
Definition: RooAbsPdf.cxx:3585
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:644
virtual RooFitResult * chi2FitTo(RooDataHist &data, const RooLinkedList &cmdList)
Calls RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList) and returns fit result.
Definition: RooAbsPdf.cxx:1794
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
Definition: RooAbsPdf.cxx:3339
Bool_t traceEvalPdf(Double_t value) const
Check that passed value is positive and not 'not-a-number'.
Definition: RooAbsPdf.cxx:457
const RooArgSet * getNormSet()
Normalization set with for above integral.
Definition: RooAbsPdf.h:366
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
Definition: RooAbsPdf.cxx:1066
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:1991
RooAbsReal * createScanCdf(const RooArgSet &iset, const RooArgSet &nset, Int_t numScanBins, Int_t intOrder)
Definition: RooAbsPdf.cxx:3481
TString _normRange
MC generator configuration specific for this object.
Definition: RooAbsPdf.h:403
void setNormRange(const char *rangeName)
Definition: RooAbsPdf.cxx:3616
const RooAbsReal * getNormIntegral(const RooArgSet &nset) const
Definition: RooAbsPdf.h:308
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:107
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:1542
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:365
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:2390
RooNumGenConfig * specialGeneratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
Definition: RooAbsPdf.cxx:3530
Bool_t canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:274
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:263
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:1964
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:437
const RooAbsReal * getCachedLastIntegral() const
Definition: RooAbsPdf.h:251
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:3167
RooAbsReal * getIntegral(RooArgSet const &set) const
Definition: RooAbsPdf.h:245
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:58
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:2000
Int_t _traceCount
Definition: RooAbsPdf.h:396
RooAbsReal * _norm
Definition: RooAbsPdf.h:364
int calcAsymptoticCorrectedCovariance(RooMinimizer &minimizer, RooAbsData const &data)
Use the asymptotically correct approach to estimate errors in the presence of weights.
Definition: RooAbsPdf.cxx:1270
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:1945
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:3498
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:551
virtual Double_t getValV(const RooArgSet *set=0) const
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:327
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:121
Int_t _errorCount
Definition: RooAbsPdf.h:392
Bool_t _selectComp
Definition: RooAbsPdf.h:399
virtual RooArgSet * getConstraints(const RooArgSet &, RooArgSet &, Bool_t) const
Definition: RooAbsPdf.h:210
@ CanBeExtended
Definition: RooAbsPdf.h:268
@ MustBeExtended
Definition: RooAbsPdf.h:268
@ CanNotBeExtended
Definition: RooAbsPdf.h:268
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:656
const char * normRange() const
Definition: RooAbsPdf.h:303
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:3404
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:516
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:682
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:2353
virtual Double_t expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.cxx:3319
Int_t _negCount
Definition: RooAbsPdf.h:397
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:1565
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:1840
Bool_t mustBeExtended() const
If true PDF must provide extended likelihood term.
Definition: RooAbsPdf.h:278
double extendedTerm(double observedEvents, const RooArgSet *nset=0) const
Return the extended likelihood term ( ) of this PDF for the given number of observed events.
Definition: RooAbsPdf.cxx:812
void setNormRangeOverride(const char *rangeName)
Definition: RooAbsPdf.cxx:3633
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of 'global observables' – for RooStats tools.
Definition: RooAbsPdf.cxx:2669
double normalizeWithNaNPacking(double rawVal, double normVal) const
Definition: RooAbsPdf.cxx:290
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:2425
const RooNumGenConfig * getGeneratorConfig() const
Return the numeric MC generator configuration used for this object.
Definition: RooAbsPdf.cxx:3558
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
Definition: RooAbsPdf.cxx:2400
void computeBatch(cudaStream_t *, double *output, size_t size, RooBatchCompute::DataMap &) const
Base function for computing multiple values of a RooAbsPdf.
Definition: RooAbsPdf.cxx:3655
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition: RooAbsPdf.h:272
RooAbsPdf()
Default constructor.
Definition: RooAbsPdf.cxx:226
virtual RooAbsGenContext * binnedGenContext(const RooArgSet &vars, Bool_t verbose=kFALSE) const
Return a binned generator context.
Definition: RooAbsPdf.cxx:1981
static RooNumGenConfig * defaultGeneratorConfig()
Returns the default numeric MC generator configuration for all RooAbsReals.
Definition: RooAbsPdf.cxx:3520
RooNumGenConfig * _specGeneratorConfig
Definition: RooAbsPdf.h:401
friend class RooAddGenContextOrig
Definition: RooAbsPdf.h:347
Double_t _rawValue
Definition: RooAbsPdf.h:363
static TString _normRangeOverride
Definition: RooAbsPdf.h:404
static Int_t _verboseEval
Definition: RooAbsPdf.h:358
virtual Bool_t redirectServersHook(const RooAbsCollection &, Bool_t, Bool_t, Bool_t)
Function that is called at the end of redirectServers().
Definition: RooAbsPdf.h:379
int calcSumW2CorrectedCovariance(RooMinimizer &minimizer, RooAbsReal const &nll) const
Apply correction to errors and covariance matrix.
Definition: RooAbsPdf.cxx:1350
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:3369
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:63
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 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...
Definition: RooAbsReal.cxx:312
Double_t getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:93
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:45
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:49
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:33
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.
std::map< DataKey, RooSpan< const double > > DataMap
Configuration struct for RooAbsPdf::minimizeNLL with all the default.
Definition: RooAbsPdf.h:166
const RooArgSet * minosSet
Definition: RooAbsPdf.h:183
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:31
static void output(int code)
Definition: gifencode.c:226