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 "RooGlobalFunc.h"
22#include "RooFit/UniqueId.h"
23
24class RooDataSet;
25class RooDataHist ;
26class RooArgSet ;
27class RooAbsGenContext ;
28class RooFitResult ;
29class RooExtendPdf ;
30class RooCategory ;
31class TPaveText;
32class TH1F;
33class TH2F;
34class TList ;
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=nullptr) ;
49 RooAbsPdf(const char *name, const char *title, double minVal, double maxVal) ;
50 // RooAbsPdf(const RooAbsPdf& other, const char* name=nullptr);
51 ~RooAbsPdf() override;
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] whatVars Set of observables to generate for each event according to this model.
58 /// \param[in] nEvents How many events to generate
59 /// \param arg1,arg2,arg3,arg4,arg5 Optional command arguments.
60 RooDataSet *generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg& arg1,
61 const RooCmdArg& arg2=RooCmdArg::none(), const RooCmdArg& arg3=RooCmdArg::none(),
62 const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none()) {
63 return generate(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5) ;
64 }
65 RooDataSet *generate(const RooArgSet &whatVars,
66 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
67 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
68 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) ;
69 RooDataSet *generate(const RooArgSet &whatVars, double nEvents = 0, bool verbose=false, bool autoBinned=true,
70 const char* binnedTag="", bool expectedData=false, bool extended = false) const;
71 RooDataSet *generate(const RooArgSet &whatVars, const RooDataSet &prototype, Int_t nEvents= 0,
72 bool verbose=false, bool randProtoOrder=false, bool resampleProto=false) const;
73
74
75 class GenSpec {
76 public:
77 virtual ~GenSpec() ;
78 GenSpec() { _genContext = nullptr ; _protoData = nullptr ; _init = false ; _extended=false, _nGen=0 ; _randProto = false ; _resampleProto=false ; }
79 private:
80 GenSpec(RooAbsGenContext* context, const RooArgSet& whatVars, RooDataSet* protoData, Int_t nGen, bool extended,
81 bool randProto, bool resampleProto, TString dsetName, bool init=false) ;
82 GenSpec(const GenSpec& other) ;
83
84 friend class RooAbsPdf ;
89 bool _extended ;
93 bool _init ;
94
95 ClassDef(GenSpec,0) // Generation specification
96 } ;
97
98 ///Prepare GenSpec configuration object for efficient generation of multiple datasets from identical specification.
99 GenSpec* prepareMultiGen(const RooArgSet &whatVars,
100 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
101 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
102 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) ;
103 ///Generate according to GenSpec obtained from prepareMultiGen().
104 RooDataSet* generate(GenSpec&) const ;
105
106
107 ////////////////////////////////////////////////////////////////////////////////
108 /// As RooAbsPdf::generateBinned(const RooArgSet&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&) const.
109 /// \param[in] whatVars set
110 /// \param[in] nEvents How many events to generate
111 /// \param arg1,arg2,arg3,arg4,arg5 ordered arguments
112 virtual RooDataHist *generateBinned(const RooArgSet &whatVars, double nEvents, const RooCmdArg& arg1,
113 const RooCmdArg& arg2=RooCmdArg::none(), const RooCmdArg& arg3=RooCmdArg::none(),
114 const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none()) const {
115 return generateBinned(whatVars,RooFit::NumEvents(nEvents),arg1,arg2,arg3,arg4,arg5);
116 }
117 virtual RooDataHist *generateBinned(const RooArgSet &whatVars,
118 const RooCmdArg& arg1=RooCmdArg::none(),const RooCmdArg& arg2=RooCmdArg::none(),
119 const RooCmdArg& arg3=RooCmdArg::none(),const RooCmdArg& arg4=RooCmdArg::none(),
120 const RooCmdArg& arg5=RooCmdArg::none(),const RooCmdArg& arg6=RooCmdArg::none()) const;
121 virtual RooDataHist *generateBinned(const RooArgSet &whatVars, double nEvents, bool expectedData=false, bool extended=false) const;
122
123 virtual RooDataSet* generateSimGlobal(const RooArgSet& whatVars, Int_t nEvents) ;
124
125 ///Helper calling plotOn(RooPlot*, RooLinkedList&) const
127 const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
128 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
129 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
130 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none(),
131 const RooCmdArg& arg9=RooCmdArg::none(), const RooCmdArg& arg10=RooCmdArg::none()
132 ) const override {
133 return RooAbsReal::plotOn(frame,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8,arg9,arg10) ;
134 }
135 RooPlot* plotOn(RooPlot* frame, RooLinkedList& cmdList) const override ;
136
137 /// Add a box with parameter values (and errors) to the specified frame
138 virtual RooPlot* paramOn(RooPlot* frame,
139 const RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
140 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
141 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
142 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
143
144 virtual RooPlot* paramOn(RooPlot* frame, const RooAbsData* data, const char *label= "", Int_t sigDigits = 2,
145 Option_t *options = "NELU", double xmin=0.65,
146 double xmax = 0.9, double ymax = 0.9) ;
147
148 // Built-in generator support
149 virtual Int_t getGenerator(const RooArgSet& directVars, RooArgSet &generateVars, bool staticInitOK=true) const;
150 virtual void initGenerator(Int_t code) ;
151 virtual void generateEvent(Int_t code);
152 virtual bool isDirectGenSafe(const RooAbsArg& arg) const ;
153
154 // Configuration of MC generators used for this pdf
155 const RooNumGenConfig* getGeneratorConfig() const ;
158 RooNumGenConfig* specialGeneratorConfig(bool createOnTheFly) ;
159 void setGeneratorConfig() ;
160 void setGeneratorConfig(const RooNumGenConfig& config) ;
161
162 // -log(L) fits to binned and unbinned data
163 virtual RooFitResult* fitTo(RooAbsData& data, const RooLinkedList& cmdList={}) ;
164 /// Takes an arbitrary number of RooCmdArg command options and calls
165 /// RooAbsPdf::fitTo(RooAbsData& data, const RooLinkedList& cmdList).
166 template <typename... Args>
167 RooFitResult* fitTo(RooAbsData& data, RooCmdArg const& arg1, Args const&... args)
168 {
169 return fitTo(data, *RooFit::Detail::createCmdList(&arg1, &args...));
170 }
171
172 /// Configuration struct for RooAbsPdf::minimizeNLL with all the default
173 //values that also should be taked as the default values for
174 //RooAbsPdf::fitTo.
176 double recoverFromNaN = 10.;
177 int optConst = 2;
178 int verbose = 0;
179 int doSave = 0;
180 int doTimer = 0;
181 int printLevel = 1;
182 int strat = 1;
183 int initHesse = 0;
184 int hesse = 1;
185 int minos = 0;
186 int numee = 10;
187 int doEEWall = 1;
188 int doWarn = 1;
189 int doSumW2 = -1;
190 int doAsymptotic = -1;
191 int maxCalls = -1;
192 int doOffset = -1;
193 int parallelize = 0;
196 bool timingAnalysis = false;
197 const RooArgSet* minosSet = nullptr;
198 std::string minType;
199 std::string minAlg = "minuit";
200 };
201 std::unique_ptr<RooFitResult> minimizeNLL(RooAbsReal & nll, RooAbsData const& data, MinimizerConfig const& cfg);
202
203 virtual RooAbsReal* createNLL(RooAbsData& data, const RooLinkedList& cmdList={}) ;
204 /// Takes an arbitrary number of RooCmdArg command options and calls
205 /// RooAbsPdf::createNLL(RooAbsData& data, const RooLinkedList& cmdList).
206 template <typename... Args>
207 RooAbsReal* createNLL(RooAbsData& data, RooCmdArg const& arg1, Args const&... args)
208 {
209 return createNLL(data, *RooFit::Detail::createCmdList(&arg1, &args...));
210 }
211
212 // Chi^2 fits to histograms
215 RooFitResult* chi2FitTo(RooDataHist& data, const RooLinkedList& cmdList) override ;
217 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
218 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) override ;
219
220 // Chi^2 fits to X-Y datasets
221 RooAbsReal* createChi2(RooDataSet& data, const RooLinkedList& cmdList) override ;
222
223
224 // Constraint management
225 virtual RooArgSet* getConstraints(const RooArgSet& /*observables*/, RooArgSet& /*constrainedParams*/,
226 bool /*stripDisconnected*/, bool /*removeConstraintsFromPdf*/=false) const
227 {
228 // Interface to retrieve constraint terms on this pdf. Default implementation returns null
229 return nullptr ;
230 }
231 RooArgSet* getAllConstraints(const RooArgSet& observables, RooArgSet& constrainedParams,
232 bool stripDisconnected=true, bool removeConstraintsFromPdf=false) const ;
233
234 // Project p.d.f into lower dimensional p.d.f
235 virtual RooAbsPdf* createProjection(const RooArgSet& iset) ;
236
237 // Create cumulative density function from p.d.f
238 RooAbsReal* createCdf(const RooArgSet& iset, const RooArgSet& nset=RooArgSet()) ;
239 RooAbsReal* createCdf(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2=RooCmdArg::none(),
240 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
241 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
242 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
243 RooAbsReal* createScanCdf(const RooArgSet& iset, const RooArgSet& nset, Int_t numScanBins, Int_t intOrder) ;
244
245 // Function evaluation support
246 double getValV(const RooArgSet* set=nullptr) const override ;
247 virtual double getLogVal(const RooArgSet* set=nullptr) const ;
248
249 RooSpan<const double> getLogValBatch(std::size_t begin, std::size_t batchSize,
250 const RooArgSet* normSet = nullptr) const;
251 RooSpan<const double> getLogProbabilities(RooBatchCompute::RunContext& evalData, const RooArgSet* normSet = nullptr) const;
252 void getLogProbabilities(RooSpan<const double> pdfValues, double * output) const;
253
254 /// \copydoc getNorm(const RooArgSet*) const
255 double getNorm(const RooArgSet& nset) const {
256 return getNorm(&nset) ;
257 }
258 virtual double getNorm(const RooArgSet* set=nullptr) const ;
259
260 virtual void resetErrorCounters(Int_t resetValue=10) ;
261 void setTraceCounter(Int_t value, bool allNodes=false) ;
262
263 double analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=nullptr) const override ;
264
265 /// Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
266 /// Always returns false, unless a PDF overrides this function.
267 virtual bool selfNormalized() const {
268 return false ;
269 }
270
271 // Support for extended maximum likelihood, switched off by default
273 /// Returns ability of PDF to provide extended likelihood terms. Possible
274 /// answers are in the enumerator RooAbsPdf::ExtendMode.
275 /// This default implementation always returns CanNotBeExtended.
276 virtual ExtendMode extendMode() const { return CanNotBeExtended; }
277 /// If true, PDF can provide extended likelihood term.
278 inline bool canBeExtended() const {
279 return (extendMode() != CanNotBeExtended) ;
280 }
281 /// If true PDF must provide extended likelihood term.
282 inline bool mustBeExtended() const {
283 return (extendMode() == MustBeExtended) ;
284 }
285 /// Return expected number of events to be used in calculation of extended
286 /// likelihood.
287 virtual double expectedEvents(const RooArgSet* nset) const ;
288 /// Return expected number of events to be used in calculation of extended
289 /// likelihood. This function should not be overridden, as it just redirects
290 /// to the actual virtual function but takes a RooArgSet reference instead of
291 /// pointer (\see expectedEvents(const RooArgSet*) const).
292 double expectedEvents(const RooArgSet& nset) const {
293 return expectedEvents(&nset) ;
294 }
295
296 // Printing interface (human readable)
297 void printValue(std::ostream& os) const override ;
298 void printMultiline(std::ostream& os, Int_t contents, bool verbose=false, TString indent="") const override ;
299
300 static void verboseEval(Int_t stat) ;
301 static int verboseEval() ;
302
303 double extendedTerm(double sumEntries, double expected, double sumEntriesW2=0.0, bool doOffset=false) const;
304 double extendedTerm(double sumEntries, RooArgSet const* nset, double sumEntriesW2=0.0, bool doOffset=false) const;
305 double extendedTerm(RooAbsData const& data, bool weightSquared, bool doOffset=false) const;
306
307 void setNormRange(const char* rangeName) ;
308 const char* normRange() const {
309 return _normRange.Length()>0 ? _normRange.Data() : nullptr ;
310 }
311 void setNormRangeOverride(const char* rangeName) ;
312
313 const RooAbsReal* getNormIntegral(const RooArgSet& nset) const { return getNormObj(nullptr,&nset,nullptr) ; }
314
315 virtual const RooAbsReal* getNormObj(const RooArgSet* set, const RooArgSet* iset, const TNamed* rangeName=nullptr) const ;
316
317 virtual RooAbsGenContext* binnedGenContext(const RooArgSet &vars, bool verbose= false) const ;
318
319 virtual RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr,
320 const RooArgSet* auxProto=nullptr, bool verbose= false) const ;
321
322 virtual RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=nullptr, const RooArgSet* auxProto=nullptr,
323 bool verbose=false, bool autoBinned=true, const char* binnedTag="") const ;
324
325 std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
326
327private:
328
329 RooDataSet *generate(RooAbsGenContext& context, const RooArgSet& whatVars, const RooDataSet* prototype,
330 double nEvents, bool verbose, bool randProtoOrder, bool resampleProto, bool skipInit=false,
331 bool extended=false) const ;
332
333 // Implementation version
334 virtual RooPlot* paramOn(RooPlot* frame, const RooArgSet& params, bool showConstants=false,
335 const char *label= "", Int_t sigDigits = 2, Option_t *options = "NELU", double xmin=0.65,
336 double xmax= 0.99,double ymax=0.95, const RooCmdArg* formatCmd=nullptr) ;
337
338 void logBatchComputationErrors(RooSpan<const double>& outputs, std::size_t begin) const;
339 bool traceEvalPdf(double value) const;
340
341 /// Setter for the _normSet member, which should never be set directly.
342 inline void setActiveNormSet(RooArgSet const* normSet) const {
343 _normSet = normSet;
344 // Also store the unique ID of the _normSet. This makes it possible to
345 // detect if the pointer was invalidated.
347 }
348
349protected:
350
351 /// Checks if `normSet` is the currently active normalization set of this
352 /// PDF, meaning is exactly the same object as the one the `_normSet` member
353 /// points to (or both are `nullptr`).
354 inline bool isActiveNormSet(RooArgSet const* normSet) const {
355 return RooFit::getUniqueId(normSet).value() == _normSetId;
356 }
357
358 double normalizeWithNaNPacking(double rawVal, double normVal) const;
359
360 RooPlot *plotOn(RooPlot *frame, PlotOpt o) const override;
361
362 friend class RooMCStudy ;
363
364 Int_t* randomizeProtoOrder(Int_t nProto,Int_t nGen,bool resample=false) const ;
365
366 // This also forces the definition of a copy ctor in derived classes
367 RooAbsPdf(const RooAbsPdf& other, const char* name = nullptr);
368
370
371 virtual bool syncNormalization(const RooArgSet* dset, bool adjustProxies=true) const ;
372
373 mutable double _rawValue ;
374 mutable RooAbsReal* _norm = nullptr; //! Normalization integral (owned by _normMgr)
375 mutable RooArgSet const* _normSet = nullptr; //! Normalization set with for above integral
376
378 public:
379 CacheElem(RooAbsReal& norm) : _norm(&norm) {} ;
380 ~CacheElem() override ;
383 } ;
384 mutable RooObjCacheManager _normMgr ; //! The cache manager
385
386 bool redirectServersHook(const RooAbsCollection & newServerList, bool mustReplaceAll,
387 bool nameChange, bool isRecursiveStep) override;
388
389
390 mutable Int_t _errorCount ; ///< Number of errors remaining to print
391 mutable Int_t _traceCount ; ///< Number of traces remaining to print
392 mutable Int_t _negCount ; ///< Number of negative probablities remaining to print
393
394 bool _selectComp ; ///< Component selection flag for RooAbsPdf::plotCompOn
395
396 std::unique_ptr<RooNumGenConfig> _specGeneratorConfig ; ///<! MC generator configuration specific for this object
397
398 TString _normRange ; ///< Normalization range
400
401private:
402 mutable RooFit::UniqueId<RooArgSet>::Value_t _normSetId = RooFit::UniqueId<RooArgSet>::nullval; ///<! Unique ID of the currently-active normalization set
403
405 int calcSumW2CorrectedCovariance(RooMinimizer& minimizer, RooAbsReal & nll) const;
406
407 ClassDefOverride(RooAbsPdf,5) // Abstract PDF with normalization support
408};
409
410
411
412
413#endif
const char Option_t
Definition: RtypesCore.h:66
#define ClassDef(name, id)
Definition: Rtypes.h:335
#define ClassDefOverride(name, id)
Definition: Rtypes.h:339
static void indent(ostringstream &buf, int indent_level)
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
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
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition: RooAbsArg.h:72
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:61
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
Normalization set with for above integral.
Definition: RooAbsPdf.h:377
~CacheElem() override
Destructor of normalization cache element.
Definition: RooAbsPdf.cxx:3287
CacheElem(RooAbsReal &norm)
Definition: RooAbsPdf.h:379
RooAbsReal * _norm
Definition: RooAbsPdf.h:382
RooArgList containedArgs(Action) override
Definition: RooAbsPdf.h:381
RooArgSet _whatVars
Definition: RooAbsPdf.h:86
RooAbsGenContext * _genContext
Definition: RooAbsPdf.h:85
GenSpec(const GenSpec &other)
TString _dsetName
Definition: RooAbsPdf.h:92
RooDataSet * _protoData
Definition: RooAbsPdf.h:87
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:2112
virtual bool syncNormalization(const RooArgSet *dset, bool adjustProxies=true) const
Verify that the normalization integral cached with this PDF is valid for given set of normalization o...
Definition: RooAbsPdf.cxx:551
int calcSumW2CorrectedCovariance(RooMinimizer &minimizer, RooAbsReal &nll) const
Apply correction to errors and covariance matrix.
Definition: RooAbsPdf.cxx:1319
double getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition: RooAbsPdf.h:255
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
Definition: RooAbsPdf.cxx:3640
RooSpan< const double > getLogValBatch(std::size_t begin, std::size_t batchSize, const RooArgSet *normSet=nullptr) const
RooObjCacheManager _normMgr
Definition: RooAbsPdf.h:384
RooFit::UniqueId< RooArgSet >::Value_t _normSetId
! Unique ID of the currently-active normalization set
Definition: RooAbsPdf.h:402
std::unique_ptr< RooNumGenConfig > _specGeneratorConfig
! MC generator configuration specific for this object
Definition: RooAbsPdf.h:396
double getValV(const RooArgSet *set=nullptr) const override
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:392
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()) override
Create a from a histogram and this function.
Definition: RooAbsPdf.cxx:1789
bool _selectComp
Component selection flag for RooAbsPdf::plotCompOn.
Definition: RooAbsPdf.h:394
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooAbsPdf.cxx:2359
double expectedEvents(const RooArgSet &nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.h:292
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:716
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:739
void setGeneratorConfig()
Remove the specialized numeric MC generator configuration associated with this object.
Definition: RooAbsPdf.cxx:3520
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:648
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
Definition: RooAbsPdf.cxx:3275
bool isActiveNormSet(RooArgSet const *normSet) const
Checks if normSet is the currently active normalization set of this PDF, meaning is exactly the same ...
Definition: RooAbsPdf.h:354
virtual double expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.cxx:3255
virtual RooAbsGenContext * binnedGenContext(const RooArgSet &vars, bool verbose=false) const
Return a binned generator context.
Definition: RooAbsPdf.cxx:1930
RooAbsReal * createScanCdf(const RooArgSet &iset, const RooArgSet &nset, Int_t numScanBins, Int_t intOrder)
Definition: RooAbsPdf.cxx:3417
TString _normRange
Normalization range.
Definition: RooAbsPdf.h:398
virtual bool isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
Definition: RooAbsPdf.cxx:2372
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList={})
Construct representation of -log(L) of PDF with given dataset.
Definition: RooAbsPdf.cxx:966
RooAbsReal * createNLL(RooAbsData &data, RooCmdArg const &arg1, Args const &... args)
Takes an arbitrary number of RooCmdArg command options and calls RooAbsPdf::createNLL(RooAbsData& dat...
Definition: RooAbsPdf.h:207
Int_t * randomizeProtoOrder(Int_t nProto, Int_t nGen, bool resample=false) const
Return lookup table with randomized order for nProto prototype events.
Definition: RooAbsPdf.cxx:2300
virtual RooFitResult * fitTo(RooAbsData &data, const RooLinkedList &cmdList={})
Fit PDF to given dataset.
Definition: RooAbsPdf.cxx:1611
void setNormRange(const char *rangeName)
Definition: RooAbsPdf.cxx:3578
const RooAbsReal * getNormIntegral(const RooArgSet &nset) const
Definition: RooAbsPdf.h:313
~RooAbsPdf() override
Destructor.
Definition: RooAbsPdf.cxx:350
bool mustBeExtended() const
If true PDF must provide extended likelihood term.
Definition: RooAbsPdf.h:282
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:375
RooFitResult * chi2FitTo(RooDataHist &data, const RooLinkedList &cmdList) override
Calls RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList) and returns fit result.
Definition: RooAbsPdf.cxx:1765
RooNumGenConfig * specialGeneratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
Definition: RooAbsPdf.cxx:3468
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:3104
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:60
virtual bool selfNormalized() const
Shows if a PDF is self-normalized, which means that no attempt is made to add a normalization term.
Definition: RooAbsPdf.h:267
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print multi line detailed information of this RooAbsPdf.
Definition: RooAbsPdf.cxx:1913
virtual RooDataHist * generateBinned(const RooArgSet &whatVars, double 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:112
Int_t _traceCount
Number of traces remaining to print.
Definition: RooAbsPdf.h:391
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:278
RooAbsReal * _norm
Definition: RooAbsPdf.h:374
int calcAsymptoticCorrectedCovariance(RooMinimizer &minimizer, RooAbsData const &data)
Use the asymptotically correct approach to estimate errors in the presence of weights.
Definition: RooAbsPdf.cxx:1239
void setTraceCounter(Int_t value, bool allNodes=false)
Reset trace counter to given value, limiting the number of future trace messages for this pdf to 'val...
Definition: RooAbsPdf.cxx:660
Int_t _errorCount
Number of errors remaining to print.
Definition: RooAbsPdf.h:390
@ CanBeExtended
Definition: RooAbsPdf.h:272
@ MustBeExtended
Definition: RooAbsPdf.h:272
@ CanNotBeExtended
Definition: RooAbsPdf.h:272
double _rawValue
Definition: RooAbsPdf.h:373
const char * normRange() const
Definition: RooAbsPdf.h:308
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:3340
Int_t _negCount
Number of negative probablities remaining to print.
Definition: RooAbsPdf.h:392
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:1372
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 ...
Definition: RooAbsPdf.cxx:516
RooFitResult * fitTo(RooAbsData &data, RooCmdArg const &arg1, Args const &... args)
Takes an arbitrary number of RooCmdArg command options and calls RooAbsPdf::fitTo(RooAbsData& data,...
Definition: RooAbsPdf.h:167
void setActiveNormSet(RooArgSet const *normSet) const
Setter for the _normSet member, which should never be set directly.
Definition: RooAbsPdf.h:342
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Analytical integral with normalization (see RooAbsReal::analyticalIntegralWN() for further informatio...
Definition: RooAbsPdf.cxx:437
void setNormRangeOverride(const char *rangeName)
Definition: RooAbsPdf.cxx:3598
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of 'global observables' – for RooStats tools.
Definition: RooAbsPdf.cxx:2615
double normalizeWithNaNPacking(double rawVal, double normVal) const
Definition: RooAbsPdf.cxx:355
virtual RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false, bool autoBinned=true, const char *binnedTag="") const
Definition: RooAbsPdf.cxx:1949
RooArgSet * getAllConstraints(const RooArgSet &observables, RooArgSet &constrainedParams, bool stripDisconnected=true, bool removeConstraintsFromPdf=false) const
This helper function finds and collects all constraints terms of all component p.d....
Definition: RooAbsPdf.cxx:3434
const RooNumGenConfig * getGeneratorConfig() const
Return the numeric MC generator configuration used for this object.
Definition: RooAbsPdf.cxx:3496
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
Definition: RooAbsPdf.cxx:2347
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition: RooAbsPdf.h:276
RooAbsPdf()
Default constructor.
Definition: RooAbsPdf.cxx:294
virtual RooArgSet * getConstraints(const RooArgSet &, RooArgSet &, bool, bool=false) const
Definition: RooAbsPdf.h:225
bool traceEvalPdf(double value) const
Check that passed value is positive and not 'not-a-number'.
Definition: RooAbsPdf.cxx:457
static RooNumGenConfig * defaultGeneratorConfig()
Returns the default numeric MC generator configuration for all RooAbsReals.
Definition: RooAbsPdf.cxx:3458
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursiveStep) override
The cache manager.
Definition: RooAbsPdf.cxx:3620
void printValue(std::ostream &os) const override
Print value of p.d.f, also print normalization integral that was last used, if any.
Definition: RooAbsPdf.cxx:1894
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false) const
Interface function to create a generator context from a p.d.f.
Definition: RooAbsPdf.cxx:1940
static TString _normRangeOverride
Definition: RooAbsPdf.h:399
static Int_t _verboseEval
Definition: RooAbsPdf.h:369
double extendedTerm(double sumEntries, double expected, double sumEntriesW2=0.0, bool doOffset=false) const
Definition: RooAbsPdf.cxx:808
virtual Int_t getGenerator(const RooArgSet &directVars, RooArgSet &generateVars, bool staticInitOK=true) const
Load generatedVars with the subset of directVars that we can generate events for, and return a code t...
Definition: RooAbsPdf.cxx:2337
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:3305
virtual double getLogVal(const RooArgSet *set=nullptr) const
Return the log of the current value with given normalization An error message is printed if the argum...
Definition: RooAbsPdf.cxx:682
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 override
Helper calling plotOn(RooPlot*, RooLinkedList&) const.
Definition: RooAbsPdf.h:126
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:60
virtual RooAbsReal * createChi2(RooDataHist &data, const RooLinkedList &cmdList)
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 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.
virtual double maxVal(Int_t code) const
Return maximum value for set of observables identified by code assigned in getMaxVal.
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:56
RooCategory is an object to represent discrete states.
Definition: RooCategory.h:28
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition: RooCmdArg.h:26
static const RooCmdArg & none()
Return reference to null argument.
Definition: RooCmdArg.cxx:51
The RooDataHist is a container class to hold N-dimensional binned data.
Definition: RooDataHist.h:39
RooDataSet is a container class to hold unbinned data.
Definition: RooDataSet.h:55
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:43
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:43
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
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:577
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:257
A doubly linked list.
Definition: TList.h:38
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.
void init()
Inspect hardware capabilities, and load the optimal library for RooFit computations.
std::unique_ptr< RooLinkedList > createCmdList(Args &&... args)
UniqueId< Class > const & getUniqueId(Class const *ptr)
A helper function to replace pointer comparisons with UniqueId comparisons.
Definition: UniqueId.h:87
Configuration struct for RooAbsPdf::minimizeNLL with all the default.
Definition: RooAbsPdf.h:175
const RooArgSet * minosSet
Definition: RooAbsPdf.h:197
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:32
A UniqueId can be added as a class member to enhance any class with a unique identifier for each inst...
Definition: UniqueId.h:39
unsigned long Value_t
Definition: UniqueId.h:41
static void output()