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 minVal, double maxVal) ;
50 // RooAbsPdf(const RooAbsPdf& other, const char* name=0);
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 = 0 ; _protoData = 0 ; _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&)
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 RooCmdArg& arg1=RooCmdArg::none(), const RooCmdArg& arg2=RooCmdArg::none(),
164 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
165 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
166 virtual RooFitResult* fitTo(RooAbsData& data, const RooLinkedList& cmdList) ;
167
168 /// Configuration struct for RooAbsPdf::minimizeNLL with all the default
169 //values that also should be taked as the default values for
170 //RooAbsPdf::fitTo.
172 double recoverFromNaN = 10.;
173 int optConst = 2;
174 int verbose = 0;
175 int doSave = 0;
176 int doTimer = 0;
177 int printLevel = 1;
178 int strat = 1;
179 int initHesse = 0;
180 int hesse = 1;
181 int minos = 0;
182 int numee = 10;
183 int doEEWall = 1;
184 int doWarn = 1;
185 int doSumW2 = -1;
186 int doAsymptotic = -1;
187 const RooArgSet* minosSet = nullptr;
188 std::string minType;
189 std::string minAlg = "minuit";
190 };
191 std::unique_ptr<RooFitResult> minimizeNLL(RooAbsReal & nll, RooAbsData const& data, MinimizerConfig const& cfg);
192
193 virtual RooAbsReal* createNLL(RooAbsData& data, const RooLinkedList& cmdList) ;
195 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
196 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
197
198 // Chi^2 fits to histograms
201 RooFitResult* chi2FitTo(RooDataHist& data, const RooLinkedList& cmdList) override ;
203 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(), const RooCmdArg& arg5=RooCmdArg::none(),
204 const RooCmdArg& arg6=RooCmdArg::none(), const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) override ;
205
206 // Chi^2 fits to X-Y datasets
207 RooAbsReal* createChi2(RooDataSet& data, const RooLinkedList& cmdList) override ;
208
209
210 // Constraint management
211 virtual RooArgSet* getConstraints(const RooArgSet& /*observables*/, RooArgSet& /*constrainedParams*/, bool /*stripDisconnected*/) const {
212 // Interface to retrieve constraint terms on this pdf. Default implementation returns null
213 return 0 ;
214 }
215 virtual RooArgSet* getAllConstraints(const RooArgSet& observables, RooArgSet& constrainedParams, bool stripDisconnected=true) const ;
216
217 // Project p.d.f into lower dimensional p.d.f
218 virtual RooAbsPdf* createProjection(const RooArgSet& iset) ;
219
220 // Create cumulative density function from p.d.f
221 RooAbsReal* createCdf(const RooArgSet& iset, const RooArgSet& nset=RooArgSet()) ;
222 RooAbsReal* createCdf(const RooArgSet& iset, const RooCmdArg& arg1, const RooCmdArg& arg2=RooCmdArg::none(),
223 const RooCmdArg& arg3=RooCmdArg::none(), const RooCmdArg& arg4=RooCmdArg::none(),
224 const RooCmdArg& arg5=RooCmdArg::none(), const RooCmdArg& arg6=RooCmdArg::none(),
225 const RooCmdArg& arg7=RooCmdArg::none(), const RooCmdArg& arg8=RooCmdArg::none()) ;
226 RooAbsReal* createScanCdf(const RooArgSet& iset, const RooArgSet& nset, Int_t numScanBins, Int_t intOrder) ;
227
228 // Function evaluation support
229 double getValV(const RooArgSet* set=0) const override ;
230 virtual double getLogVal(const RooArgSet* set=0) const ;
231
232 RooSpan<const double> getValues(RooBatchCompute::RunContext& evalData, const RooArgSet* normSet) const override;
234 RooSpan<const double> getLogValBatch(std::size_t begin, std::size_t batchSize,
235 const RooArgSet* normSet = nullptr) const;
236 RooSpan<const double> getLogProbabilities(RooBatchCompute::RunContext& evalData, const RooArgSet* normSet = nullptr) const;
237 void getLogProbabilities(RooSpan<const double> pdfValues, double * output) const;
238
239 /// \copydoc getNorm(const RooArgSet*) const
240 double getNorm(const RooArgSet& nset) const {
241 return getNorm(&nset) ;
242 }
243 virtual double getNorm(const RooArgSet* set=0) const ;
244
245 virtual void resetErrorCounters(Int_t resetValue=10) ;
246 void setTraceCounter(Int_t value, bool allNodes=false) ;
247
248 double analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=0) const override ;
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 selfNormalized() const {
253 return false ;
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 canBeExtended() const {
264 return (extendMode() != CanNotBeExtended) ;
265 }
266 /// If true PDF must provide extended likelihood term.
267 inline bool 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 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 void printValue(std::ostream& os) const override ;
283 void printMultiline(std::ostream& os, Int_t contents, bool verbose=false, TString indent="") const override ;
284
285 static void verboseEval(Int_t stat) ;
286 static int verboseEval() ;
287
288 double extendedTerm(double sumEntries, double expected, double sumEntriesW2=0.0) const;
289 double extendedTerm(double sumEntries, RooArgSet const* nset, double sumEntriesW2=0.0) const;
290 double extendedTerm(RooAbsData const& data, bool weightSquared) const;
291
292 void setNormRange(const char* rangeName) ;
293 const char* normRange() const {
294 return _normRange.Length()>0 ? _normRange.Data() : 0 ;
295 }
296 void setNormRangeOverride(const char* rangeName) ;
297
298 const RooAbsReal* getNormIntegral(const RooArgSet& nset) const { return getNormObj(0,&nset,0) ; }
299
300 virtual const RooAbsReal* getNormObj(const RooArgSet* set, const RooArgSet* iset, const TNamed* rangeName=0) const ;
301
302 virtual RooAbsGenContext* binnedGenContext(const RooArgSet &vars, bool verbose= false) const ;
303
304 virtual RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=0,
305 const RooArgSet* auxProto=0, bool verbose= false) const ;
306
307 virtual RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=0, const RooArgSet* auxProto=0,
308 bool verbose=false, bool autoBinned=true, const char* binnedTag="") const ;
309
310private:
311
312 RooDataSet *generate(RooAbsGenContext& context, const RooArgSet& whatVars, const RooDataSet* prototype,
313 double nEvents, bool verbose, bool randProtoOrder, bool resampleProto, bool skipInit=false,
314 bool extended=false) const ;
315
316 // Implementation version
317 virtual RooPlot* paramOn(RooPlot* frame, const RooArgSet& params, bool showConstants=false,
318 const char *label= "", Int_t sigDigits = 2, Option_t *options = "NELU", double xmin=0.65,
319 double xmax= 0.99,double ymax=0.95, const RooCmdArg* formatCmd=0) ;
320
321 void logBatchComputationErrors(RooSpan<const double>& outputs, std::size_t begin) const;
322 bool traceEvalPdf(double value) const;
323
324
325protected:
326 double normalizeWithNaNPacking(double rawVal, double normVal) const;
327
328 RooPlot *plotOn(RooPlot *frame, PlotOpt o) const override;
329
330 friend class RooEffGenContext ;
331 friend class RooAddGenContext ;
332 friend class RooProdGenContext ;
333 friend class RooSimGenContext ;
335 friend class RooConvGenContext ;
336 friend class RooSimultaneous ;
337 friend class RooAddGenContextOrig ;
338 friend class RooProdPdf ;
339 friend class RooMCStudy ;
340
341 Int_t* randomizeProtoOrder(Int_t nProto,Int_t nGen,bool resample=false) const ;
342
343 friend class RooExtendPdf ;
344 // This also forces the definition of a copy ctor in derived classes
345 RooAbsPdf(const RooAbsPdf& other, const char* name = 0);
346
347 friend class RooRealIntegral ;
349
350 virtual bool syncNormalization(const RooArgSet* dset, bool adjustProxies=true) const ;
351
352 friend class RooAbsAnaConvPdf ;
353 mutable double _rawValue ;
354 mutable RooAbsReal* _norm = nullptr; //! Normalization integral (owned by _normMgr)
355 mutable RooArgSet const* _normSet = nullptr; //! Normalization set with for above integral
356 inline const RooArgSet* getNormSet() { return _normSet; }
357
359 public:
360 CacheElem(RooAbsReal& norm) : _norm(&norm) {} ;
361 ~CacheElem() override ;
364 } ;
365 mutable RooObjCacheManager _normMgr ; //! The cache manager
366
367 friend class CacheElem ; // Cache needs to be able to clear _norm pointer
368
369 bool redirectServersHook(const RooAbsCollection&, bool, bool, bool) override {
370 // Hook function intercepting redirectServer calls. Discard current normalization
371 // object if any server is redirected
372
373 // Object is own by _normCacheManager that will delete object as soon as cache
374 // is sterilized by server redirect
375 _norm = nullptr ;
376
377 // Similar to the situation with the normalization integral above: if a
378 // server is redirected, the cached normalization set might not point to
379 // the right observables anymore. We need to reset it.
380 _normSet = nullptr ;
381 return false ;
382 } ;
383
384
385 mutable Int_t _errorCount ; ///< Number of errors remaining to print
386 mutable Int_t _traceCount ; ///< Number of traces remaining to print
387 mutable Int_t _negCount ; ///< Number of negative probablities remaining to print
388
389 bool _selectComp ; ///< Component selection flag for RooAbsPdf::plotCompOn
390
391 RooNumGenConfig* _specGeneratorConfig ; ///<! MC generator configuration specific for this object
392
393 TString _normRange ; ///< Normalization range
395
396private:
398 int calcSumW2CorrectedCovariance(RooMinimizer& minimizer, RooAbsReal & nll) const;
399
400 ClassDefOverride(RooAbsPdf,5) // Abstract PDF with normalization support
401};
402
403
404
405
406#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
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:645
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:62
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
~CacheElem() override
Destructor of normalization cache element.
Definition: RooAbsPdf.cxx:3247
CacheElem(RooAbsReal &norm)
Definition: RooAbsPdf.h:360
RooAbsReal * _norm
Definition: RooAbsPdf.h:363
RooArgList containedArgs(Action) override
Definition: RooAbsPdf.h:362
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:2059
virtual double 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:667
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:540
int calcSumW2CorrectedCovariance(RooMinimizer &minimizer, RooAbsReal &nll) const
Apply correction to errors and covariance matrix.
Definition: RooAbsPdf.cxx:1282
double getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition: RooAbsPdf.h:240
RooSpan< const double > getLogValBatch(std::size_t begin, std::size_t batchSize, const RooArgSet *normSet=nullptr) const
RooObjCacheManager _normMgr
Definition: RooAbsPdf.h:365
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:1736
bool _selectComp
Component selection flag for RooAbsPdf::plotCompOn.
Definition: RooAbsPdf.h:389
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
Definition: RooAbsPdf.cxx:2308
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:701
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:724
void setGeneratorConfig()
Remove the specialized numeric MC generator configuration associated with this object.
Definition: RooAbsPdf.cxx:3481
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:633
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
Definition: RooAbsPdf.cxx:3235
const RooArgSet * getNormSet()
Normalization set with for above integral.
Definition: RooAbsPdf.h:356
virtual double expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
Definition: RooAbsPdf.cxx:3215
virtual RooAbsReal * createNLL(RooAbsData &data, const RooLinkedList &cmdList)
Construct representation of -log(L) of PDFwith given dataset.
Definition: RooAbsPdf.cxx:996
virtual RooAbsGenContext * binnedGenContext(const RooArgSet &vars, bool verbose=false) const
Return a binned generator context.
Definition: RooAbsPdf.cxx:1877
RooAbsReal * createScanCdf(const RooArgSet &iset, const RooArgSet &nset, Int_t numScanBins, Int_t intOrder)
Definition: RooAbsPdf.cxx:3377
TString _normRange
Normalization range.
Definition: RooAbsPdf.h:393
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:2321
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:2249
void setNormRange(const char *rangeName)
Definition: RooAbsPdf.cxx:3512
const RooAbsReal * getNormIntegral(const RooArgSet &nset) const
Definition: RooAbsPdf.h:298
~RooAbsPdf() override
Destructor.
Definition: RooAbsPdf.cxx:312
bool mustBeExtended() const
If true PDF must provide extended likelihood term.
Definition: RooAbsPdf.h:267
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:1458
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition: RooAbsPdf.h:355
RooFitResult * chi2FitTo(RooDataHist &data, const RooLinkedList &cmdList) override
Calls RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList) and returns fit result.
Definition: RooAbsPdf.cxx:1690
RooNumGenConfig * specialGeneratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
Definition: RooAbsPdf.cxx:3426
virtual RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, bool verbose=false, bool autoBinned=true, const char *binnedTag="") const
Definition: RooAbsPdf.cxx:1896
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:3063
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=0) const override
Analytical integral with normalization (see RooAbsReal::analyticalIntegralWN() for further informatio...
Definition: RooAbsPdf.cxx:426
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:252
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:1860
RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet) const override
Compute batch of values for given input data, and normalise by integrating over the observables in no...
Definition: RooAbsPdf.cxx:407
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:386
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition: RooAbsPdf.h:263
RooAbsReal * _norm
Definition: RooAbsPdf.h:354
int calcAsymptoticCorrectedCovariance(RooMinimizer &minimizer, RooAbsData const &data)
Use the asymptotically correct approach to estimate errors in the presence of weights.
Definition: RooAbsPdf.cxx:1202
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:645
Int_t _errorCount
Number of errors remaining to print.
Definition: RooAbsPdf.h:385
@ CanBeExtended
Definition: RooAbsPdf.h:257
@ MustBeExtended
Definition: RooAbsPdf.h:257
@ CanNotBeExtended
Definition: RooAbsPdf.h:257
double _rawValue
Definition: RooAbsPdf.h:353
const char * normRange() const
Definition: RooAbsPdf.h:293
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:3300
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:505
Int_t _negCount
Number of negative probablities remaining to print.
Definition: RooAbsPdf.h:387
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:1481
void setNormRangeOverride(const char *rangeName)
Definition: RooAbsPdf.cxx:3529
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of 'global observables' – for RooStats tools.
Definition: RooAbsPdf.cxx:2565
double normalizeWithNaNPacking(double rawVal, double normVal) const
Definition: RooAbsPdf.cxx:318
double getValV(const RooArgSet *set=0) const override
Return current value, normalized by integrating over the observables in nset.
Definition: RooAbsPdf.cxx:355
virtual RooArgSet * getAllConstraints(const RooArgSet &observables, RooArgSet &constrainedParams, bool stripDisconnected=true) const
This helper function finds and collects all constraints terms of all component p.d....
Definition: RooAbsPdf.cxx:3394
const RooNumGenConfig * getGeneratorConfig() const
Return the numeric MC generator configuration used for this object.
Definition: RooAbsPdf.cxx:3454
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
Definition: RooAbsPdf.cxx:2296
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition: RooAbsPdf.h:261
RooAbsPdf()
Default constructor.
Definition: RooAbsPdf.cxx:254
bool traceEvalPdf(double value) const
Check that passed value is positive and not 'not-a-number'.
Definition: RooAbsPdf.cxx:446
static RooNumGenConfig * defaultGeneratorConfig()
Returns the default numeric MC generator configuration for all RooAbsReals.
Definition: RooAbsPdf.cxx:3416
RooNumGenConfig * _specGeneratorConfig
! MC generator configuration specific for this object
Definition: RooAbsPdf.h:391
friend class RooAddGenContextOrig
Definition: RooAbsPdf.h:337
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:1841
virtual RooArgSet * getConstraints(const RooArgSet &, RooArgSet &, bool) const
Definition: RooAbsPdf.h:211
static TString _normRangeOverride
Definition: RooAbsPdf.h:394
static Int_t _verboseEval
Definition: RooAbsPdf.h:348
virtual RooAbsGenContext * genContext(const RooArgSet &vars, const RooDataSet *prototype=0, const RooArgSet *auxProto=0, bool verbose=false) const
Interface function to create a generator context from a p.d.f.
Definition: RooAbsPdf.cxx:1887
bool redirectServersHook(const RooAbsCollection &, bool, bool, bool) override
Function that is called at the end of redirectServers().
Definition: RooAbsPdf.h:369
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:2286
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:3265
double extendedTerm(double sumEntries, double expected, double sumEntriesW2=0.0) const
Definition: RooAbsPdf.cxx:786
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:64
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 RooSpan< const double > getValues(RooBatchCompute::RunContext &evalData, const RooArgSet *normSet=nullptr) const
Definition: RooAbsReal.cxx:311
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.
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:57
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
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:55
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: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
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:574
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.
Configuration struct for RooAbsPdf::minimizeNLL with all the default.
Definition: RooAbsPdf.h:171
const RooArgSet * minosSet
Definition: RooAbsPdf.h:187
This struct enables passing computation data around between elements of a computation graph.
Definition: RunContext.h:32
static void output()