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