Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
RooAbsPdf.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitCore *
4 * @(#)root/roofitcore:$Id$
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
17//////////////////////////////////////////////////////////////////////////////
18/** \class RooAbsPdf
19 \ingroup Roofitcore
20
21## RooAbsPdf, the base class of all PDFs
22
23RooAbsPdf is the abstract interface for all probability density
24functions. The class provides hybrid analytical/numerical
25normalization for its implementations, error tracing, and a Monte Carlo
26generator interface.
27
28### A Minimal PDF Implementation
29
30A minimal implementation of a PDF class derived from RooAbsPdf
31should override the `evaluate()` function. This function should
32return the PDF's value (which does not need to be normalised).
33
34
35#### Normalization/Integration
36
37Although the normalization of a PDF is an integral part of a
38probability density function, normalization is treated separately
39in RooAbsPdf. The reason is that a RooAbsPdf object is more than a
40PDF: it can be a building block for a more complex composite PDF
41if any of its variables are functions instead of variables. In
42such cases, the normalization of the composite PDF may not simply be
43integral over the dependents of the top-level PDF: these are
44functions with potentially non-trivial Jacobian terms themselves.
45\note Therefore, no explicit attempt should be made to normalize the
46function output in evaluate(). In particular, normalisation constants
47can be omitted to speed up the function evaluations, and included later
48in the integration of the PDF (see below), which is rarely called in
49comparison to the `evaluate()` function.
50
51In addition, RooAbsPdf objects do not have a static concept of what
52variables are parameters, and what variables are dependents (which
53need to be integrated over for a correct PDF normalization).
54Instead, the choice of normalization is always specified each time a
55normalized value is requested from the PDF via the getVal()
56method.
57
58RooAbsPdf manages the entire normalization logic of each PDF with
59the help of a RooRealIntegral object, which coordinates the integration
60of a given choice of normalization. By default, RooRealIntegral will
61perform an entirely numeric integration of all dependents. However,
62PDFs can advertise one or more (partial) analytical integrals of
63their function, and these will be used by RooRealIntegral, if it
64determines that this is safe (i.e., no hidden Jacobian terms,
65multiplication with other PDFs that have one or more dependents in
66commen etc).
67
68#### Implementing analytical integrals
69To implement analytical integrals, two functions must be implemented. First,
70
71```
72Int_t getAnalyticalIntegral(const RooArgSet& integSet, RooArgSet& anaIntSet)
73```
74should return the analytical integrals that are supported. `integSet`
75is the set of dependents for which integration is requested. The
76function should copy the subset of dependents it can analytically
77integrate to `anaIntSet`, and return a unique identification code for
78this integration configuration. If no integration can be
79performed, zero should be returned. Second,
80
81```
82double analyticalIntegral(Int_t code)
83```
84
85implements the actual analytical integral(s) advertised by
86`getAnalyticalIntegral()`. This function will only be called with
87codes returned by `getAnalyticalIntegral()`, except code zero.
88
89The integration range for each dependent to be integrated can
90be obtained from the dependent's proxy functions `min()` and
91`max()`. Never call these proxy functions for any proxy not known to
92be a dependent via the integration code. Doing so may be
93ill-defined, e.g., in case the proxy holds a function, and will
94trigger an assert. Integrated category dependents should always be
95summed over all of their states.
96
97
98
99### Direct generation of observables
100
101Distributions for any PDF can be generated with the accept/reject method,
102but for certain PDFs, more efficient methods may be implemented. To
103implement direct generation of one or more observables, two
104functions need to be implemented, similar to those for analytical
105integrals:
106
107```
108Int_t getGenerator(const RooArgSet& generateVars, RooArgSet& directVars)
109```
110and
111```
112void generateEvent(Int_t code)
113```
114
115The first function advertises observables, for which distributions can be generated,
116similar to the way analytical integrals are advertised. The second
117function implements the actual generator for the advertised observables.
118
119The generated dependent values should be stored in the proxy
120objects. For this, the assignment operator can be used (i.e. `xProxy
121= 3.0` ). Never call assign to any proxy not known to be a dependent
122via the generation code. Doing so may be ill-defined, e.g. in case
123the proxy holds a function, and will trigger an assert.
124
125
126### Batched function evaluations (Advanced usage)
127
128To speed up computations with large numbers of data events in unbinned fits,
129it is beneficial to override `computeBatch()`. Like this, large spans of
130computations can be done, without having to call `evaluate()` for each single data event.
131`computeBatch()` should execute the same computation as `evaluate()`, but it
132may choose an implementation that is capable of SIMD computations.
133If computeBatch is not implemented, the classic and slower `evaluate()` will be
134called for each data event.
136
137#include "RooAbsPdf.h"
138
139#include "RooNormalizedPdf.h"
140#include "RooMsgService.h"
141#include "RooDataSet.h"
142#include "RooArgSet.h"
143#include "RooArgProxy.h"
144#include "RooRealProxy.h"
145#include "RooRealVar.h"
146#include "RooGenContext.h"
147#include "RooBinnedGenContext.h"
148#include "RooPlot.h"
149#include "RooCurve.h"
150#include "RooNLLVar.h"
151#include "RooCategory.h"
152#include "RooNameReg.h"
153#include "RooCmdConfig.h"
154#include "RooGlobalFunc.h"
155#include "RooAddition.h"
156#include "RooRandom.h"
157#include "RooNumIntConfig.h"
158#include "RooProjectedPdf.h"
159#include "RooCustomizer.h"
160#include "RooParamBinning.h"
161#include "RooNumCdf.h"
162#include "RooFitResult.h"
163#include "RooNumGenConfig.h"
164#include "RooCachedReal.h"
165#include "RooXYChi2Var.h"
166#include "RooChi2Var.h"
167#include "RooMinimizer.h"
168#include "RooRealIntegral.h"
169#include "RooWorkspace.h"
170#include "RooNaNPacker.h"
171#include "RooHelpers.h"
172#include "RooFormulaVar.h"
173#include "RooDerivative.h"
175#include "RooVDTHeaders.h"
178#include "RunContext.h"
179#include "ConstraintHelpers.h"
180#include "RooFitDriver.h"
181#include "RooSimultaneous.h"
182
183#include "ROOT/StringUtils.hxx"
184#include "TMath.h"
185#include "TPaveText.h"
186#include "TMatrixD.h"
187#include "TMatrixDSym.h"
188#include "Math/CholeskyDecomp.h"
189
190#include <algorithm>
191#include <iostream>
192#include <string>
193#include <cmath>
194#include <stdexcept>
195
196namespace {
197
198bool interpretExtendedCmdArg(RooAbsPdf const& pdf, int extendedCmdArg) {
199 // Process automatic extended option
200 if (extendedCmdArg == 2) {
202 if (ext) {
203 oocoutI(&pdf, Minimization)
204 << "p.d.f. provides expected number of events, including extended term in likelihood." << std::endl;
205 }
206 return ext;
207 }
208 return extendedCmdArg;
209}
210
211inline double getLog(double prob, RooAbsReal const *caller)
212{
213
214 if (std::abs(prob) > 1e6) {
215 oocoutW(caller, Eval) << "RooAbsPdf::getLogVal(" << caller->GetName()
216 << ") WARNING: top-level pdf has a large value: " << prob << std::endl;
217 }
218
219 if (prob < 0) {
220 caller->logEvalError("getLogVal() top-level p.d.f evaluates to a negative number");
221 return RooNaNPacker::packFloatIntoNaN(-prob);
222 }
223
224 if (prob == 0) {
225 caller->logEvalError("getLogVal() top-level p.d.f evaluates to zero");
226
227 return -std::numeric_limits<double>::infinity();
228 }
229
230 if (TMath::IsNaN(prob)) {
231 caller->logEvalError("getLogVal() top-level p.d.f evaluates to NaN");
232
233 return prob;
234 }
235
236 return std::log(prob);
237}
238
239/// To set the fitrange attribute of the PDF and custom ranges for the
240/// observables so that RooPlot can automatically plot the fitting range.
241void resetFitrangeAttributes(RooAbsArg& pdf, RooAbsData const& data, std::string const& baseName,
242 const char* rangeName, bool splitRange)
243{
244 // Clear possible range attributes from previous fits.
245 pdf.removeStringAttribute("fitrange");
246
247 // No fitrange was speficied, so we do nothing. Or "SplitRange" is used, and
248 // then there are no uniquely defined ranges for the observables (as they
249 // are different in each category).
250 if(!rangeName || splitRange) return;
251
252 RooArgSet observables;
253 pdf.getObservables(data.get(), observables) ;
254
255 std::string fitrangeValue;
256 auto subranges = ROOT::Split(rangeName, ",");
257 for (auto const &subrange : subranges) {
258 if (subrange.empty())
259 continue;
260 std::string fitrangeValueSubrange = std::string("fit_") + baseName;
261 if (subranges.size() > 1) {
262 fitrangeValueSubrange += "_" + subrange;
263 }
264 fitrangeValue += fitrangeValueSubrange + ",";
265 for (RooAbsArg *arg : observables) {
266
267 if(arg->isCategory()) continue;
268 auto& observable = static_cast<RooRealVar&>(*arg);
269
270 observable.setRange(fitrangeValueSubrange.c_str(), observable.getMin(subrange.c_str()),
271 observable.getMax(subrange.c_str()));
272 }
273 }
274 pdf.setStringAttribute("fitrange", fitrangeValue.substr(0, fitrangeValue.size() - 1).c_str());
275}
276
277
278} // namespace
279
280using namespace std;
281
283
285
287
288
291
292////////////////////////////////////////////////////////////////////////////////
293/// Default constructor
294
295RooAbsPdf::RooAbsPdf() :_normMgr(this,10)
296{
297 _errorCount = 0 ;
298 _negCount = 0 ;
299 _rawValue = 0 ;
300 _selectComp = false ;
301 _traceCount = 0 ;
302}
303
304
305
306////////////////////////////////////////////////////////////////////////////////
307/// Constructor with name and title only
308
309RooAbsPdf::RooAbsPdf(const char *name, const char *title) :
310 RooAbsReal(name,title), _normMgr(this,10), _selectComp(true)
311{
313 setTraceCounter(0) ;
314}
315
316
317
318////////////////////////////////////////////////////////////////////////////////
319/// Constructor with name, title, and plot range
320
321RooAbsPdf::RooAbsPdf(const char *name, const char *title,
322 double plotMin, double plotMax) :
323 RooAbsReal(name,title,plotMin,plotMax), _normMgr(this,10), _selectComp(true)
324{
326 setTraceCounter(0) ;
327}
328
329
330
331////////////////////////////////////////////////////////////////////////////////
332/// Copy constructor
333
334RooAbsPdf::RooAbsPdf(const RooAbsPdf& other, const char* name) :
335 RooAbsReal(other,name),
336 _normMgr(other._normMgr,this), _selectComp(other._selectComp), _normRange(other._normRange)
337{
340
341 if (other._specGeneratorConfig) {
342 _specGeneratorConfig = std::make_unique<RooNumGenConfig>(*other._specGeneratorConfig);
343 }
344}
345
346
347
348////////////////////////////////////////////////////////////////////////////////
349/// Destructor
350
352{
354
355
356double RooAbsPdf::normalizeWithNaNPacking(double rawVal, double normVal) const {
357
358 if (normVal < 0. || (normVal == 0. && rawVal != 0)) {
359 //Unreasonable normalisations. A zero integral can be tolerated if the function vanishes, though.
360 const std::string msg = "p.d.f normalization integral is zero or negative: " + std::to_string(normVal);
361 logEvalError(msg.c_str());
363 return RooNaNPacker::packFloatIntoNaN(-normVal + (rawVal < 0. ? -rawVal : 0.));
364 }
365
366 if (rawVal < 0.) {
367 logEvalError(Form("p.d.f value is less than zero (%f), trying to recover", rawVal));
369 return RooNaNPacker::packFloatIntoNaN(-rawVal);
370 }
371
372 if (TMath::IsNaN(rawVal)) {
373 logEvalError("p.d.f value is Not-a-Number");
375 return rawVal;
376 }
377
378 return (rawVal == 0. && normVal == 0.) ? 0. : rawVal / normVal;
379}
380
381
382////////////////////////////////////////////////////////////////////////////////
383/// Return current value, normalized by integrating over
384/// the observables in `nset`. If `nset` is 0, the unnormalized value
385/// is returned. All elements of `nset` must be lvalues.
386///
387/// Unnormalized values are not cached.
388/// Doing so would be complicated as `_norm->getVal()` could
389/// spoil the cache and interfere with returning the cached
390/// return value. Since unnormalized calls are typically
391/// done in integration calls, there is no performance hit.
392
393double RooAbsPdf::getValV(const RooArgSet* nset) const
394{
395
396 // Special handling of case without normalization set (used in numeric integration of pdfs)
397 if (!nset) {
398 RooArgSet const* tmp = _normSet ;
399 _normSet = nullptr ;
400 double val = evaluate() ;
401 _normSet = tmp ;
402
403 return TMath::IsNaN(val) ? 0. : val;
404 }
405
406
407 // Process change in last data set used
408 bool nintChanged(false) ;
409 if (!isActiveNormSet(nset) || _norm==0) {
410 nintChanged = syncNormalization(nset) ;
411 }
412
413 // Return value of object. Calculated if dirty, otherwise cached value is returned.
414 if (isValueDirty() || nintChanged || _norm->isValueDirty()) {
415
416 // Evaluate numerator
417 const double rawVal = evaluate();
418
419 // Evaluate denominator
420 const double normVal = _norm->getVal();
421
422 _value = normalizeWithNaNPacking(rawVal, normVal);
423
425 }
426
427 return _value ;
428}
429
430
431////////////////////////////////////////////////////////////////////////////////
432/// Analytical integral with normalization (see RooAbsReal::analyticalIntegralWN() for further information).
433///
434/// This function applies the normalization specified by `normSet` to the integral returned
435/// by RooAbsReal::analyticalIntegral(). The passthrough scenario (code=0) is also changed
436/// to return a normalized answer.
437
438double RooAbsPdf::analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName) const
439{
440 cxcoutD(Eval) << "RooAbsPdf::analyticalIntegralWN(" << GetName() << ") code = " << code << " normset = " << (normSet?*normSet:RooArgSet()) << endl ;
441
442
443 if (code==0) return getVal(normSet) ;
444 if (normSet) {
445 return analyticalIntegral(code,rangeName) / getNorm(normSet) ;
446 } else {
447 return analyticalIntegral(code,rangeName) ;
448 }
449}
450
451
452
453////////////////////////////////////////////////////////////////////////////////
454/// Check that passed value is positive and not 'not-a-number'. If
455/// not, print an error, until the error counter reaches its set
456/// maximum.
457
459{
460 // check for a math error or negative value
461 bool error(false) ;
462 if (TMath::IsNaN(value)) {
463 logEvalError(Form("p.d.f value is Not-a-Number (%f), forcing value to zero",value)) ;
464 error=true ;
465 }
466 if (value<0) {
467 logEvalError(Form("p.d.f value is less than zero (%f), forcing value to zero",value)) ;
468 error=true ;
469 }
470
471 // do nothing if we are no longer tracing evaluations and there was no error
472 if(!error) return error ;
473
474 // otherwise, print out this evaluations input values and result
475 if(++_errorCount <= 10) {
476 cxcoutD(Tracing) << "*** Evaluation Error " << _errorCount << " ";
477 if(_errorCount == 10) cxcoutD(Tracing) << "(no more will be printed) ";
478 }
479 else {
480 return error ;
481 }
482
483 Print() ;
484 return error ;
485}
486
487
488////////////////////////////////////////////////////////////////////////////////
489/// Get normalisation term needed to normalise the raw values returned by
490/// getVal(). Note that `getVal(normalisationVariables)` will automatically
491/// apply the normalisation term returned here.
492/// \param nset Set of variables to normalise over.
493double RooAbsPdf::getNorm(const RooArgSet* nset) const
494{
495 if (!nset) return 1 ;
496
497 syncNormalization(nset,true) ;
498 if (_verboseEval>1) cxcoutD(Tracing) << ClassName() << "::getNorm(" << GetName() << "): norm(" << _norm << ") = " << _norm->getVal() << endl ;
499
500 double ret = _norm->getVal() ;
501 if (ret==0.) {
502 if(++_errorCount <= 10) {
503 coutW(Eval) << "RooAbsPdf::getNorm(" << GetName() << ":: WARNING normalization is zero, nset = " ; nset->Print("1") ;
504 if(_errorCount == 10) coutW(Eval) << "RooAbsPdf::getNorm(" << GetName() << ") INFO: no more messages will be printed " << endl ;
505 }
506 }
507
508 return ret ;
509}
510
511
512
513////////////////////////////////////////////////////////////////////////////////
514/// Return pointer to RooAbsReal object that implements calculation of integral over observables iset in range
515/// rangeName, optionally taking the integrand normalized over observables nset
516
517const RooAbsReal* RooAbsPdf::getNormObj(const RooArgSet* nset, const RooArgSet* iset, const TNamed* rangeName) const
518{
519
520 // Check normalization is already stored
521 CacheElem* cache = (CacheElem*) _normMgr.getObj(nset,iset,0,rangeName) ;
522 if (cache) {
523 return cache->_norm ;
524 }
525
526 // If not create it now
527 RooArgSet depList;
528 getObservables(iset, depList);
529 RooAbsReal* norm = std::unique_ptr<RooAbsReal>{createIntegral(depList,*nset, *getIntegratorConfig(), RooNameReg::str(rangeName))}.release();
530
531 // Store it in the cache
532 cache = new CacheElem(*norm) ;
533 _normMgr.setObj(nset,iset,cache,rangeName) ;
534
535 // And return the newly created integral
536 return norm ;
537}
538
539
540
541////////////////////////////////////////////////////////////////////////////////
542/// Verify that the normalization integral cached with this PDF
543/// is valid for given set of normalization observables.
544///
545/// If not, the cached normalization integral (if any) is deleted
546/// and a new integral is constructed for use with 'nset'.
547/// Elements in 'nset' can be discrete and real, but must be lvalues.
548///
549/// For functions that declare to be self-normalized by overloading the
550/// selfNormalized() function, a unit normalization is always constructed.
551
552bool RooAbsPdf::syncNormalization(const RooArgSet* nset, bool adjustProxies) const
553{
554 setActiveNormSet(nset);
555
556 // Check if data sets are identical
557 CacheElem* cache = (CacheElem*) _normMgr.getObj(nset) ;
558 if (cache) {
559
560 bool nintChanged = (_norm!=cache->_norm) ;
561 _norm = cache->_norm ;
562
563 // In the past, this condition read `if (nintChanged && adjustProxies)`.
564 // However, the cache checks if the nset was already cached **by content**,
565 // and not by RooArgSet instance! So it can happen that the normalization
566 // set object is different, but the integral object is the same, in which
567 // case it would be wrong to not adjust the proxies. They always have to be
568 // adjusted when the nset changed, which is always the case when
569 // `syncNormalization()` is called.
570 if (adjustProxies) {
571 // Update dataset pointers of proxies
572 ((RooAbsPdf*) this)->setProxyNormSet(nset) ;
573 }
574
575 return nintChanged ;
576 }
577
578 // Update dataset pointers of proxies
579 if (adjustProxies) {
580 ((RooAbsPdf*) this)->setProxyNormSet(nset) ;
581 }
582
583 RooArgSet depList;
584 getObservables(nset, depList);
585
586 if (_verboseEval>0) {
587 if (!selfNormalized()) {
588 cxcoutD(Tracing) << ClassName() << "::syncNormalization(" << GetName()
589 << ") recreating normalization integral " << endl ;
590 depList.printStream(ccoutD(Tracing),kName|kValue|kArgs,kSingleLine) ;
591 } else {
592 cxcoutD(Tracing) << ClassName() << "::syncNormalization(" << GetName() << ") selfNormalized, creating unit norm" << endl;
593 }
594 }
595
596 // Destroy old normalization & create new
597 if (selfNormalized() || !dependsOn(depList)) {
598 auto ntitle = std::string(GetTitle()) + " Unit Normalization";
599 auto nname = std::string(GetName()) + "_UnitNorm";
600 _norm = new RooRealVar(nname.c_str(),ntitle.c_str(),1) ;
601 } else {
602 const char* nr = (_normRangeOverride.Length()>0 ? _normRangeOverride.Data() : (_normRange.Length()>0 ? _normRange.Data() : 0)) ;
603
604// cout << "RooAbsPdf::syncNormalization(" << GetName() << ") rangeName for normalization is " << (nr?nr:"<null>") << endl ;
605 RooAbsReal* normInt = std::unique_ptr<RooAbsReal>{createIntegral(depList,*getIntegratorConfig(),nr)}.release();
606 normInt->getVal() ;
607// cout << "resulting normInt = " << normInt->GetName() << endl ;
608
609 const char* cacheParamsStr = getStringAttribute("CACHEPARAMINT") ;
610 if (cacheParamsStr && strlen(cacheParamsStr)) {
611
612 std::unique_ptr<RooArgSet> intParams{normInt->getVariables()} ;
613
614 RooArgSet cacheParams = RooHelpers::selectFromArgSet(*intParams, cacheParamsStr);
615
616 if (!cacheParams.empty()) {
617 cxcoutD(Caching) << "RooAbsReal::createIntObj(" << GetName() << ") INFO: constructing " << cacheParams.getSize()
618 << "-dim value cache for integral over " << depList << " as a function of " << cacheParams << " in range " << (nr?nr:"<default>") << endl ;
619 string name = Form("%s_CACHE_[%s]",normInt->GetName(),cacheParams.contentsString().c_str()) ;
620 RooCachedReal* cachedIntegral = new RooCachedReal(name.c_str(),name.c_str(),*normInt,cacheParams) ;
621 cachedIntegral->setInterpolationOrder(2) ;
622 cachedIntegral->addOwnedComponents(*normInt) ;
623 cachedIntegral->setCacheSource(true) ;
624 if (normInt->operMode()==ADirty) {
625 cachedIntegral->setOperMode(ADirty) ;
626 }
627 normInt= cachedIntegral ;
628 }
629
630 }
631 _norm = normInt ;
632 }
633
634 // Register new normalization with manager (takes ownership)
635 cache = new CacheElem(*_norm) ;
636 _normMgr.setObj(nset,cache) ;
637
638// cout << "making new object " << _norm->GetName() << endl ;
639
640 return true ;
641}
642
643
644
645////////////////////////////////////////////////////////////////////////////////
646/// Reset error counter to given value, limiting the number
647/// of future error messages for this pdf to 'resetValue'
648
650{
651 _errorCount = resetValue ;
652 _negCount = resetValue ;
653}
654
655
656
657////////////////////////////////////////////////////////////////////////////////
658/// Reset trace counter to given value, limiting the
659/// number of future trace messages for this pdf to 'value'
660
662{
663 if (!allNodes) {
665 return ;
666 } else {
667 RooArgList branchList ;
668 branchNodeServerList(&branchList) ;
669 for(auto * pdf : dynamic_range_cast<RooAbsPdf*>(branchList)) {
670 if (pdf) pdf->setTraceCounter(value,false) ;
671 }
672 }
673
674}
675
676
677
678
679////////////////////////////////////////////////////////////////////////////////
680/// Return the log of the current value with given normalization
681/// An error message is printed if the argument of the log is negative.
682
683double RooAbsPdf::getLogVal(const RooArgSet* nset) const
684{
685 return getLog(getVal(nset), this);
686}
687
688
689////////////////////////////////////////////////////////////////////////////////
690/// Check for infinity or NaN.
691/// \param[in] inputs Array to check
692/// \return True if either infinity or NaN were found.
693namespace {
694template<class T>
695bool checkInfNaNNeg(const T& inputs) {
696 // check for a math error or negative value
697 bool inf = false;
698 bool nan = false;
699 bool neg = false;
700
701 for (double val : inputs) { //CHECK_VECTORISE
702 inf |= !std::isfinite(val);
703 nan |= TMath::IsNaN(val); // Works also during fast math
704 neg |= val < 0;
705 }
706
707 return inf || nan || neg;
708}
709}
710
711
712////////////////////////////////////////////////////////////////////////////////
713/// Scan through outputs and fix+log all nans and negative values.
714/// \param[in,out] outputs Array to be scanned & fixed.
715/// \param[in] begin Begin of event range. Only needed to print the correct event number
716/// where the error occurred.
717void RooAbsPdf::logBatchComputationErrors(RooSpan<const double>& outputs, std::size_t begin) const {
718 for (unsigned int i=0; i<outputs.size(); ++i) {
719 const double value = outputs[i];
720 if (TMath::IsNaN(outputs[i])) {
721 logEvalError(Form("p.d.f value of (%s) is Not-a-Number (%f) for entry %zu",
722 GetName(), value, begin+i));
723 } else if (!std::isfinite(outputs[i])){
724 logEvalError(Form("p.d.f value of (%s) is (%f) for entry %zu",
725 GetName(), value, begin+i));
726 } else if (outputs[i] < 0.) {
727 logEvalError(Form("p.d.f value of (%s) is less than zero (%f) for entry %zu",
728 GetName(), value, begin+i));
729 }
730 }
731}
732
733
735 for (std::size_t i = 0; i < pdfValues.size(); ++i) {
736 output[i] = getLog(pdfValues[i], this);
737 }
738}
739
740////////////////////////////////////////////////////////////////////////////////
741/// Return the extended likelihood term (\f$ N_\mathrm{expect} - N_\mathrm{observed} \cdot \log(N_\mathrm{expect} \f$)
742/// of this PDF for the given number of observed events.
743///
744/// For successful operation, the PDF implementation must indicate that
745/// it is extendable by overloading `canBeExtended()`, and must
746/// implement the `expectedEvents()` function.
747///
748/// \param[in] observed The number of observed events.
749/// \param[in] nset The normalization set when asking the pdf for the expected
750/// number of events.
751/// \param[in] observedSumW2 The number of observed events when weighting with
752/// squared weights. If non-zero, the weight-squared error
753/// correction is applied to the extended term.
754/// \param[in] doOffset Offset the extended term by a counterterm where the
755/// expected number of events equals the observed number of events.
756/// This constant shift results in a term closer to zero that is
757/// approximately chi-square distributed. It is useful to do this
758/// also when summing multiple NLL terms to avoid numeric precision
759/// loss that happens if you sum multiple terms of different orders
760/// of magnitude.
761///
762/// The weight-squared error correction works as follows:
763/// adjust poisson such that
764/// estimate of \f$N_\mathrm{expect}\f$ stays at the same value, but has a different variance, rescale
765/// both the observed and expected count of the Poisson with a factor \f$ \sum w_{i} / \sum w_{i}^2 \f$
766/// (the effective weight of the Poisson term),
767/// i.e., change \f$\mathrm{Poisson}(N_\mathrm{observed} = \sum w_{i} | N_\mathrm{expect} )\f$
768/// to \f$ \mathrm{Poisson}(\sum w_{i} \cdot \sum w_{i} / \sum w_{i}^2 | N_\mathrm{expect} \cdot \sum w_{i} / \sum w_{i}^2 ) \f$,
769/// weighted by the effective weight \f$ \sum w_{i}^2 / \sum w_{i} \f$ in the likelihood.
770/// Since here we compute the likelihood with the weight square, we need to multiply by the
771/// square of the effective weight:
772/// - \f$ W_\mathrm{expect} = N_\mathrm{expect} \cdot \sum w_{i} / \sum w_{i}^2 \f$ : effective expected entrie
773/// - \f$ W_\mathrm{observed} = \sum w_{i} \cdot \sum w_{i} / \sum w_{i}^2 \f$ : effective observed entries
774///
775/// The extended term for the likelihood weighted by the square of the weight will be then:
776///
777/// \f$ \left(\sum w_{i}^2 / \sum w_{i}\right)^2 \cdot W_\mathrm{expect} - (\sum w_{i}^2 / \sum w_{i})^2 \cdot W_\mathrm{observed} \cdot \log{W_\mathrm{expect}} \f$
778///
779/// aund this is using the previous expressions for \f$ W_\mathrm{expect} \f$ and \f$ W_\mathrm{observed} \f$:
780///
781/// \f$ \sum w_{i}^2 / \sum w_{i} \cdot N_\mathrm{expect} - \sum w_{i}^2 \cdot \log{W_\mathrm{expect}} \f$
782///
783/// Since the weights are constants in the likelihood we can use \f$\log{N_\mathrm{expect}}\f$ instead of \f$\log{W_\mathrm{expect}}\f$.
784///
785/// See also RooAbsPdf::extendedTerm(RooAbsData const& data, bool weightSquared, bool doOffset),
786/// which takes a dataset to extract \f$N_\mathrm{observed}\f$ and the
787/// normalization set.
788double RooAbsPdf::extendedTerm(double sumEntries, RooArgSet const* nset, double sumEntriesW2, bool doOffset) const
789{
790 return extendedTerm(sumEntries, expectedEvents(nset), sumEntriesW2, doOffset);
791}
792
793double RooAbsPdf::extendedTerm(double sumEntries, double expected, double sumEntriesW2, bool doOffset) const
794{
795 // check if this PDF supports extended maximum likelihood fits
796 if(!canBeExtended()) {
797 coutE(InputArguments) << GetName() << ": this PDF does not support extended maximum likelihood"
798 << std::endl;
799 return 0.0;
800 }
801
802 if(expected < 0.0) {
803 coutE(InputArguments) << GetName() << ": calculated negative expected events: " << expected
804 << std::endl;
805 logEvalError("extendedTerm #expected events is <0 return a NaN");
806 return TMath::QuietNaN();
807 }
808
809
810 // Explicitly handle case Nobs=Nexp=0
811 if (std::abs(expected)<1e-10 && std::abs(sumEntries)<1e-10) {
812 return 0.0;
813 }
814
815 // Check for errors in Nexpected
816 if (TMath::IsNaN(expected)) {
817 logEvalError("extendedTerm #expected events is a NaN") ;
818 return TMath::QuietNaN() ;
819 }
820
821 double extra = doOffset
822 ? (expected - sumEntries) - sumEntries * (std::log(expected) - std::log(sumEntries))
823 : expected - sumEntries * std::log(expected);
824
825 if(sumEntriesW2 != 0.0) {
826 extra *= sumEntriesW2 / sumEntries;
827 }
828
829 return extra;
830}
831
832////////////////////////////////////////////////////////////////////////////////
833/// Return the extended likelihood term (\f$ N_\mathrm{expect} - N_\mathrm{observed} \cdot \log(N_\mathrm{expect} \f$)
834/// of this PDF for the given number of observed events.
835///
836/// This function is a wrapper around
837/// RooAbsPdf::extendedTerm(double, const RooArgSet*, double, bool), where the
838/// number of observed events and observables to be used as the normalization
839/// set for the pdf is extracted from a RooAbsData.
840///
841/// For successful operation, the PDF implementation must indicate that
842/// it is extendable by overloading `canBeExtended()`, and must
843/// implement the `expectedEvents()` function.
844///
845/// \param[in] data The RooAbsData to retrieve the set of observables and
846/// number of expected events.
847/// \param[in] weightSquared If set to `true`, the extended term will be scaled by
848/// the ratio of squared event weights over event weights:
849/// \f$ \sum w_{i}^2 / \sum w_{i} \f$.
850/// Indended to be used by fits with the `SumW2Error()` option that
851/// can be passed to RooAbsPdf::fitTo(RooAbsData&, const RooLinkedList&)
852/// (see the documentation of said function to learn more about the
853/// interpretation of fits with squared weights).
854/// \param[in] doOffset See RooAbsPdf::extendedTerm(double, RooArgSet const*, double bool).
855
856double RooAbsPdf::extendedTerm(RooAbsData const& data, bool weightSquared, bool doOffset) const {
857 double sumW = data.sumEntries();
858 double sumW2 = 0.0;
859 if (weightSquared) {
860 sumW2 = data.sumEntriesW2();
861 }
862 return extendedTerm(sumW, data.get(), sumW2, doOffset);
863}
864
865
866////////////////////////////////////////////////////////////////////////////////
867/// Construct representation of -log(L) of PDF with given dataset. If dataset is unbinned, an unbinned likelihood is constructed. If the dataset
868/// is binned, a binned likelihood is constructed.
869///
870/// The following named arguments are supported
871///
872/// <table>
873/// <tr><th> Type of CmdArg <th> Effect on nll
874/// <tr><td> `ConditionalObservables(Args_t &&... argsOrArgSet)` <td> Do not normalize PDF over listed observables.
875// Arguments can either be multiple RooRealVar or a single RooArgSet containing them.
876/// <tr><td> `Extended(bool flag)` <td> Add extended likelihood term, off by default
877/// <tr><td> `Range(const char* name)` <td> Fit only data inside range with given name
878/// <tr><td> `Range(double lo, double hi)` <td> Fit only data inside given range. A range named "fit" is created on the fly on all observables.
879/// Multiple comma separated range names can be specified.
880/// <tr><td> `SumCoefRange(const char* name)` <td> Set the range in which to interpret the coefficients of RooAddPdf components
881/// <tr><td> `NumCPU(int num, int strat)` <td> Parallelize NLL calculation on num CPUs
882/// <table>
883/// <tr><th> Strategy <th> Effect
884/// <tr><td> 0 = RooFit::BulkPartition (Default) <td> Divide events in N equal chunks
885/// <tr><td> 1 = RooFit::Interleave <td> Process event i%N in process N. Recommended for binned data with
886/// a substantial number of zero-bins, which will be distributed across processes more equitably in this strategy
887/// <tr><td> 2 = RooFit::SimComponents <td> Process each component likelihood of a RooSimultaneous fully in a single process
888/// and distribute components over processes. This approach can be benificial if normalization calculation time
889/// dominates the total computation time of a component (since the normalization calculation must be performed
890/// in each process in strategies 0 and 1. However beware that if the RooSimultaneous components do not share many
891/// parameters this strategy is inefficient: as most minuit-induced likelihood calculations involve changing
892/// a single parameter, only 1 of the N processes will be active most of the time if RooSimultaneous components
893/// do not share many parameters
894/// <tr><td> 3 = RooFit::Hybrid <td> Follow strategy 0 for all RooSimultaneous components, except those with less than
895/// 30 dataset entries, for which strategy 2 is followed.
896/// </table>
897/// <tr><td> `BatchMode(bool on)` <td> Batch evaluation mode. See fitTo().
898/// <tr><td> `Optimize(bool flag)` <td> Activate constant term optimization (on by default)
899/// <tr><td> `SplitRange(bool flag)` <td> Use separate fit ranges in a simultaneous fit. Actual range name for each subsample is assumed to
900/// be `rangeName_indexState`, where `indexState` is the state of the master index category of the simultaneous fit.
901/// Using `Range("range"), SplitRange()` as switches, different ranges could be set like this:
902/// ```
903/// myVariable.setRange("range_pi0", 135, 210);
904/// myVariable.setRange("range_gamma", 50, 210);
905/// ```
906/// <tr><td> `Constrain(const RooArgSet&pars)` <td> For p.d.f.s that contain internal parameter constraint terms (that is usually product PDFs, where one
907/// term of the product depends on parameters but not on the observable(s),), only apply constraints to the given subset of parameters.
908/// <tr><td> `ExternalConstraints(const RooArgSet& )` <td> Include given external constraints to likelihood by multiplying them with the original likelihood.
909/// <tr><td> `GlobalObservables(const RooArgSet&)` <td> Define the set of normalization observables to be used for the constraint terms.
910/// If none are specified the constrained parameters are used.
911/// <tr><td> `GlobalObservablesSource(const char* sourceName)` <td> Which source to prioritize for global observable values.
912/// Can be either:
913/// - `data`: to take the values from the dataset,
914/// falling back to the pdf value if a given global observable is not available.
915/// If no `GlobalObservables` or `GlobalObservablesTag` command argument is given, the set
916/// of global observables will be automatically defined to be the set stored in the data.
917/// - `model`: to take all values from the pdf and completely ignore the set of global observables stored in the data
918/// (not even using it to automatically define the set of global observables
919/// if the `GlobalObservables` or `GlobalObservablesTag` command arguments are not given).
920/// The default option is `data`.
921/// <tr><td> `GlobalObservablesTag(const char* tagName)` <td> Define the set of normalization observables to be used for the constraint terms by
922/// a string attribute associated with pdf observables that match the given tagName.
923/// <tr><td> `Verbose(bool flag)` <td> Controls RooFit informational messages in likelihood construction
924/// <tr><td> `CloneData(bool flag)` <td> Use clone of dataset in NLL (default is true).
925/// \warning Deprecated option that is ignored. It is up to the implementation of the NLL creation method if the data is cloned or not.
926/// <tr><td> `Offset(std::string const& mode)` <td> Likelihood offsetting mode. Can be either:
927/// - `"none"` (default): no offsetting
928/// - `"initial"`: Offset likelihood by initial value (so that starting value of FCN in minuit is zero).
929/// This can improve numeric stability in simultaneous fits with components with large likelihood values.
930/// - `"bin"`: Offset by the likelihood bin-by-bin with a template histogram model based on the obersved data.
931/// This results in per-bin values that are all in the same order of magnitude, which reduces precision loss in the sum,
932/// which can drastically improve numeric stability.
933/// Furthermore, 2 * NLL defined like this is approximately chi-square distributed, allowing for goodness-of-fit tests.
934/// <tr><td> `IntegrateBins(double precision)` <td> In binned fits, integrate the PDF over the bins instead of using the probability density at the bin centre.
935/// This can reduce the bias observed when fitting functions with high curvature to binned data.
936/// - precision > 0: Activate bin integration everywhere. Use precision between 0.01 and 1.E-6, depending on binning.
937/// Note that a low precision such as 0.01 might yield identical results to 1.E-4, since the integrator might reach 1.E-4 already in its first
938/// integration step. If lower precision is desired (more speed), a RooBinSamplingPdf has to be created manually, and its integrator
939/// has to be manipulated directly.
940/// - precision = 0: Activate bin integration only for continuous PDFs fit to a RooDataHist.
941/// - precision < 0: Deactivate.
942/// \see RooBinSamplingPdf
943/// <tr><td> `ModularL(bool flag)` <td> Enable or disable modular likelihoods, which will become the default in a future release.
944/// This does not change any user-facing code, but only enables a different likelihood class in the back-end. Note that this
945/// should be set to true for parallel minimization of likelihoods!
946/// Note that it is currently not recommended to use Modular likelihoods without any parallelization enabled in the minimization, since
947/// some features such as offsetting might not yet work in this case.
948/// </table>
949///
950///
951
953{
954 auto baseName = std::string("nll_") + GetName() + "_" + data.GetName();
955
956 // Figure out the integer that corresponds to the default BatchMode option.
957 const int defaultBatchModeInt = RooFit::BatchMode(RooFit::Experimental::defaultBatchMode()).getInt(0);
958
959 // Select the pdf-specific commands
960 RooCmdConfig pc(Form("RooAbsPdf::createNLL(%s)",GetName())) ;
961
962 pc.defineString("rangeName","RangeWithName",0,"",true) ;
963 pc.defineString("addCoefRange","SumCoefRange",0,"") ;
964 pc.defineString("globstag","GlobalObservablesTag",0,"") ;
965 pc.defineString("globssource","GlobalObservablesSource",0,"data") ;
966 pc.defineDouble("rangeLo","Range",0,-999.) ;
967 pc.defineDouble("rangeHi","Range",1,-999.) ;
968 pc.defineInt("splitRange","SplitRange",0,0) ;
969 pc.defineInt("ext","Extended",0,2) ;
970 pc.defineInt("numcpu","NumCPU",0,1) ;
971 pc.defineInt("interleave","NumCPU",1,0) ;
972 pc.defineInt("verbose","Verbose",0,0) ;
973 pc.defineInt("optConst","Optimize",0,0) ;
974 pc.defineInt("cloneData","CloneData", 0, 2);
975 pc.defineSet("projDepSet","ProjectedObservables",0,0) ;
976 pc.defineSet("cPars","Constrain",0,0) ;
977 pc.defineSet("glObs","GlobalObservables",0,0) ;
978 pc.defineInt("doOffset","OffsetLikelihood",0,0) ;
979 pc.defineSet("extCons","ExternalConstraints",0,0) ;
980 pc.defineInt("BatchMode", "BatchMode", 0, defaultBatchModeInt);
981 pc.defineDouble("IntegrateBins", "IntegrateBins", 0, -1.);
982 pc.defineMutex("Range","RangeWithName") ;
983 pc.defineMutex("GlobalObservables","GlobalObservablesTag") ;
984 pc.defineInt("ModularL", "ModularL", 0, 0);
985
986 // New style likelihoods define parallelization through Parallelize(...) on fitTo or attributes on RooMinimizer::Config.
987 pc.defineMutex("ModularL", "NumCPU");
988
989 // New style likelihoods define offsetting on minimizer, not on likelihood
990 pc.defineMutex("ModularL", "OffsetLikelihood");
991
992 // Process and check varargs
993 pc.process(cmdList) ;
994 if (!pc.ok(true)) {
995 return 0 ;
996 }
997
998 if (pc.getInt("ModularL")) {
999 int lut[3] = {2, 1, 0};
1001 static_cast<RooFit::TestStatistics::RooAbsL::Extended>(lut[pc.getInt("ext")])};
1002
1003 RooArgSet cParsSet;
1004 RooArgSet extConsSet;
1005 RooArgSet glObsSet;
1006
1007 if (auto tmp = pc.getSet("cPars"))
1008 cParsSet.add(*tmp);
1009
1010 if (auto tmp = pc.getSet("extCons"))
1011 extConsSet.add(*tmp);
1012
1013 if (auto tmp = pc.getSet("glObs"))
1014 glObsSet.add(*tmp);
1015
1016 const std::string rangeName = pc.getString("globstag", "", false);
1017
1019 builder.Extended(ext)
1020 .ConstrainedParameters(cParsSet)
1021 .ExternalConstraints(extConsSet)
1022 .GlobalObservables(glObsSet)
1023 .GlobalObservablesTag(rangeName.c_str());
1024
1025 return RooFit::Detail::owningPtr(std::make_unique<RooFit::TestStatistics::RooRealL>("likelihood", "", builder.build()));
1026 }
1027
1028 // Decode command line arguments
1029 const char* rangeName = pc.getString("rangeName",0,true) ;
1030 const char* addCoefRangeName = pc.getString("addCoefRange",0,true) ;
1031 const bool ext = interpretExtendedCmdArg(*this, pc.getInt("ext")) ;
1032 Int_t numcpu = pc.getInt("numcpu") ;
1033 Int_t numcpu_strategy = pc.getInt("interleave");
1034 // strategy 3 works only for RooSimultaneous.
1035 if (numcpu_strategy==3 && !this->InheritsFrom("RooSimultaneous") ) {
1036 coutW(Minimization) << "Cannot use a NumCpu Strategy = 3 when the pdf is not a RooSimultaneous, "
1037 "falling back to default strategy = 0" << endl;
1038 numcpu_strategy = 0;
1039 }
1040 RooFit::MPSplit interl = (RooFit::MPSplit) numcpu_strategy;
1041
1042 Int_t splitRange = pc.getInt("splitRange") ;
1043 bool verbose = pc.getInt("verbose") ;
1044 Int_t optConst = pc.getInt("optConst") ;
1045 Int_t cloneData = pc.getInt("cloneData") ;
1046 auto offset = static_cast<RooFit::OffsetMode>(pc.getInt("doOffset"));
1047
1048 if(offset == RooFit::OffsetMode::Bin && dynamic_cast<RooDataSet*>(&data)) {
1049 coutE(Minimization) << "The Offset(\"bin\") option doesn't suppot fits to RooDataSet yet, only to RooDataHist."
1050 " Falling back to no offsetting." << endl;
1052 }
1053
1054 // If no explicit cloneData command is specified, cloneData is set to true if optimization is activated
1055 if (cloneData==2) {
1056 cloneData = optConst ;
1057 }
1058
1059 if (pc.hasProcessed("Range")) {
1060 double rangeLo = pc.getDouble("rangeLo") ;
1061 double rangeHi = pc.getDouble("rangeHi") ;
1062
1063 // Create range with name 'fit' with above limits on all observables
1064 RooArgSet obs;
1065 getObservables(data.get(), obs) ;
1066 for (auto arg : obs) {
1067 RooRealVar* rrv = dynamic_cast<RooRealVar*>(arg) ;
1068 if (rrv) rrv->setRange("fit",rangeLo,rangeHi) ;
1069 }
1070
1071 // Set range name to be fitted to "fit"
1072 rangeName = "fit" ;
1073 }
1074
1075 // Set the fitrange attribute of th PDF, add observables ranges for plotting
1076 resetFitrangeAttributes(*this, data, baseName, rangeName, splitRange);
1077
1078 RooArgSet projDeps ;
1079 auto tmp = pc.getSet("projDepSet");
1080 if (tmp) {
1081 projDeps.add(*tmp) ;
1082 }
1083
1084 const std::string globalObservablesSource = pc.getString("globssource","data",false);
1085 if(globalObservablesSource != "data" && globalObservablesSource != "model") {
1086 std::string errMsg = "RooAbsPdf::fitTo: GlobalObservablesSource can only be \"data\" or \"model\"!";
1087 coutE(InputArguments) << errMsg << std::endl;
1088 throw std::invalid_argument(errMsg);
1089 }
1090 const bool takeGlobalObservablesFromData = globalObservablesSource == "data";
1091
1092 // Lambda function to create the correct constraint term for a PDF. In old
1093 // RooFit, we use this PDF itself as the argument, for the new BatchMode
1094 // we're passing a clone.
1095 auto createConstr = [&](RooAbsPdf const& pdf, bool removeConstraintsFromPdf=false) -> std::unique_ptr<RooAbsReal> {
1096 return createConstraintTerm(
1097 baseName + "_constr", // name
1098 pdf, // pdf
1099 data, // data
1100 pc.getSet("cPars"), // Constrain RooCmdArg
1101 pc.getSet("extCons"), // ExternalConstraints RooCmdArg
1102 pc.getSet("glObs"), // GlobalObservables RooCmdArg
1103 pc.getString("globstag",0,true), // GlobalObservablesTag RooCmdArg
1104 takeGlobalObservablesFromData, // From GlobalObservablesSource RooCmdArg
1105 removeConstraintsFromPdf
1106 );
1107 };
1108
1109 auto batchMode = static_cast<RooFit::BatchModeOption>(pc.getInt("BatchMode"));
1110
1111 // Construct BatchModeNLL if requested
1112 if (batchMode != RooFit::BatchModeOption::Off) {
1113
1114 // Set the normalization range. We need to do it now, because it will be
1115 // considered in `compileForNormSet`.
1116 TString oldNormRange = _normRange;
1117 setNormRange(rangeName);
1118
1119 auto normSet = new RooArgSet; // INTENTIONAL LEAK FOR NOW!
1120 getObservables(data.get(), *normSet);
1121 normSet->remove(projDeps, true, true);
1122
1123 this->setAttribute("SplitRange", splitRange);
1124 this->setStringAttribute("RangeName", rangeName);
1125
1126 std::unique_ptr<RooAbsPdf> pdfClone = RooFit::Detail::compileForNormSet(*this, *normSet);
1127
1128 // reset attributes
1129 this->setAttribute("SplitRange", false);
1130 this->setStringAttribute("RangeName", nullptr);
1131
1132 // Reset the normalization range
1133 _normRange = oldNormRange;
1134
1135 if (addCoefRangeName) {
1136 cxcoutI(Fitting) << "RooAbsPdf::fitTo(" << GetName()
1137 << ") fixing interpretation of coefficients of any component to range "
1138 << addCoefRangeName << "\n";
1139 pdfClone->fixAddCoefRange(addCoefRangeName, false);
1140 }
1141
1142 std::unique_ptr<RooAbsReal> compiledConstr;
1143 if(std::unique_ptr<RooAbsReal> constr = createConstr(*this)) {
1144 compiledConstr = RooFit::Detail::compileForNormSet(*constr, *data.get());
1145 compiledConstr->addOwnedComponents(std::move(constr));
1146 }
1147
1149 *pdfClone,
1150 data,
1151 std::move(compiledConstr),
1152 rangeName ? rangeName : "",
1153 projDeps,
1154 ext,
1155 pc.getDouble("IntegrateBins"),
1156 offset);
1157
1158 auto driver = std::make_unique<ROOT::Experimental::RooFitDriver>(*nll, batchMode);
1159
1160 auto driverWrapper = std::make_unique<ROOT::Experimental::RooAbsRealWrapper>(
1161 std::move(driver), rangeName ? rangeName : "", dynamic_cast<RooSimultaneous *>(pdfClone.get()), takeGlobalObservablesFromData);
1162 driverWrapper->setData(data, false);
1163 driverWrapper->addOwnedComponents(std::move(nll));
1164 driverWrapper->addOwnedComponents(std::move(pdfClone));
1165
1166 return RooFit::Detail::owningPtr<RooAbsReal>(std::move(driverWrapper));
1167 }
1168
1169 // Construct NLL
1171 std::unique_ptr<RooAbsReal> nll ;
1173 cfg.addCoefRangeName = addCoefRangeName ? addCoefRangeName : "";
1174 cfg.nCPU = numcpu;
1175 cfg.interleave = interl;
1176 cfg.verbose = verbose;
1177 cfg.splitCutRange = static_cast<bool>(splitRange);
1178 cfg.cloneInputData = static_cast<bool>(cloneData);
1179 cfg.integrateOverBinsPrecision = pc.getDouble("IntegrateBins");
1180 cfg.binnedL = false;
1181 cfg.takeGlobalObservablesFromData = takeGlobalObservablesFromData;
1182 cfg.rangeName = rangeName ? rangeName : "";
1183 auto nllVar = std::make_unique<RooNLLVar>(baseName.c_str(),"-log(likelihood)",*this,data,projDeps, ext, cfg);
1184 nllVar->enableBinOffsetting(offset == RooFit::OffsetMode::Bin);
1185 nll = std::move(nllVar);
1187
1188 // Include constraints, if any, in likelihood
1189 if (std::unique_ptr<RooAbsReal> constraintTerm = createConstr(*this)) {
1190
1191 // Even though it is technically only required when the computation graph
1192 // is changed because global observables are taken from data, it is safer
1193 // to clone the constraint model in general to reset the normalization
1194 // integral caches and avoid ASAN build failures (the PDF of the main
1195 // measurement is cloned too anyway, so not much overhead). This can be
1196 // reconsidered after the caching of normalization sets by pointer is changed
1197 // to a more memory-safe solution.
1198 constraintTerm = RooHelpers::cloneTreeWithSameParameters(*constraintTerm, data.get());
1199
1200 // Redirect the global observables to the ones from the dataset if applicable.
1201 constraintTerm->setData(data, false);
1202
1203 // The computation graph for the constraints is very small, no need to do
1204 // the tracking of clean and dirty nodes here.
1205 constraintTerm->setOperMode(RooAbsArg::ADirty);
1206
1207 auto orignll = std::move(nll) ;
1208 nll = std::make_unique<RooAddition>(Form("%s_with_constr",baseName.c_str()),"nllWithCons",RooArgSet(*orignll,*constraintTerm)) ;
1209 nll->addOwnedComponents(std::move(orignll),std::move(constraintTerm)) ;
1210 }
1211
1212 if (optConst) {
1213 nll->constOptimizeTestStatistic(RooAbsArg::Activate,optConst>1) ;
1214 }
1215
1217 nll->enableOffsetting(true) ;
1218 }
1219
1220 return RooFit::Detail::owningPtr(std::move(nll));
1221}
1222
1223
1224////////////////////////////////////////////////////////////////////////////////
1225/// Use the asymptotically correct approach to estimate errors in the presence of weights.
1226/// This is slower but more accurate than `SumW2Error`. See also https://arxiv.org/abs/1911.01303).
1227/// Applies the calculated covaraince matrix to the RooMinimizer and returns
1228/// the quality of the covariance matrix.
1229/// See also the documentation of RooAbsPdf::fitTo(), where this function is used.
1230/// \param[in] minimizer The RooMinimizer to get the fit result from. The state
1231/// of the minimizer will be altered by this function: the covariance
1232/// matrix caltulated here will be applied to it via
1233/// RooMinimizer::applyCovarianceMatrix().
1234/// \param[in] data The dataset that was used for the fit.
1236{
1237 // Calculated corrected errors for weighted likelihood fits
1238 std::unique_ptr<RooFitResult> rw(minimizer.save());
1239 // Weighted inverse Hessian matrix
1240 const TMatrixDSym &matV = rw->covarianceMatrix();
1241 coutI(Fitting)
1242 << "RooAbsPdf::fitTo(" << this->GetName()
1243 << ") Calculating covariance matrix according to the asymptotically correct approach. If you find this "
1244 "method useful please consider citing https://arxiv.org/abs/1911.01303."
1245 << endl;
1246
1247 // Initialise matrix containing first derivatives
1248 auto nFloatPars = rw->floatParsFinal().getSize();
1249 TMatrixDSym num(nFloatPars);
1250 for (int k = 0; k < nFloatPars; k++) {
1251 for (int l = 0; l < nFloatPars; l++) {
1252 num(k, l) = 0.0;
1253 }
1254 }
1255 RooArgSet obs;
1256 this->getObservables(data.get(), obs);
1257 // Create derivative objects
1258 std::vector<std::unique_ptr<RooDerivative>> derivatives;
1259 const RooArgList &floated = rw->floatParsFinal();
1260 std::unique_ptr<RooArgSet> floatingparams{
1261 static_cast<RooArgSet *>(this->getParameters(data)->selectByAttrib("Constant", false))};
1262 for (const auto paramresult : floated) {
1263 auto paraminternal = static_cast<RooRealVar *>(floatingparams->find(*paramresult));
1264 assert(floatingparams->find(*paramresult)->IsA() == RooRealVar::Class());
1265 derivatives.emplace_back(this->derivative(*paraminternal, obs, 1));
1266 }
1267
1268 // Loop over data
1269 for (int j = 0; j < data.numEntries(); j++) {
1270 // Sets obs to current data point, this is where the pdf will be evaluated
1271 obs.assign(*data.get(j));
1272 // Determine first derivatives
1273 std::vector<double> diffs(floated.getSize(), 0.0);
1274 for (int k = 0; k < floated.getSize(); k++) {
1275 const auto paramresult = static_cast<RooRealVar *>(floated.at(k));
1276 auto paraminternal = static_cast<RooRealVar *>(floatingparams->find(*paramresult));
1277 // first derivative to parameter k at best estimate point for this measurement
1278 double diff = derivatives[k]->getVal();
1279 // need to reset to best fit point after differentiation
1280 *paraminternal = paramresult->getVal();
1281 diffs[k] = diff;
1282 }
1283 // Fill numerator matrix
1284 double prob = getVal(&obs);
1285 for (int k = 0; k < floated.getSize(); k++) {
1286 for (int l = 0; l < floated.getSize(); l++) {
1287 num(k, l) += data.weightSquared() * diffs[k] * diffs[l] / (prob * prob);
1288 }
1289 }
1290 }
1291 num.Similarity(matV);
1292
1293 // Propagate corrected errors to parameters objects
1294 minimizer.applyCovarianceMatrix(num);
1295
1296 // The derivatives are found in RooFit and not with the minimizer (e.g.
1297 // minuit), so the quality of the corrected covariance matrix corresponds to
1298 // the quality of the original covariance matrix
1299 return rw->covQual();
1300}
1301
1302
1303////////////////////////////////////////////////////////////////////////////////
1304/// Apply correction to errors and covariance matrix. This uses two covariance
1305/// matrices, one with the weights, the other with squared weights, to obtain
1306/// the correct errors for weighted likelihood fits.
1307/// Applies the calculated covaraince matrix to the RooMinimizer and returns
1308/// the quality of the covariance matrix.
1309/// See also the documentation of RooAbsPdf::fitTo(), where this function is used.
1310/// \param[in] minimizer The RooMinimizer to get the fit result from. The state
1311/// of the minimizer will be altered by this function: the covariance
1312/// matrix caltulated here will be applied to it via
1313/// RooMinimizer::applyCovarianceMatrix().
1314/// \param[in] nll The NLL object that was used for the fit.
1316{
1317 // Calculated corrected errors for weighted likelihood fits
1318 std::unique_ptr<RooFitResult> rw{minimizer.save()};
1319 nll.applyWeightSquared(true);
1320 coutI(Fitting) << "RooAbsPdf::fitTo(" << this->GetName()
1321 << ") Calculating sum-of-weights-squared correction matrix for covariance matrix"
1322 << std::endl;
1323 minimizer.hesse();
1324 std::unique_ptr<RooFitResult> rw2{minimizer.save()};
1325 nll.applyWeightSquared(false);
1326
1327 // Apply correction matrix
1328 const TMatrixDSym &matV = rw->covarianceMatrix();
1329 TMatrixDSym matC = rw2->covarianceMatrix();
1331 if (!decomp) {
1332 coutE(Fitting) << "RooAbsPdf::fitTo(" << this->GetName()
1333 << ") ERROR: Cannot apply sum-of-weights correction to covariance matrix: correction "
1334 "matrix calculated with weight-squared is singular"
1335 << std::endl;
1336 return -1;
1337 }
1338
1339 // replace C by its inverse
1340 decomp.Invert(matC);
1341 // the class lies about the matrix being symmetric, so fill in the
1342 // part above the diagonal
1343 for (int i = 0; i < matC.GetNrows(); ++i) {
1344 for (int j = 0; j < i; ++j) {
1345 matC(j, i) = matC(i, j);
1346 }
1347 }
1348 matC.Similarity(matV);
1349 // C now contiains V C^-1 V
1350 // Propagate corrected errors to parameters objects
1351 minimizer.applyCovarianceMatrix(matC);
1352
1353 return std::min(rw->covQual(), rw2->covQual());
1354}
1355
1356
1357////////////////////////////////////////////////////////////////////////////////
1358/// Minimizes a given NLL variable by finding the optimal parameters with the
1359/// RooMinimzer. The NLL variable can be created with RooAbsPdf::createNLL.
1360/// If you are looking for a function that combines likelihood creation with
1361/// fitting, see RooAbsPdf::fitTo.
1362/// \param[in] nll The negative log-likelihood variable to minimize.
1363/// \param[in] data The dataset that was als used for the NLL. It's a necessary
1364/// parameter because it is used in the asymptotic error correction.
1365/// \param[in] cfg Configuration struct with all the configuration options for
1366/// the RooMinimizer. These are a subset of the options that you can
1367/// also pass to RooAbsPdf::fitTo via the RooFit command arguments.
1368std::unique_ptr<RooFitResult> RooAbsPdf::minimizeNLL(RooAbsReal & nll,
1369 RooAbsData const& data, MinimizerConfig const& cfg) {
1370
1371 // Determine if the dataset has weights
1372 bool weightedData = data.isNonPoissonWeighted();
1373
1374 // Warn user that a method to determine parameter uncertainties should be provided if weighted data is offered
1375 if (weightedData && cfg.doSumW2==-1 && cfg.doAsymptotic==-1) {
1376 coutW(InputArguments) << "RooAbsPdf::fitTo(" << GetName() << ") WARNING: a likelihood fit is requested of what appears to be weighted data.\n"
1377 << " While the estimated values of the parameters will always be calculated taking the weights into account,\n"
1378 << " there are multiple ways to estimate the errors of the parameters. You are advised to make an \n"
1379 << " explicit choice for the error calculation:\n"
1380 << " - Either provide SumW2Error(true), to calculate a sum-of-weights-corrected HESSE error matrix\n"
1381 << " (error will be proportional to the number of events in MC).\n"
1382 << " - Or provide SumW2Error(false), to return errors from original HESSE error matrix\n"
1383 << " (which will be proportional to the sum of the weights, i.e., a dataset with <sum of weights> events).\n"
1384 << " - Or provide AsymptoticError(true), to use the asymptotically correct expression\n"
1385 << " (for details see https://arxiv.org/abs/1911.01303)."
1386 << endl ;
1387 }
1388
1389 if (cfg.minos && (cfg.doSumW2==1 || cfg.doAsymptotic == 1)) {
1390 coutE(InputArguments) << "RooAbsPdf::fitTo(" << GetName() << "): sum-of-weights and asymptotic error correction do not work with MINOS errors. Not fitting." << endl;
1391 return nullptr;
1392 }
1393 if (cfg.doAsymptotic==1 && cfg.minos) {
1394 coutW(InputArguments) << "RooAbsPdf::fitTo(" << GetName() << ") WARNING: asymptotic correction does not apply to MINOS errors" << endl ;
1395 }
1396
1397 //avoid setting both SumW2 and Asymptotic for uncertainty correction
1398 if (cfg.doSumW2==1 && cfg.doAsymptotic==1) {
1399 coutE(InputArguments) << "RooAbsPdf::fitTo(" << GetName() << ") ERROR: Cannot compute both asymptotically correct and SumW2 errors." << endl ;
1400 return nullptr;
1401 }
1402
1403 // Instantiate RooMinimizer
1404 RooMinimizer::Config minimizerConfig;
1405 minimizerConfig.enableParallelGradient = cfg.enableParallelGradient;
1406 minimizerConfig.enableParallelDescent = cfg.enableParallelDescent;
1407 minimizerConfig.parallelize = cfg.parallelize;
1408 minimizerConfig.timingAnalysis = cfg.timingAnalysis;
1409 minimizerConfig.offsetting = cfg.doOffset;
1410 RooMinimizer m(nll, minimizerConfig);
1411
1412 m.setMinimizerType(cfg.minType.c_str());
1413 m.setEvalErrorWall(cfg.doEEWall);
1414 m.setRecoverFromNaNStrength(cfg.recoverFromNaN);
1415 m.setPrintEvalErrors(cfg.numee);
1416 if (cfg.maxCalls > 0) m.setMaxFunctionCalls(cfg.maxCalls);
1417 if (cfg.printLevel!=1) m.setPrintLevel(cfg.printLevel);
1418 if (cfg.optConst) m.optimizeConst(cfg.optConst); // Activate constant term optimization
1419 if (cfg.verbose) m.setVerbose(1); // Activate verbose options
1420 if (cfg.doTimer) m.setProfile(1); // Activate timer options
1421 if (cfg.strat!=1) m.setStrategy(cfg.strat); // Modify fit strategy
1422 if (cfg.initHesse) m.hesse(); // Initialize errors with hesse
1423 m.minimize(cfg.minType.c_str(), cfg.minAlg.c_str()); // Minimize using chosen algorithm
1424 if (cfg.hesse) m.hesse(); // Evaluate errors with Hesse
1425
1426 int corrCovQual = -1;
1427
1428 if (m.getNPar()>0) {
1429 if (cfg.doAsymptotic == 1) corrCovQual = calcAsymptoticCorrectedCovariance(m, data); // Asymptotically correct
1430 if (cfg.doSumW2 == 1) corrCovQual = calcSumW2CorrectedCovariance(m, nll);
1431 }
1432
1433 if (cfg.minos) cfg.minosSet ? m.minos(*cfg.minosSet) : m.minos(); // Evaluate errs with Minos
1434
1435 // Optionally return fit result
1436 std::unique_ptr<RooFitResult> ret;
1437 if (cfg.doSave) {
1438 auto name = std::string("fitresult_") + GetName() + "_" + data.GetName();
1439 auto title = std::string("Result of fit of p.d.f. ") + GetName() + " to dataset " + data.GetName();
1440 ret.reset(m.save(name.c_str(),title.c_str()));
1441 if((cfg.doSumW2==1 || cfg.doAsymptotic==1) && m.getNPar()>0) ret->setCovQual(corrCovQual);
1442 }
1443
1444 if (cfg.optConst) m.optimizeConst(0) ;
1445 return ret ;
1446}
1447
1448
1449
1450////////////////////////////////////////////////////////////////////////////////
1451/// Fit PDF to given dataset. If dataset is unbinned, an unbinned maximum likelihood is performed. If the dataset
1452/// is binned, a binned maximum likelihood is performed. By default the fit is executed through the MINUIT
1453/// commands MIGRAD, HESSE in succession.
1454/// \param[in] data Data to fit the PDF to
1455/// \param[in] arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8 One or more arguments to control the behaviour of the fit
1456/// \return RooFitResult with fit status and parameters if option Save() is used, `nullptr` otherwise. The user takes ownership of the fit result.
1457///
1458/// The following named arguments are supported
1459///
1460/// <table>
1461/// <tr><th> Type of CmdArg <th> Options to control construction of -log(L)
1462/// <tr><td> `ConditionalObservables(Args_t &&... argsOrArgSet)` <td> Do not normalize PDF over listed observables.
1463// Arguments can either be multiple RooRealVar or a single RooArgSet containing them.
1464/// <tr><td> `Extended(bool flag)` <td> Add extended likelihood term, off by default
1465/// <tr><td> `Range(const char* name)` <td> Fit only data inside range with given name. Multiple comma-separated range names can be specified.
1466/// In this case, the unnormalized PDF \f$f(x)\f$ is normalized by the integral over all ranges \f$r_i\f$:
1467/// \f[
1468/// p(x) = \frac{f(x)}{\sum_i \int_{r_i} f(x) dx}.
1469/// \f]
1470/// <tr><td> `Range(double lo, double hi)` <td> Fit only data inside given range. A range named "fit" is created on the fly on all observables.
1471/// <tr><td> `SumCoefRange(const char* name)` <td> Set the range in which to interpret the coefficients of RooAddPdf components
1472/// <tr><td> `NumCPU(int num, int strat)` <td> Parallelize NLL calculation on `num` CPUs
1473/// <table>
1474/// <tr><th> Strategy <th> Effect
1475/// <tr><td> 0 = RooFit::BulkPartition (Default) <td> Divide events in N equal chunks
1476/// <tr><td> 1 = RooFit::Interleave <td> Process event i%N in process N. Recommended for binned data with
1477/// a substantial number of zero-bins, which will be distributed across processes more equitably in this strategy
1478/// <tr><td> 2 = RooFit::SimComponents <td> Process each component likelihood of a RooSimultaneous fully in a single process
1479/// and distribute components over processes. This approach can be benificial if normalization calculation time
1480/// dominates the total computation time of a component (since the normalization calculation must be performed
1481/// in each process in strategies 0 and 1. However beware that if the RooSimultaneous components do not share many
1482/// parameters this strategy is inefficient: as most minuit-induced likelihood calculations involve changing
1483/// a single parameter, only 1 of the N processes will be active most of the time if RooSimultaneous components
1484/// do not share many parameters
1485/// <tr><td> 3 = RooFit::Hybrid <td> Follow strategy 0 for all RooSimultaneous components, except those with less than
1486/// 30 dataset entries, for which strategy 2 is followed.
1487/// </table>
1488/// <tr><td> `SplitRange(bool flag)` <td> Use separate fit ranges in a simultaneous fit. Actual range name for each subsample is assumed
1489/// to by `rangeName_indexState` where indexState is the state of the master index category of the simultaneous fit.
1490/// Using `Range("range"), SplitRange()` as switches, different ranges could be set like this:
1491/// ```
1492/// myVariable.setRange("range_pi0", 135, 210);
1493/// myVariable.setRange("range_gamma", 50, 210);
1494/// ```
1495/// <tr><td> `Constrain(const RooArgSet&pars)` <td> For p.d.f.s that contain internal parameter constraint terms (that is usually product PDFs, where one
1496/// term of the product depends on parameters but not on the observable(s),), only apply constraints to the given subset of parameters.
1497/// <tr><td> `ExternalConstraints(const RooArgSet& )` <td> Include given external constraints to likelihood by multiplying them with the original likelihood.
1498/// <tr><td> `GlobalObservables(const RooArgSet&)` <td> Define the set of normalization observables to be used for the constraint terms.
1499/// If none are specified the constrained parameters are used.
1500/// <tr><td> `Offset(std::string const& mode)` <td> Likelihood offsetting mode. See createNLL(RooAbsData&, RooLinkedList const&).
1501/// <tr><td> `BatchMode(bool on)` <td> **Experimental** batch evaluation mode. This computes a batch of likelihood values at a time,
1502/// uses faster math functions and possibly auto vectorisation (this depends on the compiler flags).
1503/// Depending on hardware capabilities, the compiler flags and whether a batch evaluation function was
1504/// implemented for the PDFs of the model, likelihood computations are 2x to 10x faster.
1505/// The relative difference of the single log-likelihoods w.r.t. the legacy mode is usually better than 1.E-12,
1506/// and fit parameters usually agree to better than 1.E-6.
1507/// <tr><td> `IntegrateBins(double precision)` <td> In binned fits, integrate the PDF over the bins instead of using the probability density at the bin centre.
1508/// This can reduce the bias observed when fitting functions with high curvature to binned data.
1509/// - precision > 0: Activate bin integration everywhere. Use precision between 0.01 and 1.E-6, depending on binning.
1510/// Note that a low precision such as 0.01 might yield identical results to 1.E-4, since the integrator might reach 1.E-4 already in its first
1511/// integration step. If lower precision is desired (more speed), a RooBinSamplingPdf has to be created manually, and its integrator
1512/// has to be manipulated directly.
1513/// - precision = 0: Activate bin integration only for continuous PDFs fit to a RooDataHist.
1514/// - precision < 0: Deactivate.
1515/// \see RooBinSamplingPdf
1516///
1517/// <tr><th><th> Options to control flow of fit procedure
1518/// <tr><td> `Minimizer("<type>", "<algo>")` <td> Choose minimization package and optionally the algorithm to use. Default is MINUIT/MIGRAD through the RooMinimizer interface,
1519/// but others can be specified (through RooMinimizer interface).
1520/// <table>
1521/// <tr><th> Type <th> Algorithm
1522/// <tr><td> Minuit <td> migrad, simplex, minimize (=migrad+simplex), migradimproved (=migrad+improve)
1523/// <tr><td> Minuit2 <td> migrad, simplex, minimize, scan
1524/// <tr><td> GSLMultiMin <td> conjugatefr, conjugatepr, bfgs, bfgs2, steepestdescent
1525/// <tr><td> GSLSimAn <td> -
1526/// </table>
1527///
1528/// <tr><td> `InitialHesse(bool flag)` <td> Flag controls if HESSE before MIGRAD as well, off by default
1529/// <tr><td> `Optimize(bool flag)` <td> Activate constant term optimization of test statistic during minimization (on by default)
1530/// <tr><td> `Hesse(bool flag)` <td> Flag controls if HESSE is run after MIGRAD, on by default
1531/// <tr><td> `Minos(bool flag)` <td> Flag controls if MINOS is run after HESSE, off by default
1532/// <tr><td> `Minos(const RooArgSet& set)` <td> Only run MINOS on given subset of arguments
1533/// <tr><td> `Save(bool flag)` <td> Flag controls if RooFitResult object is produced and returned, off by default
1534/// <tr><td> `Strategy(Int_t flag)` <td> Set Minuit strategy (0 to 2, default is 1)
1535/// <tr><td> `MaxCalls(int n)` <td> Change maximum number of likelihood function calss from MINUIT (if `n <= 0`, the default of 500 * #%parameters is used)
1536/// <tr><td> `EvalErrorWall(bool flag=true)` <td> When parameters are in disallowed regions (e.g. PDF is negative), return very high value to fitter
1537/// to force it out of that region. This can, however, mean that the fitter gets lost in this region. If
1538/// this happens, try switching it off.
1539/// <tr><td> `RecoverFromUndefinedRegions(double strength)` <td> When PDF is invalid (e.g. parameter in undefined region), try to direct minimiser away from that region.
1540/// `strength` controls the magnitude of the penalty term. Leaving out this argument defaults to 10. Switch off with `strength = 0.`.
1541///
1542/// <tr><td> `SumW2Error(bool flag)` <td> Apply correction to errors and covariance matrix.
1543/// This uses two covariance matrices, one with the weights, the other with squared weights,
1544/// to obtain the correct errors for weighted likelihood fits. If this option is activated, the
1545/// corrected covariance matrix is calculated as \f$ V_\mathrm{corr} = V C^{-1} V \f$, where \f$ V \f$ is the original
1546/// covariance matrix and \f$ C \f$ is the inverse of the covariance matrix calculated using the
1547/// squared weights. This allows to switch between two interpretations of errors:
1548/// <table>
1549/// <tr><th> SumW2Error <th> Interpretation
1550/// <tr><td> true <td> The errors reflect the uncertainty of the Monte Carlo simulation.
1551/// Use this if you want to know how much accuracy you can get from the available Monte Carlo statistics.
1552///
1553/// **Example**: Simulation with 1000 events, the average weight is 0.1.
1554/// The errors are as big as if one fitted to 1000 events.
1555/// <tr><td> false <td> The errors reflect the errors of a dataset, which is as big as the sum of weights.
1556/// Use this if you want to know what statistical errors you would get if you had a dataset with as many
1557/// events as the (weighted) Monte Carlo simulation represents.
1558///
1559/// **Example** (Data as above):
1560/// The errors are as big as if one fitted to 100 events.
1561/// </table>
1562/// \note If the `SumW2Error` correction is enabled, the covariance matrix quality stored in the RooFitResult
1563/// object will be the minimum of the original covariance matrix quality and the quality of the covariance
1564/// matrix calculated with the squared weights.
1565/// <tr><td> `AsymptoticError()` <td> Use the asymptotically correct approach to estimate errors in the presence of weights.
1566/// This is slower but more accurate than `SumW2Error`. See also https://arxiv.org/abs/1911.01303).
1567/// <tr><td> `PrefitDataFraction(double fraction)`
1568/// <td> Runs a prefit on a small dataset of size fraction*(actual data size). This can speed up fits
1569/// by finding good starting values for the parameters for the actual fit.
1570/// \warning Prefitting may give bad results when used in binned analysis.
1571///
1572/// <tr><th><th> Options to control informational output
1573/// <tr><td> `Verbose(bool flag)` <td> Flag controls if verbose output is printed (NLL, parameter changes during fit).
1574/// <tr><td> `Timer(bool flag)` <td> Time CPU and wall clock consumption of fit steps, off by default.
1575/// <tr><td> `PrintLevel(Int_t level)` <td> Set Minuit print level (-1 to 3, default is 1). At -1 all RooFit informational messages are suppressed as well.
1576/// See RooMinimizer::PrintLevel for the meaning of the levels.
1577/// <tr><td> `Warnings(bool flag)` <td> Enable or disable MINUIT warnings (enabled by default)
1578/// <tr><td> `PrintEvalErrors(Int_t numErr)` <td> Control number of p.d.f evaluation errors printed per likelihood evaluation.
1579/// A negative value suppresses output completely, a zero value will only print the error count per p.d.f component,
1580/// a positive value will print details of each error up to `numErr` messages per p.d.f component.
1581/// <tr><td> `ModularL(bool flag)` <td> Enable or disable modular likelihoods, which will become the default in a future release.
1582/// This does not change any user-facing code, but only enables a different likelihood class in the back-end.
1583/// This option is forced to true in case parallelization is requested with the `Parallelize` named argument.
1584/// Note that it is currently not recommended to use modular likelihoods without any parallelization enabled, since
1585/// some features such as offsetting might not yet work in this case.
1586/// <tr><td> `Parallelize(Int_t nWorkers)` <td> Control global parallelization settings. Arguments 1 and above enable the use of RooFit's parallel minimization
1587/// backend and uses the number given as the number of workers to use in the parallelization. -1 also enables
1588/// RooFit's parallel minimization backend, and sets the number of workers to the number of available processes.
1589/// 0 disables this feature.
1590/// <tr><td> `ParallelGradientOptions(bool enable=true, int orderStrategy=0, int chainFactor=1)` <td> **Experimental** Control gradient parallelization settings. The first argument
1591/// only disables or enables gradient parallelization, this is on by default.
1592/// The second argument determines the internal partial derivative calculation
1593/// ordering strategy. The third argument determines the number of partial
1594/// derivatives that are executed per task package on each worker.
1595/// <tr><td> `ParallelDescentOptions(bool enable=false, int splitStrategy=0, int numSplits=4)` <td> **Experimental** Control settings related to the parallelization of likelihoods
1596/// outside of the gradient calculation but in the minimization, most prominently
1597/// in the linesearch step. The first argument this disables or enables likelihood
1598/// parallelization. The second argument determines whether to split the task batches
1599/// per event or per likelihood component. And the third argument how many events or
1600/// respectively components to include in each batch.
1601/// <tr><td> `TimingAnalysis(bool flag)` <td> **Experimental** log timings. This feature logs timings with NewStyle likelihoods on multiple processes simultaneously
1602/// and outputs the timings at the end of a run to json log files, which can be analyzed with the
1603/// `RooFit::MultiProcess::HeatmapAnalyzer`. Only works with simultaneous likelihoods.
1604/// </table>
1605///
1606
1608{
1609 // Select the pdf-specific commands
1610 RooCmdConfig pc(Form("RooAbsPdf::fitTo(%s)",GetName())) ;
1611
1612 RooLinkedList fitCmdList(cmdList) ;
1613 std::string nllCmdListString = "ProjectedObservables,Extended,Range,"
1614 "RangeWithName,SumCoefRange,NumCPU,SplitRange,Constrained,Constrain,ExternalConstraints,"
1615 "CloneData,GlobalObservables,GlobalObservablesSource,GlobalObservablesTag,"
1616 "BatchMode,IntegrateBins,ModularL";
1617
1618 if (((RooCmdArg*)cmdList.FindObject("ModularL")) && !((RooCmdArg*)cmdList.FindObject("ModularL"))->getInt(0))
1619 nllCmdListString += ",OffsetLikelihood";
1620
1621 RooLinkedList nllCmdList = pc.filterCmdList(fitCmdList, nllCmdListString.c_str());
1622
1623 // Default-initialized instance of MinimizerConfig to get the default
1624 // minimizer parameter values.
1625 MinimizerConfig minimizerDefaults;
1626
1627 pc.defineDouble("prefit", "Prefit",0,0);
1628 pc.defineDouble("RecoverFromUndefinedRegions", "RecoverFromUndefinedRegions",0,minimizerDefaults.recoverFromNaN);
1629 pc.defineInt("optConst","Optimize",0,minimizerDefaults.optConst) ;
1630 pc.defineInt("verbose","Verbose",0,minimizerDefaults.verbose) ;
1631 pc.defineInt("doSave","Save",0,minimizerDefaults.doSave) ;
1632 pc.defineInt("doTimer","Timer",0,minimizerDefaults.doTimer) ;
1633 pc.defineInt("printLevel","PrintLevel",0,minimizerDefaults.printLevel) ;
1634 pc.defineInt("strat","Strategy",0,minimizerDefaults.strat) ;
1635 pc.defineInt("initHesse","InitialHesse",0,minimizerDefaults.initHesse) ;
1636 pc.defineInt("hesse","Hesse",0,minimizerDefaults.hesse) ;
1637 pc.defineInt("minos","Minos",0,minimizerDefaults.minos) ;
1638 pc.defineInt("numee","PrintEvalErrors",0,minimizerDefaults.numee) ;
1639 pc.defineInt("doEEWall","EvalErrorWall",0,minimizerDefaults.doEEWall) ;
1640 pc.defineInt("doWarn","Warnings",0,minimizerDefaults.doWarn) ;
1641 pc.defineInt("doSumW2","SumW2Error",0,minimizerDefaults.doSumW2) ;
1642 pc.defineInt("doAsymptoticError","AsymptoticError",0,minimizerDefaults.doAsymptotic) ;
1643 pc.defineInt("maxCalls","MaxCalls",0,minimizerDefaults.maxCalls);
1644 pc.defineInt("doOffset","OffsetLikelihood",0,0) ;
1645 pc.defineInt("parallelize", "Parallelize", 0, 0); // Three parallelize arguments
1646 pc.defineInt("enableParallelGradient", "ParallelGradientOptions", 0, 0);
1647 pc.defineInt("enableParallelDescent", "ParallelDescentOptions", 0, 0);
1648 pc.defineInt("timingAnalysis", "TimingAnalysis", 0, 0);
1649 pc.defineString("mintype","Minimizer",0,minimizerDefaults.minType.c_str()) ;
1650 pc.defineString("minalg","Minimizer",1,minimizerDefaults.minAlg.c_str()) ;
1651 pc.defineSet("minosSet","Minos",0,minimizerDefaults.minosSet) ;
1652 pc.defineMutex("Range","RangeWithName") ;
1653
1654 // Process and check varargs
1655 pc.process(fitCmdList) ;
1656 if (!pc.ok(true)) {
1657 return 0 ;
1658 }
1659
1660 // TimingAnalysis works only for RooSimultaneous.
1661 if (pc.getInt("timingAnalysis") && !this->InheritsFrom("RooSimultaneous") ) {
1662 coutW(Minimization) << "The timingAnalysis feature was built for minimization with RooSimulteneous "
1663 "and is not implemented for other PDF's. Please create a RooSimultenous to "
1664 "enable this feature." << endl;
1665 }
1666
1667 // Decode command line arguments
1668 double prefit = pc.getDouble("prefit");
1669 Int_t optConst = pc.getInt("optConst") ;
1670
1671 if (optConst > 1) {
1672 // optConst >= 2 is pre-computating values, which are never used when
1673 // the batchMode is on. This just wastes time.
1674
1675 RooCmdConfig conf(Form("RooAbsPdf::fitTo(%s)", GetName()));
1676 conf.defineInt("BatchMode","BatchMode",0,0);
1677 conf.allowUndefined(true);
1678 conf.process(nllCmdList);
1679 if (conf.getInt("BatchMode") != 0) {
1680 optConst = 1;
1681 }
1682 }
1683
1684 if (prefit != 0) {
1685 size_t nEvents = static_cast<size_t>(prefit*data.numEntries());
1686 if (prefit > 0.5 || nEvents < 100) {
1687 coutW(InputArguments) << "PrefitDataFraction should be in suitable range."
1688 << "With the current PrefitDataFraction=" << prefit
1689 << ", the number of events would be " << nEvents<< " out of "
1690 << data.numEntries() << ". Skipping prefit..." << endl;
1691 }
1692 else {
1693 size_t step = data.numEntries()/nEvents;
1694
1695 RooDataSet tiny("tiny", "tiny", *data.get(),
1696 data.isWeighted() ? RooFit::WeightVar() : RooCmdArg());
1697
1698 for (int i=0; i<data.numEntries(); i+=step)
1699 {
1700 const RooArgSet *event = data.get(i);
1701 tiny.add(*event, data.weight());
1702 }
1703 RooLinkedList tinyCmdList(cmdList) ;
1704 pc.filterCmdList(tinyCmdList,"Prefit,Hesse,Minos,Verbose,Save,Timer");
1705 RooCmdArg hesse_option = RooFit::Hesse(false);
1706 RooCmdArg print_option = RooFit::PrintLevel(-1);
1707
1708 tinyCmdList.Add(&hesse_option);
1709 tinyCmdList.Add(&print_option);
1710
1711 fitTo(tiny,tinyCmdList);
1712 }
1713 }
1714
1715 RooCmdArg modularL_option;
1716 if (pc.getInt("parallelize") != 0 || pc.getInt("enableParallelGradient") || pc.getInt("enableParallelDescent")) {
1717 // Set to new style likelihood if parallelization is requested
1718 modularL_option = RooFit::ModularL(true);
1719 nllCmdList.Add(&modularL_option);
1720 }
1721
1722 std::unique_ptr<RooAbsReal> nll{createNLL(data,nllCmdList)};
1723
1724 MinimizerConfig cfg;
1725 cfg.recoverFromNaN = pc.getDouble("RecoverFromUndefinedRegions");
1726 cfg.optConst = optConst;
1727 cfg.verbose = pc.getInt("verbose");
1728 cfg.doSave = pc.getInt("doSave");
1729 cfg.doTimer = pc.getInt("doTimer");
1730 cfg.printLevel = pc.getInt("printLevel");
1731 cfg.strat = pc.getInt("strat");
1732 cfg.initHesse = pc.getInt("initHesse");
1733 cfg.hesse = pc.getInt("hesse");
1734 cfg.minos = pc.getInt("minos");
1735 cfg.numee = pc.getInt("numee");
1736 cfg.doEEWall = pc.getInt("doEEWall");
1737 cfg.doWarn = pc.getInt("doWarn");
1738 cfg.doSumW2 = pc.getInt("doSumW2");
1739 cfg.doAsymptotic = pc.getInt("doAsymptoticError");
1740 cfg.maxCalls = pc.getInt("maxCalls");
1741 cfg.minosSet = pc.getSet("minosSet");
1742 cfg.minType = pc.getString("mintype","");
1743 cfg.minAlg = pc.getString("minalg","minuit");
1744 cfg.doOffset = pc.getInt("doOffset");
1745 cfg.parallelize = pc.getInt("parallelize");
1746 cfg.enableParallelGradient = pc.getInt("enableParallelGradient");
1747 cfg.enableParallelDescent = pc.getInt("enableParallelDescent");
1748 cfg.timingAnalysis = pc.getInt("timingAnalysis");
1749 return RooFit::Detail::owningPtr(minimizeNLL(*nll, data, cfg));
1750}
1751
1752
1753
1754////////////////////////////////////////////////////////////////////////////////
1755/// Calls RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList) and returns fit result.
1756
1758{
1759 // Select the pdf-specific commands
1760 RooCmdConfig pc(Form("RooAbsPdf::chi2FitTo(%s)",GetName())) ;
1761
1762 // Pull arguments to be passed to chi2 construction from list
1763 RooLinkedList fitCmdList(cmdList) ;
1764 RooLinkedList chi2CmdList = pc.filterCmdList(fitCmdList,"Range,RangeWithName,NumCPU,Optimize,ProjectedObservables,AddCoefRange,SplitRange,DataError,Extended,IntegrateBins") ;
1765
1766 std::unique_ptr<RooAbsReal> chi2{createChi2(data,chi2CmdList)};
1767 return chi2FitDriver(*chi2,fitCmdList) ;
1768}
1769
1770
1771
1772
1773////////////////////////////////////////////////////////////////////////////////
1774/// Create a \f$ \chi^2 \f$ from a histogram and this function.
1775///
1776/// Options to control construction of the chi-square
1777/// ------------------------------------------
1778/// The list of supported command arguments is given in the documentation for
1779/// RooChi2Var::RooChi2Var(const char *name, const char* title, RooAbsReal& func, RooDataHist& hdata, const RooCmdArg&,const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&,const RooCmdArg&, const RooCmdArg&,const RooCmdArg&,const RooCmdArg&).
1780
1782 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
1783 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
1784{
1785 RooLinkedList cmdList ;
1786 cmdList.Add((TObject*)&arg1) ; cmdList.Add((TObject*)&arg2) ;
1787 cmdList.Add((TObject*)&arg3) ; cmdList.Add((TObject*)&arg4) ;
1788 cmdList.Add((TObject*)&arg5) ; cmdList.Add((TObject*)&arg6) ;
1789 cmdList.Add((TObject*)&arg7) ; cmdList.Add((TObject*)&arg8) ;
1790
1791 RooCmdConfig pc(Form("RooAbsPdf::createChi2(%s)",GetName())) ;
1792 pc.defineString("rangeName","RangeWithName",0,"",true) ;
1793 pc.allowUndefined(true) ;
1794 pc.process(cmdList) ;
1795 if (!pc.ok(true)) {
1796 return 0 ;
1797 }
1798 const char* rangeName = pc.getString("rangeName",0,true) ;
1799
1800 // Construct Chi2
1802 RooAbsReal* chi2 ;
1803 string baseName = Form("chi2_%s_%s",GetName(),data.GetName()) ;
1804
1805 // Clear possible range attributes from previous fits.
1806 removeStringAttribute("fitrange");
1807
1808 if (!rangeName || strchr(rangeName,',')==0) {
1809 // Simple case: default range, or single restricted range
1810
1811 chi2 = new RooChi2Var(baseName.c_str(),baseName.c_str(),*this,data,arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8) ;
1812
1813 } else {
1814
1815 // Find which argument is RangeWithName
1816 const RooCmdArg* rarg(0) ;
1817 string rcmd = "RangeWithName" ;
1818 if (arg1.GetName()==rcmd) rarg = &arg1 ;
1819 if (arg2.GetName()==rcmd) rarg = &arg2 ;
1820 if (arg3.GetName()==rcmd) rarg = &arg3 ;
1821 if (arg4.GetName()==rcmd) rarg = &arg4 ;
1822 if (arg5.GetName()==rcmd) rarg = &arg5 ;
1823 if (arg6.GetName()==rcmd) rarg = &arg6 ;
1824 if (arg7.GetName()==rcmd) rarg = &arg7 ;
1825 if (arg8.GetName()==rcmd) rarg = &arg8 ;
1826
1827 // Composite case: multiple ranges
1828 RooArgList chi2List ;
1829 for (std::string& token : ROOT::Split(rangeName, ",")) {
1830 RooCmdArg subRangeCmd = RooFit::Range(token.c_str()) ;
1831 // Construct chi2 while substituting original RangeWithName argument with subrange argument created above
1832 RooAbsReal* chi2Comp = new RooChi2Var(Form("%s_%s", baseName.c_str(), token.c_str()), "chi^2", *this, data,
1833 &arg1==rarg?subRangeCmd:arg1,&arg2==rarg?subRangeCmd:arg2,
1834 &arg3==rarg?subRangeCmd:arg3,&arg4==rarg?subRangeCmd:arg4,
1835 &arg5==rarg?subRangeCmd:arg5,&arg6==rarg?subRangeCmd:arg6,
1836 &arg7==rarg?subRangeCmd:arg7,&arg8==rarg?subRangeCmd:arg8) ;
1837 chi2List.add(*chi2Comp) ;
1838 }
1839 chi2 = new RooAddition(baseName.c_str(),"chi^2",chi2List,true) ;
1840 }
1842
1843
1844 return chi2 ;
1845}
1846
1847
1848
1849
1850////////////////////////////////////////////////////////////////////////////////
1851/// Argument-list version of RooAbsPdf::createChi2()
1852
1854{
1855 // Select the pdf-specific commands
1856 RooCmdConfig pc(Form("RooAbsPdf::createChi2(%s)",GetName())) ;
1857
1858 pc.defineInt("integrate","Integrate",0,0) ;
1859 pc.defineObject("yvar","YVar",0,0) ;
1860
1861 // Process and check varargs
1862 pc.process(cmdList) ;
1863 if (!pc.ok(true)) {
1864 return 0 ;
1865 }
1866
1867 // Decode command line arguments
1868 bool integrate = pc.getInt("integrate") ;
1869 RooRealVar* yvar = (RooRealVar*) pc.getObject("yvar") ;
1870
1871 string name = Form("chi2_%s_%s",GetName(),data.GetName()) ;
1872
1873 if (yvar) {
1874 return new RooXYChi2Var(name.c_str(),name.c_str(),*this,data,*yvar,integrate) ;
1875 } else {
1876 return new RooXYChi2Var(name.c_str(),name.c_str(),*this,data,integrate) ;
1877 }
1878}
1879
1880
1881
1882
1883////////////////////////////////////////////////////////////////////////////////
1884/// Print value of p.d.f, also print normalization integral that was last used, if any
1885
1886void RooAbsPdf::printValue(ostream& os) const
1887{
1888 // silent warning messages coming when evaluating a RooAddPdf without a normalization set
1890
1891 getVal() ;
1892
1893 if (_norm) {
1894 os << getVal() << "/" << _norm->getVal() ;
1895 } else {
1896 os << getVal();
1897 }
1898}
1899
1900
1901
1902////////////////////////////////////////////////////////////////////////////////
1903/// Print multi line detailed information of this RooAbsPdf
1904
1905void RooAbsPdf::printMultiline(ostream& os, Int_t contents, bool verbose, TString indent) const
1906{
1907 RooAbsReal::printMultiline(os,contents,verbose,indent);
1908 os << indent << "--- RooAbsPdf ---" << endl;
1909 os << indent << "Cached value = " << _value << endl ;
1910 if (_norm) {
1911 os << indent << " Normalization integral: " << endl ;
1912 auto moreIndent = std::string(indent.Data()) + " " ;
1913 _norm->printStream(os,kName|kAddress|kTitle|kValue|kArgs,kSingleLine,moreIndent.c_str()) ;
1914 }
1915}
1916
1917
1918
1919////////////////////////////////////////////////////////////////////////////////
1920/// Return a binned generator context
1921
1923{
1924 return new RooBinnedGenContext(*this,vars,0,0,verbose) ;
1925}
1926
1927
1928////////////////////////////////////////////////////////////////////////////////
1929/// Interface function to create a generator context from a p.d.f. This default
1930/// implementation returns a 'standard' context that works for any p.d.f
1931
1933 const RooArgSet* auxProto, bool verbose) const
1934{
1935 return new RooGenContext(*this,vars,prototype,auxProto,verbose) ;
1936}
1937
1938
1939////////////////////////////////////////////////////////////////////////////////
1940
1941RooAbsGenContext* RooAbsPdf::autoGenContext(const RooArgSet &vars, const RooDataSet* prototype, const RooArgSet* auxProto,
1942 bool verbose, bool autoBinned, const char* binnedTag) const
1943{
1944 if (prototype || (auxProto && auxProto->getSize()>0)) {
1945 return genContext(vars,prototype,auxProto,verbose);
1946 }
1947
1948 RooAbsGenContext *context(0) ;
1949 if ( (autoBinned && isBinnedDistribution(vars)) || ( binnedTag && strlen(binnedTag) && (getAttribute(binnedTag)||string(binnedTag)=="*"))) {
1950 context = binnedGenContext(vars,verbose) ;
1951 } else {
1952 context= genContext(vars,0,0,verbose);
1953 }
1954 return context ;
1955}
1956
1957
1958
1959////////////////////////////////////////////////////////////////////////////////
1960/// Generate a new dataset containing the specified variables with events sampled from our distribution.
1961/// Generate the specified number of events or expectedEvents() if not specified.
1962/// \param[in] whatVars Choose variables in which to generate events. Variables not listed here will remain
1963/// constant and not be used for event generation.
1964/// \param[in] arg1,arg2,arg3,arg4,arg5,arg6 Optional RooCmdArg() to change behaviour of generate().
1965/// \return RooDataSet *, owned by caller.
1966///
1967/// Any variables of this PDF that are not in whatVars will use their
1968/// current values and be treated as fixed parameters. Returns zero
1969/// in case of an error.
1970///
1971/// <table>
1972/// <tr><th> Type of CmdArg <th> Effect on generate
1973/// <tr><td> `Name(const char* name)` <td> Name of the output dataset
1974/// <tr><td> `Verbose(bool flag)` <td> Print informational messages during event generation
1975/// <tr><td> `NumEvents(int nevt)` <td> Generate specified number of events
1976/// <tr><td> `Extended()` <td> If no number of events to be generated is given,
1977/// use expected number of events from extended likelihood term.
1978/// This evidently only works for extended PDFs.
1979/// <tr><td> `GenBinned(const char* tag)` <td> Use binned generation for all component pdfs that have 'setAttribute(tag)' set
1980/// <tr><td> `AutoBinned(bool flag)` <td> Automatically deploy binned generation for binned distributions (e.g. RooHistPdf, sums and products of
1981/// RooHistPdfs etc)
1982/// \note Datasets that are generated in binned mode are returned as weighted unbinned datasets. This means that
1983/// for each bin, there will be one event in the dataset with a weight corresponding to the (possibly randomised) bin content.
1984///
1985///
1986/// <tr><td> `AllBinned()` <td> As above, but for all components.
1987/// \note The notion of components is only meaningful for simultaneous PDFs
1988/// as binned generation is always executed at the top-level node for a regular
1989/// PDF, so for those it only mattes that the top-level node is tagged.
1990///
1991/// <tr><td> ProtoData(const RooDataSet& data, bool randOrder)
1992/// <td> Use specified dataset as prototype dataset. If randOrder in ProtoData() is set to true,
1993/// the order of the events in the dataset will be read in a random order if the requested
1994/// number of events to be generated does not match the number of events in the prototype dataset.
1995/// \note If ProtoData() is used, the specified existing dataset as a prototype: the new dataset will contain
1996/// the same number of events as the prototype (unless otherwise specified), and any prototype variables not in
1997/// whatVars will be copied into the new dataset for each generated event and also used to set our PDF parameters.
1998/// The user can specify a number of events to generate that will override the default. The result is a
1999/// copy of the prototype dataset with only variables in whatVars randomized. Variables in whatVars that
2000/// are not in the prototype will be added as new columns to the generated dataset.
2001///
2002/// </table>
2003///
2004/// #### Accessing the underlying event generator
2005/// Depending on the fit model (if it is difficult to sample), it may be necessary to change generator settings.
2006/// For the default generator (RooFoamGenerator), the number of samples or cells could be increased by e.g. using
2007/// myPdf->specialGeneratorConfig()->getConfigSection("RooFoamGenerator").setRealValue("nSample",1e4);
2008///
2009/// The foam generator e.g. has the following config options:
2010/// - nCell[123N]D
2011/// - nSample
2012/// - chatLevel
2013/// \see rf902_numgenconfig.C
2014
2015RooDataSet *RooAbsPdf::generate(const RooArgSet& whatVars, const RooCmdArg& arg1,const RooCmdArg& arg2,
2016 const RooCmdArg& arg3,const RooCmdArg& arg4, const RooCmdArg& arg5,const RooCmdArg& arg6)
2017{
2018 // Select the pdf-specific commands
2019 RooCmdConfig pc(Form("RooAbsPdf::generate(%s)",GetName())) ;
2020 pc.defineObject("proto","PrototypeData",0,0) ;
2021 pc.defineString("dsetName","Name",0,"") ;
2022 pc.defineInt("randProto","PrototypeData",0,0) ;
2023 pc.defineInt("resampleProto","PrototypeData",1,0) ;
2024 pc.defineInt("verbose","Verbose",0,0) ;
2025 pc.defineInt("extended","Extended",0,0) ;
2026 pc.defineInt("nEvents","NumEvents",0,0) ;
2027 pc.defineInt("autoBinned","AutoBinned",0,1) ;
2028 pc.defineInt("expectedData","ExpectedData",0,0) ;
2029 pc.defineDouble("nEventsD","NumEventsD",0,-1.) ;
2030 pc.defineString("binnedTag","GenBinned",0,"") ;
2031 pc.defineMutex("GenBinned","ProtoData") ;
2032 pc.defineMutex("Extended", "NumEvents");
2033
2034 // Process and check varargs
2035 pc.process(arg1,arg2,arg3,arg4,arg5,arg6) ;
2036 if (!pc.ok(true)) {
2037 return 0 ;
2038 }
2039
2040 // Decode command line arguments
2041 RooDataSet* protoData = static_cast<RooDataSet*>(pc.getObject("proto",0)) ;
2042 const char* dsetName = pc.getString("dsetName") ;
2043 bool verbose = pc.getInt("verbose") ;
2044 bool randProto = pc.getInt("randProto") ;
2045 bool resampleProto = pc.getInt("resampleProto") ;
2046 bool extended = pc.getInt("extended") ;
2047 bool autoBinned = pc.getInt("autoBinned") ;
2048 const char* binnedTag = pc.getString("binnedTag") ;
2049 Int_t nEventsI = pc.getInt("nEvents") ;
2050 double nEventsD = pc.getInt("nEventsD") ;
2051 //bool verbose = pc.getInt("verbose") ;
2052 bool expectedData = pc.getInt("expectedData") ;
2053
2054 double nEvents = (nEventsD>0) ? nEventsD : double(nEventsI);
2055
2056 // Force binned mode for expected data mode
2057 if (expectedData) {
2058 binnedTag="*" ;
2059 }
2060
2061 if (extended) {
2062 if (nEvents == 0) nEvents = expectedEvents(&whatVars);
2063 } else if (nEvents==0) {
2064 cxcoutI(Generation) << "No number of events specified , number of events generated is "
2065 << GetName() << "::expectedEvents() = " << expectedEvents(&whatVars)<< endl ;
2066 }
2067
2068 if (extended && protoData && !randProto) {
2069 cxcoutI(Generation) << "WARNING Using generator option Extended() (Poisson distribution of #events) together "
2070 << "with a prototype dataset implies incomplete sampling or oversampling of proto data. "
2071 << "Set randomize flag in ProtoData() option to randomize prototype dataset order and thus "
2072 << "to randomize the set of over/undersampled prototype events for each generation cycle." << endl ;
2073 }
2074
2075
2076 // Forward to appropriate implementation
2077 RooDataSet* data ;
2078 if (protoData) {
2079 data = generate(whatVars,*protoData,Int_t(nEvents),verbose,randProto,resampleProto) ;
2080 } else {
2081 data = generate(whatVars,nEvents,verbose,autoBinned,binnedTag,expectedData, extended) ;
2082 }
2083
2084 // Rename dataset to given name if supplied
2085 if (dsetName && strlen(dsetName)>0) {
2086 data->SetName(dsetName) ;
2087 }
2088
2089 return data ;
2090}
2091
2092
2093
2094
2095
2096
2097////////////////////////////////////////////////////////////////////////////////
2098/// \note This method does not perform any generation. To generate according to generations specification call RooAbsPdf::generate(RooAbsPdf::GenSpec&) const
2099///
2100/// Details copied from RooAbsPdf::generate():
2101/// --------------------------------------------
2102/// \copydetails RooAbsPdf::generate(const RooArgSet&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&)
2103
2105 const RooCmdArg& arg1,const RooCmdArg& arg2,
2106 const RooCmdArg& arg3,const RooCmdArg& arg4,
2107 const RooCmdArg& arg5,const RooCmdArg& arg6)
2108{
2109
2110 // Select the pdf-specific commands
2111 RooCmdConfig pc(Form("RooAbsPdf::generate(%s)",GetName())) ;
2112 pc.defineObject("proto","PrototypeData",0,0) ;
2113 pc.defineString("dsetName","Name",0,"") ;
2114 pc.defineInt("randProto","PrototypeData",0,0) ;
2115 pc.defineInt("resampleProto","PrototypeData",1,0) ;
2116 pc.defineInt("verbose","Verbose",0,0) ;
2117 pc.defineInt("extended","Extended",0,0) ;
2118 pc.defineInt("nEvents","NumEvents",0,0) ;
2119 pc.defineInt("autoBinned","AutoBinned",0,1) ;
2120 pc.defineString("binnedTag","GenBinned",0,"") ;
2121 pc.defineMutex("GenBinned","ProtoData") ;
2122
2123
2124 // Process and check varargs
2125 pc.process(arg1,arg2,arg3,arg4,arg5,arg6) ;
2126 if (!pc.ok(true)) {
2127 return 0 ;
2128 }
2129
2130 // Decode command line arguments
2131 RooDataSet* protoData = static_cast<RooDataSet*>(pc.getObject("proto",0)) ;
2132 const char* dsetName = pc.getString("dsetName") ;
2133 Int_t nEvents = pc.getInt("nEvents") ;
2134 bool verbose = pc.getInt("verbose") ;
2135 bool randProto = pc.getInt("randProto") ;
2136 bool resampleProto = pc.getInt("resampleProto") ;
2137 bool extended = pc.getInt("extended") ;
2138 bool autoBinned = pc.getInt("autoBinned") ;
2139 const char* binnedTag = pc.getString("binnedTag") ;
2140
2141 RooAbsGenContext* cx = autoGenContext(whatVars,protoData,0,verbose,autoBinned,binnedTag) ;
2142
2143 return new GenSpec(cx,whatVars,protoData,nEvents,extended,randProto,resampleProto,dsetName) ;
2144}
2145
2146
2147////////////////////////////////////////////////////////////////////////////////
2148/// If many identical generation requests
2149/// are needed, e.g. in toy MC studies, it is more efficient to use the prepareMultiGen()/generate()
2150/// combination than calling the standard generate() multiple times as
2151/// initialization overhead is only incurred once.
2152
2154{
2155 //Int_t nEvt = spec._extended ? RooRandom::randomGenerator()->Poisson(spec._nGen) : spec._nGen ;
2156 //Int_t nEvt = spec._extended ? RooRandom::randomGenerator()->Poisson(spec._nGen==0?expectedEvents(spec._whatVars):spec._nGen) : spec._nGen ;
2157 //Int_t nEvt = spec._nGen == 0 ? RooRandom::randomGenerator()->Poisson(expectedEvents(spec._whatVars)) : spec._nGen;
2158
2159 double nEvt = spec._nGen == 0 ? expectedEvents(spec._whatVars) : spec._nGen;
2160
2161 RooDataSet* ret = generate(*spec._genContext,spec._whatVars,spec._protoData, nEvt,false,spec._randProto,spec._resampleProto,
2162 spec._init,spec._extended) ;
2163 spec._init = true ;
2164 return ret ;
2165}
2166
2167
2168
2169
2170
2171////////////////////////////////////////////////////////////////////////////////
2172/// Generate a new dataset containing the specified variables with
2173/// events sampled from our distribution.
2174///
2175/// \param[in] whatVars Generate a dataset with the variables (and categories) in this set.
2176/// Any variables of this PDF that are not in `whatVars` will use their
2177/// current values and be treated as fixed parameters.
2178/// \param[in] nEvents Generate the specified number of events or else try to use
2179/// expectedEvents() if nEvents <= 0 (default).
2180/// \param[in] verbose Show which generator strategies are being used.
2181/// \param[in] autoBinned If original distribution is binned, return bin centers and randomise weights
2182/// instead of generating single events.
2183/// \param[in] binnedTag
2184/// \param[in] expectedData Call setExpectedData on the genContext.
2185/// \param[in] extended Randomise number of events generated according to Poisson(nEvents). Only useful
2186/// if PDF is extended.
2187/// \return New dataset. Returns zero in case of an error. The caller takes ownership of the returned
2188/// dataset.
2189
2190RooDataSet *RooAbsPdf::generate(const RooArgSet &whatVars, double nEvents, bool verbose, bool autoBinned, const char* binnedTag, bool expectedData, bool extended) const
2191{
2192 if (nEvents==0 && extendMode()==CanNotBeExtended) {
2193 return new RooDataSet("emptyData","emptyData",whatVars) ;
2194 }
2195
2196 // Request for binned generation
2197 std::unique_ptr<RooAbsGenContext> context{autoGenContext(whatVars,0,0,verbose,autoBinned,binnedTag)};
2198 if (expectedData) {
2199 context->setExpectedData(true) ;
2200 }
2201
2202 RooDataSet *generated = 0;
2203 if(0 != context && context->isValid()) {
2204 generated= context->generate(nEvents, false, extended);
2205 }
2206 else {
2207 coutE(Generation) << "RooAbsPdf::generate(" << GetName() << ") cannot create a valid context" << endl;
2208 }
2209 return generated;
2210}
2211
2212
2213
2214
2215////////////////////////////////////////////////////////////////////////////////
2216/// Internal method
2217
2218RooDataSet *RooAbsPdf::generate(RooAbsGenContext& context, const RooArgSet &whatVars, const RooDataSet *prototype,
2219 double nEvents, bool /*verbose*/, bool randProtoOrder, bool resampleProto,
2220 bool skipInit, bool extended) const
2221{
2222 if (nEvents==0 && (prototype==0 || prototype->numEntries()==0)) {
2223 return new RooDataSet("emptyData","emptyData",whatVars) ;
2224 }
2225
2226 RooDataSet *generated = 0;
2227
2228 // Resampling implies reshuffling in the implementation
2229 if (resampleProto) {
2230 randProtoOrder=true ;
2231 }
2232
2233 if (randProtoOrder && prototype && prototype->numEntries()!=nEvents) {
2234 coutI(Generation) << "RooAbsPdf::generate (Re)randomizing event order in prototype dataset (Nevt=" << nEvents << ")" << endl ;
2235 Int_t* newOrder = randomizeProtoOrder(prototype->numEntries(),Int_t(nEvents),resampleProto) ;
2236 context.setProtoDataOrder(newOrder) ;
2237 delete[] newOrder ;
2238 }
2239
2240 if(context.isValid()) {
2241 generated= context.generate(nEvents,skipInit,extended);
2242 }
2243 else {
2244 coutE(Generation) << "RooAbsPdf::generate(" << GetName() << ") do not have a valid generator context" << endl;
2245 }
2246 return generated;
2247}
2248
2249
2250
2251
2252////////////////////////////////////////////////////////////////////////////////
2253/// Generate a new dataset using a prototype dataset as a model,
2254/// with values of the variables in `whatVars` sampled from our distribution.
2255///
2256/// \param[in] whatVars Generate for these variables.
2257/// \param[in] prototype Use this dataset
2258/// as a prototype: the new dataset will contain the same number of
2259/// events as the prototype (by default), and any prototype variables not in
2260/// whatVars will be copied into the new dataset for each generated
2261/// event and also used to set our PDF parameters. The user can specify a
2262/// number of events to generate that will override the default. The result is a
2263/// copy of the prototype dataset with only variables in whatVars
2264/// randomized. Variables in whatVars that are not in the prototype
2265/// will be added as new columns to the generated dataset.
2266/// \param[in] nEvents Number of events to generate. Defaults to 0, which means number
2267/// of event in prototype dataset.
2268/// \param[in] verbose Show which generator strategies are being used.
2269/// \param[in] randProtoOrder Randomise order of retrieval of events from proto dataset.
2270/// \param[in] resampleProto Resample from the proto dataset.
2271/// \return The new dataset. Returns zero in case of an error. The caller takes ownership of the
2272/// returned dataset.
2273
2274RooDataSet *RooAbsPdf::generate(const RooArgSet &whatVars, const RooDataSet& prototype,
2275 Int_t nEvents, bool verbose, bool randProtoOrder, bool resampleProto) const
2276{
2277 std::unique_ptr<RooAbsGenContext> context{genContext(whatVars,&prototype,0,verbose)};
2278 if (context) {
2279 RooDataSet* data = generate(*context,whatVars,&prototype,nEvents,verbose,randProtoOrder,resampleProto) ;
2280 return data ;
2281 } else {
2282 coutE(Generation) << "RooAbsPdf::generate(" << GetName() << ") ERROR creating generator context" << endl ;
2283 return nullptr;
2284 }
2285}
2286
2287
2288
2289////////////////////////////////////////////////////////////////////////////////
2290/// Return lookup table with randomized order for nProto prototype events.
2291
2292Int_t* RooAbsPdf::randomizeProtoOrder(Int_t nProto, Int_t, bool resampleProto) const
2293{
2294 // Make output list
2295 Int_t* lut = new Int_t[nProto] ;
2296
2297 // Randomly sample input list into output list
2298 if (!resampleProto) {
2299 // In this mode, randomization is a strict reshuffle of the order
2300 std::iota(lut, lut + nProto, 0); // fill the vector with 0 to nProto - 1
2301 // Shuffle code taken from https://en.cppreference.com/w/cpp/algorithm/random_shuffle.
2302 // The std::random_shuffle function was deprecated in C++17. We could have
2303 // used std::shuffle instead, but this is not straight-forward to use with
2304 // RooRandom::integer() and we didn't want to change the random number
2305 // generator. It might cause unwanted effects like reproducibility problems.
2306 for (int i = nProto-1; i > 0; --i) {
2307 std::swap(lut[i], lut[RooRandom::integer(i+1)]);
2308 }
2309 } else {
2310 // In this mode, we resample, i.e. events can be used more than once
2311 std::generate(lut, lut + nProto, [&]{ return RooRandom::integer(nProto); });
2312 }
2313
2314
2315 return lut ;
2316}
2317
2318
2319
2320////////////////////////////////////////////////////////////////////////////////
2321/// Load generatedVars with the subset of directVars that we can generate events for,
2322/// and return a code that specifies the generator algorithm we will use. A code of
2323/// zero indicates that we cannot generate any of the directVars (in this case, nothing
2324/// should be added to generatedVars). Any non-zero codes will be passed to our generateEvent()
2325/// implementation, but otherwise its value is arbitrary. The default implemetation of
2326/// this method returns zero. Subclasses will usually implement this method using the
2327/// matchArgs() methods to advertise the algorithms they provide.
2328
2329Int_t RooAbsPdf::getGenerator(const RooArgSet &/*directVars*/, RooArgSet &/*generatedVars*/, bool /*staticInitOK*/) const
2330{
2331 return 0 ;
2332}
2333
2334
2335
2336////////////////////////////////////////////////////////////////////////////////
2337/// Interface for one-time initialization to setup the generator for the specified code.
2338
2340{
2341}
2342
2343
2344
2345////////////////////////////////////////////////////////////////////////////////
2346/// Interface for generation of an event using the algorithm
2347/// corresponding to the specified code. The meaning of each code is
2348/// defined by the getGenerator() implementation. The default
2349/// implementation does nothing.
2350
2352{
2353}
2354
2355
2356
2357////////////////////////////////////////////////////////////////////////////////
2358/// Check if given observable can be safely generated using the
2359/// pdfs internal generator mechanism (if that existsP). Observables
2360/// on which a PDF depends via more than route are not safe
2361/// for use with internal generators because they introduce
2362/// correlations not known to the internal generator
2363
2365{
2366 // Arg must be direct server of self
2367 if (!findServer(arg.GetName())) return false ;
2368
2369 // There must be no other dependency routes
2370 for (const auto server : _serverList) {
2371 if(server == &arg) continue;
2372 if(server->dependsOn(arg)) {
2373 return false ;
2374 }
2375 }
2376
2377 return true ;
2378}
2379
2380
2381////////////////////////////////////////////////////////////////////////////////
2382/// Generate a new dataset containing the specified variables with events sampled from our distribution.
2383/// \param[in] whatVars Choose variables in which to generate events. Variables not listed here will remain
2384/// constant and not be used for event generation
2385/// \param[in] arg1,arg2,arg3,arg4,arg5,arg6 Optional RooCmdArg to change behaviour of generateBinned()
2386/// \return RooDataHist *, to be managed by caller.
2387///
2388/// Generate the specified number of events or expectedEvents() if not specified.
2389///
2390/// Any variables of this PDF that are not in whatVars will use their
2391/// current values and be treated as fixed parameters. Returns zero
2392/// in case of an error. The caller takes ownership of the returned
2393/// dataset.
2394///
2395/// The following named arguments are supported
2396/// | Type of CmdArg | Effect on generation
2397/// |---------------------------|-----------------------
2398/// | `Name(const char* name)` | Name of the output dataset
2399/// | `Verbose(bool flag)` | Print informational messages during event generation
2400/// | `NumEvents(int nevt)` | Generate specified number of events
2401/// | `Extended()` | The actual number of events generated will be sampled from a Poisson distribution with mu=nevt. This can be *much* faster for peaked PDFs, but the number of events is not exactly what was requested.
2402/// | `ExpectedData()` | Return a binned dataset _without_ statistical fluctuations (also aliased as Asimov())
2403///
2404
2405RooDataHist *RooAbsPdf::generateBinned(const RooArgSet& whatVars, const RooCmdArg& arg1,const RooCmdArg& arg2,
2406 const RooCmdArg& arg3,const RooCmdArg& arg4, const RooCmdArg& arg5,const RooCmdArg& arg6) const
2407{
2408
2409 // Select the pdf-specific commands
2410 RooCmdConfig pc(Form("RooAbsPdf::generate(%s)",GetName())) ;
2411 pc.defineString("dsetName","Name",0,"") ;
2412 pc.defineInt("verbose","Verbose",0,0) ;
2413 pc.defineInt("extended","Extended",0,0) ;
2414 pc.defineInt("nEvents","NumEvents",0,0) ;
2415 pc.defineDouble("nEventsD","NumEventsD",0,-1.) ;
2416 pc.defineInt("expectedData","ExpectedData",0,0) ;
2417
2418 // Process and check varargs
2419 pc.process(arg1,arg2,arg3,arg4,arg5,arg6) ;
2420 if (!pc.ok(true)) {
2421 return 0 ;
2422 }
2423
2424 // Decode command line arguments
2425 double nEvents = pc.getDouble("nEventsD") ;
2426 if (nEvents<0) {
2427 nEvents = pc.getInt("nEvents") ;
2428 }
2429 //bool verbose = pc.getInt("verbose") ;
2430 bool extended = pc.getInt("extended") ;
2431 bool expectedData = pc.getInt("expectedData") ;
2432 const char* dsetName = pc.getString("dsetName") ;
2433
2434 if (extended) {
2435 //nEvents = (nEvents==0?Int_t(expectedEvents(&whatVars)+0.5):nEvents) ;
2436 nEvents = (nEvents==0 ? expectedEvents(&whatVars) :nEvents) ;
2437 cxcoutI(Generation) << " Extended mode active, number of events generated (" << nEvents << ") is Poisson fluctuation on "
2438 << GetName() << "::expectedEvents() = " << nEvents << endl ;
2439 // If Poisson fluctuation results in zero events, stop here
2440 if (nEvents==0) {
2441 return 0 ;
2442 }
2443 } else if (nEvents==0) {
2444 cxcoutI(Generation) << "No number of events specified , number of events generated is "
2445 << GetName() << "::expectedEvents() = " << expectedEvents(&whatVars)<< endl ;
2446 }
2447
2448 // Forward to appropriate implementation
2449 RooDataHist* data = generateBinned(whatVars,nEvents,expectedData,extended) ;
2450
2451 // Rename dataset to given name if supplied
2452 if (dsetName && strlen(dsetName)>0) {
2453 data->SetName(dsetName) ;
2454 }
2455
2456 return data ;
2457}
2458
2459
2460
2461
2462////////////////////////////////////////////////////////////////////////////////
2463/// Generate a new dataset containing the specified variables with
2464/// events sampled from our distribution.
2465///
2466/// \param[in] whatVars Variables that values should be generated for.
2467/// \param[in] nEvents How many events to generate. If `nEvents <=0`, use the value returned by expectedEvents() as target.
2468/// \param[in] expectedData If set to true (false by default), the returned histogram returns the 'expected'
2469/// data sample, i.e. no statistical fluctuations are present.
2470/// \param[in] extended For each bin, generate Poisson(x, mu) events, where `mu` is chosen such that *on average*,
2471/// one would obtain `nEvents` events. This means that the true number of events will fluctuate around the desired value,
2472/// but the generation happens a lot faster.
2473/// Especially if the PDF is sharply peaked, the multinomial event generation necessary to generate *exactly* `nEvents` events can
2474/// be very slow.
2475///
2476/// The binning used for generation of events is the currently set binning for the variables.
2477/// It can e.g. be changed using
2478/// ```
2479/// x.setBins(15);
2480/// x.setRange(-5., 5.);
2481/// pdf.generateBinned(RooArgSet(x), 1000);
2482/// ```
2483///
2484/// Any variables of this PDF that are not in `whatVars` will use their
2485/// current values and be treated as fixed parameters.
2486/// \return RooDataHist* owned by the caller. Returns `nullptr` in case of an error.
2487RooDataHist *RooAbsPdf::generateBinned(const RooArgSet &whatVars, double nEvents, bool expectedData, bool extended) const
2488{
2489 // Create empty RooDataHist
2490 RooDataHist* hist = new RooDataHist("genData","genData",whatVars) ;
2491
2492 // Scale to number of events and introduce Poisson fluctuations
2493 if (nEvents<=0) {
2494 if (!canBeExtended()) {
2495 coutE(InputArguments) << "RooAbsPdf::generateBinned(" << GetName() << ") ERROR: No event count provided and p.d.f does not provide expected number of events" << endl ;
2496 delete hist ;
2497 return nullptr;
2498 } else {
2499
2500 // Don't round in expectedData or extended mode
2501 if (expectedData || extended) {
2502 nEvents = expectedEvents(&whatVars) ;
2503 } else {
2504 nEvents = std::round(expectedEvents(&whatVars));
2505 }
2506 }
2507 }
2508
2509 // Sample p.d.f. distribution
2510 fillDataHist(hist,&whatVars,1,true) ;
2511
2512 vector<int> histOut(hist->numEntries()) ;
2513 double histMax(-1) ;
2514 Int_t histOutSum(0) ;
2515 for (int i=0 ; i<hist->numEntries() ; i++) {
2516 hist->get(i) ;
2517 if (expectedData) {
2518
2519 // Expected data, multiply p.d.f by nEvents
2520 double w=hist->weight()*nEvents ;
2521 hist->set(i, w, sqrt(w));
2522
2523 } else if (extended) {
2524
2525 // Extended mode, set contents to Poisson(pdf*nEvents)
2526 double w = RooRandom::randomGenerator()->Poisson(hist->weight()*nEvents) ;
2527 hist->set(w,sqrt(w)) ;
2528
2529 } else {
2530
2531 // Regular mode, fill array of weights with Poisson(pdf*nEvents), but to not fill
2532 // histogram yet.
2533 if (hist->weight()>histMax) {
2534 histMax = hist->weight() ;
2535 }
2536 histOut[i] = RooRandom::randomGenerator()->Poisson(hist->weight()*nEvents) ;
2537 histOutSum += histOut[i] ;
2538 }
2539 }
2540
2541
2542 if (!expectedData && !extended) {
2543
2544 // Second pass for regular mode - Trim/Extend dataset to exact number of entries
2545
2546 // Calculate difference between what is generated so far and what is requested
2547 Int_t nEvtExtra = std::abs(Int_t(nEvents)-histOutSum) ;
2548 Int_t wgt = (histOutSum>nEvents) ? -1 : 1 ;
2549
2550 // Perform simple binned accept/reject procedure to get to exact event count
2551 std::size_t counter = 0;
2552 bool havePrintedInfo = false;
2553 while(nEvtExtra>0) {
2554
2555 Int_t ibinRand = RooRandom::randomGenerator()->Integer(hist->numEntries()) ;
2556 hist->get(ibinRand) ;
2557 double ranY = RooRandom::randomGenerator()->Uniform(histMax) ;
2558
2559 if (ranY<hist->weight()) {
2560 if (wgt==1) {
2561 histOut[ibinRand]++ ;
2562 } else {
2563 // If weight is negative, prior bin content must be at least 1
2564 if (histOut[ibinRand]>0) {
2565 histOut[ibinRand]-- ;
2566 } else {
2567 continue ;
2568 }
2569 }
2570 nEvtExtra-- ;
2571 }
2572
2573 if ((counter++ > 10*nEvents || nEvents > 1.E7) && !havePrintedInfo) {
2574 havePrintedInfo = true;
2575 coutP(Generation) << "RooAbsPdf::generateBinned(" << GetName() << ") Performing costly accept/reject sampling. If this takes too long, use "
2576 << "extended mode to speed up the process." << std::endl;
2577 }
2578 }
2579
2580 // Transfer working array to histogram
2581 for (int i=0 ; i<hist->numEntries() ; i++) {
2582 hist->get(i) ;
2583 hist->set(histOut[i],sqrt(1.0*histOut[i])) ;
2584 }
2585
2586 } else if (expectedData) {
2587
2588 // Second pass for expectedData mode -- Normalize to exact number of requested events
2589 // Minor difference may be present in first round due to difference between
2590 // bin average and bin integral in sampling bins
2591 double corr = nEvents/hist->sumEntries() ;
2592 for (int i=0 ; i<hist->numEntries() ; i++) {
2593 hist->get(i) ;
2594 hist->set(hist->weight()*corr,sqrt(hist->weight()*corr)) ;
2595 }
2596
2597 }
2598
2599 return hist;
2600}
2601
2602
2603
2604////////////////////////////////////////////////////////////////////////////////
2605/// Special generator interface for generation of 'global observables' -- for RooStats tools
2606
2608{
2609 return generate(whatVars,nEvents) ;
2610}
2611
2612namespace {
2613void removeRangeOverlap(std::vector<std::pair<double, double>>& ranges) {
2614 //Sort from left to right
2615 std::sort(ranges.begin(), ranges.end());
2616
2617 for (auto it = ranges.begin(); it != ranges.end(); ++it) {
2618 double& startL = it->first;
2619 double& endL = it->second;
2620
2621 for (auto innerIt = it+1; innerIt != ranges.end(); ++innerIt) {
2622 const double startR = innerIt->first;
2623 const double endR = innerIt->second;
2624
2625 if (startL <= startR && startR <= endL) {
2626 //Overlapping ranges, extend left one
2627 endL = std::max(endL, endR);
2628 *innerIt = make_pair(0., 0.);
2629 }
2630 }
2631 }
2632
2633 auto newEnd = std::remove_if(ranges.begin(), ranges.end(),
2634 [](const std::pair<double,double>& input){
2635 return input.first == input.second;
2636 });
2637 ranges.erase(newEnd, ranges.end());
2638}
2639}
2640
2641
2642////////////////////////////////////////////////////////////////////////////////
2643/// Plot (project) PDF on specified frame.
2644/// - If a PDF is plotted in an empty frame, it
2645/// will show a unit-normalized curve in the frame variable. When projecting a multi-
2646/// dimensional PDF onto the frame axis, hidden parameters are taken are taken at
2647/// their current value.
2648/// - If a PDF is plotted in a frame in which a dataset has already been plotted, it will
2649/// show a projection integrated over all variables that were present in the shown
2650/// dataset (except for the one on the x-axis). The normalization of the curve will
2651/// be adjusted to the event count of the plotted dataset. An informational message
2652/// will be printed for each projection step that is performed.
2653/// - If a PDF is plotted in a frame showing a dataset *after* a fit, the above happens,
2654/// but the PDF will be drawn and normalised only in the fit range. If this is not desired,
2655/// plotting and normalisation range can be overridden using Range() and NormRange() as
2656/// documented in the table below.
2657///
2658/// This function takes the following named arguments (for more arguments, see also
2659/// RooAbsReal::plotOn(RooPlot*,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,
2660/// const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,const RooCmdArg&,
2661/// const RooCmdArg&) const )
2662///
2663///
2664/// <table>
2665/// <tr><th> Type of argument <th> Controlling normalisation
2666/// <tr><td> `NormRange(const char* name)` <td> Calculate curve normalization w.r.t. specified range[s].
2667/// See the tutorial rf212_plottingInRanges_blinding.C
2668/// \note Setting a Range() by default also sets a NormRange() on the same range, meaning that the
2669/// PDF is plotted and normalised in the same range. Overriding this can be useful if the PDF was fit
2670/// in limited range[s] such as side bands, `NormRange("sidebandLeft,sidebandRight")`, but the PDF
2671/// should be drawn in the full range, `Range("")`.
2672///
2673/// <tr><td> `Normalization(double scale, ScaleType code)` <td> Adjust normalization by given scale factor.
2674/// Interpretation of number depends on code:
2675/// `RooAbsReal::Relative`: relative adjustment factor
2676/// `RooAbsReal::NumEvent`: scale to match given number of events.
2677///
2678/// <tr><th> Type of argument <th> Misc control
2679/// <tr><td> `Name(const chat* name)` <td> Give curve specified name in frame. Useful if curve is to be referenced later
2680/// <tr><td> `Asymmetry(const RooCategory& c)` <td> Show the asymmetry of the PDF in given two-state category
2681/// \f$ \frac{F(+)-F(-)}{F(+)+F(-)} \f$ rather than the PDF projection. Category must have two
2682/// states with indices -1 and +1 or three states with indeces -1,0 and +1.
2683/// <tr><td> `ShiftToZero(bool flag)` <td> Shift entire curve such that lowest visible point is at exactly zero.
2684/// Mostly useful when plotting -log(L) or \f$ \chi^2 \f$ distributions
2685/// <tr><td> `AddTo(const char* name, double_t wgtSelf, double_t wgtOther)` <td> Create a projection of this PDF onto the x-axis, but
2686/// instead of plotting it directly, add it to an existing curve with given name (and relative weight factors).
2687/// <tr><td> `Components(const char* names)` <td> When plotting sums of PDFs, plot only the named components (*e.g.* only
2688/// the signal of a signal+background model).
2689/// <tr><td> `Components(const RooArgSet& compSet)` <td> As above, but pass a RooArgSet of the components themselves.
2690///
2691/// <tr><th> Type of argument <th> Projection control
2692/// <tr><td> `Slice(const RooArgSet& set)` <td> Override default projection behaviour by omitting observables listed
2693/// in set from the projection, i.e. by not integrating over these.
2694/// Slicing is usually only sensible in discrete observables, by e.g. creating a slice
2695/// of the PDF at the current value of the category observable.
2696/// <tr><td> `Slice(RooCategory& cat, const char* label)` <td> Override default projection behaviour by omitting the specified category
2697/// observable from the projection, i.e., by not integrating over all states of this category.
2698/// The slice is positioned at the given label value. Multiple Slice() commands can be given to specify slices
2699/// in multiple observables, e.g.
2700/// ```{.cpp}
2701/// pdf.plotOn(frame, Slice(tagCategory, "2tag"), Slice(jetCategory, "3jet"));
2702/// ```
2703/// <tr><td> `Project(const RooArgSet& set)` <td> Override default projection behaviour by projecting
2704/// over observables given in set, completely ignoring the default projection behavior. Advanced use only.
2705/// <tr><td> `ProjWData(const RooAbsData& d)` <td> Override default projection _technique_ (integration). For observables
2706/// present in given dataset projection of PDF is achieved by constructing an average over all observable
2707/// values in given set. Consult RooFit plotting tutorial for further explanation of meaning & use of this technique
2708/// <tr><td> `ProjWData(const RooArgSet& s, const RooAbsData& d)` <td> As above but only consider subset 's' of
2709/// observables in dataset 'd' for projection through data averaging
2710/// <tr><td> `ProjectionRange(const char* rn)` <td> When projecting the PDF onto the plot axis, it is usually integrated
2711/// over the full range of the invisible variables. The ProjectionRange overrides this.
2712/// This is useful if the PDF was fitted in a limited range in y, but it is now projected onto x. If
2713/// `ProjectionRange("<name of fit range>")` is passed, the projection is normalised correctly.
2714///
2715/// <tr><th> Type of argument <th> Plotting control
2716/// <tr><td> `LineStyle(Int_t style)` <td> Select line style by ROOT line style code, default is solid
2717/// <tr><td> `LineColor(Int_t color)` <td> Select line color by ROOT color code, default is blue
2718/// <tr><td> `LineWidth(Int_t width)` <td> Select line with in pixels, default is 3
2719/// <tr><td> `FillStyle(Int_t style)` <td> Select fill style, default is not filled. If a filled style is selected,
2720/// also use VLines() to add vertical downward lines at end of curve to ensure proper closure
2721/// <tr><td> `FillColor(Int_t color)` <td> Select fill color by ROOT color code
2722/// <tr><td> `Range(const char* name)` <td> Only draw curve in range defined by given name. Multiple comma-separated ranges can be given.
2723/// An empty string "" or `nullptr` means to use the default range of the variable.
2724/// <tr><td> `Range(double lo, double hi)` <td> Only draw curve in specified range
2725/// <tr><td> `VLines()` <td> Add vertical lines to y=0 at end points of curve
2726/// <tr><td> `Precision(double eps)` <td> Control precision of drawn curve w.r.t to scale of plot, default is 1e-3. A higher precision will
2727/// result in more and more densely spaced curve points. A negative precision value will disable
2728/// adaptive point spacing and restrict sampling to the grid point of points defined by the binning
2729/// of the plotted observable (recommended for expensive functions such as profile likelihoods)
2730/// <tr><td> `Invisible(bool flag)` <td> Add curve to frame, but do not display. Useful in combination AddTo()
2731/// <tr><td> `VisualizeError(const RooFitResult& fitres, double Z=1, bool linearMethod=true)`
2732/// <td> Visualize the uncertainty on the parameters, as given in fitres, at 'Z' sigma.
2733/// The linear method is fast but may not be accurate in the presence of strong correlations (~>0.9) and at Z>2 due to linear and Gaussian approximations made.
2734/// Intervals from the sampling method can be asymmetric, and may perform better in the presence of strong correlations, but may take (much) longer to calculate
2735/// \note To include the uncertainty from the expected number of events,
2736/// the Normalization() argument with `ScaleType` `RooAbsReal::RelativeExpected` has to be passed, e.g.
2737/// ```{.cpp}
2738/// pdf.plotOn(frame, VisualizeError(fitResult), Normalization(1.0, RooAbsReal::RelativeExpected));
2739/// ```
2740///
2741/// <tr><td> `VisualizeError(const RooFitResult& fitres, const RooArgSet& param, double Z=1, bool linearMethod=true)`
2742/// <td> Visualize the uncertainty on the subset of parameters 'param', as given in fitres, at 'Z' sigma
2743/// </table>
2744
2746{
2747
2748 // Pre-processing if p.d.f. contains a fit range and there is no command specifying one,
2749 // add a fit range as default range
2750 std::unique_ptr<RooCmdArg> plotRange;
2751 std::unique_ptr<RooCmdArg> normRange2;
2752 if (getStringAttribute("fitrange") && !cmdList.FindObject("Range") &&
2753 !cmdList.FindObject("RangeWithName")) {
2754 plotRange.reset(static_cast<RooCmdArg*>(RooFit::Range(getStringAttribute("fitrange")).Clone()));
2755 cmdList.Add(plotRange.get());
2756 }
2757
2758 if (getStringAttribute("fitrange") && !cmdList.FindObject("NormRange")) {
2759 normRange2.reset(static_cast<RooCmdArg*>(RooFit::NormRange(getStringAttribute("fitrange")).Clone()));
2760 cmdList.Add(normRange2.get());
2761 }
2762
2763 if (plotRange || normRange2) {
2764 coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") p.d.f was fitted in a subrange and no explicit "
2765 << (plotRange?"Range()":"") << ((plotRange&&normRange2)?" and ":"")
2766 << (normRange2?"NormRange()":"") << " was specified. Plotting / normalising in fit range. To override, do one of the following"
2767 << "\n\t- Clear the automatic fit range attribute: <pdf>.removeStringAttribute(\"fitrange\");"
2768 << "\n\t- Explicitly specify the plotting range: Range(\"<rangeName>\")."
2769 << "\n\t- Explicitly specify where to compute the normalisation: NormRange(\"<rangeName>\")."
2770 << "\n\tThe default (full) range can be denoted with Range(\"\") / NormRange(\"\")."<< endl ;
2771 }
2772
2773 // Sanity checks
2774 if (plotSanityChecks(frame)) return frame ;
2775
2776 // Select the pdf-specific commands
2777 RooCmdConfig pc(Form("RooAbsPdf::plotOn(%s)",GetName())) ;
2778 pc.defineDouble("scaleFactor","Normalization",0,1.0) ;
2779 pc.defineInt("scaleType","Normalization",0,Relative) ;
2780 pc.defineSet("compSet","SelectCompSet",0) ;
2781 pc.defineString("compSpec","SelectCompSpec",0) ;
2782 pc.defineObject("asymCat","Asymmetry",0) ;
2783 pc.defineDouble("rangeLo","Range",0,-999.) ;
2784 pc.defineDouble("rangeHi","Range",1,-999.) ;
2785 pc.defineString("rangeName","RangeWithName",0,"") ;
2786 pc.defineString("normRangeName","NormRange",0,"") ;
2787 pc.defineInt("rangeAdjustNorm","Range",0,0) ;
2788 pc.defineInt("rangeWNAdjustNorm","RangeWithName",0,0) ;
2789 pc.defineMutex("SelectCompSet","SelectCompSpec") ;
2790 pc.defineMutex("Range","RangeWithName") ;
2791 pc.allowUndefined() ; // unknowns may be handled by RooAbsReal
2792
2793 // Process and check varargs
2794 pc.process(cmdList) ;
2795 if (!pc.ok(true)) {
2796 return frame ;
2797 }
2798
2799 // Decode command line arguments
2800 ScaleType stype = (ScaleType) pc.getInt("scaleType") ;
2801 double scaleFactor = pc.getDouble("scaleFactor") ;
2802 const RooAbsCategoryLValue* asymCat = (const RooAbsCategoryLValue*) pc.getObject("asymCat") ;
2803 const char* compSpec = pc.getString("compSpec") ;
2804 const RooArgSet* compSet = pc.getSet("compSet");
2805 bool haveCompSel = ( (compSpec && strlen(compSpec)>0) || compSet) ;
2806
2807 // Suffix for curve name
2808 std::string nameSuffix ;
2809 if (compSpec && strlen(compSpec)>0) {
2810 nameSuffix.append("_Comp[") ;
2811 nameSuffix.append(compSpec) ;
2812 nameSuffix.append("]") ;
2813 } else if (compSet) {
2814 nameSuffix.append("_Comp[") ;
2815 nameSuffix.append(compSet->contentsString().c_str()) ;
2816 nameSuffix.append("]") ;
2817 }
2818
2819 // Remove PDF-only commands from command list
2820 RooCmdConfig::stripCmdList(cmdList,"SelectCompSet,SelectCompSpec") ;
2821
2822 // Adjust normalization, if so requested
2823 if (asymCat) {
2824 RooCmdArg cnsuffix("CurveNameSuffix",0,0,0,0,nameSuffix.c_str(),0,0,0) ;
2825 cmdList.Add(&cnsuffix);
2826 return RooAbsReal::plotOn(frame,cmdList) ;
2827 }
2828
2829 // More sanity checks
2830 double nExpected(1) ;
2831 if (stype==RelativeExpected) {
2832 if (!canBeExtended()) {
2833 coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName()
2834 << "): ERROR the 'Expected' scale option can only be used on extendable PDFs" << endl ;
2835 return frame ;
2836 }
2837 frame->updateNormVars(*frame->getPlotVar()) ;
2838 nExpected = expectedEvents(frame->getNormVars()) ;
2839 }
2840
2841 if (stype != Raw) {
2842
2843 if (frame->getFitRangeNEvt() && stype==Relative) {
2844
2845 bool hasCustomRange(false), adjustNorm(false) ;
2846
2847 std::vector<pair<double,double> > rangeLim;
2848
2849 // Retrieve plot range to be able to adjust normalization to data
2850 if (pc.hasProcessed("Range")) {
2851
2852 double rangeLo = pc.getDouble("rangeLo") ;
2853 double rangeHi = pc.getDouble("rangeHi") ;
2854 rangeLim.push_back(make_pair(rangeLo,rangeHi)) ;
2855 adjustNorm = pc.getInt("rangeAdjustNorm") ;
2856 hasCustomRange = true ;
2857
2858 coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") only plotting range ["
2859 << rangeLo << "," << rangeHi << "]" ;
2860 if (!pc.hasProcessed("NormRange")) {
2861 ccoutI(Plotting) << ", curve is normalized to data in " << (adjustNorm?"given":"full") << " range" << endl ;
2862 } else {
2863 ccoutI(Plotting) << endl ;
2864 }
2865
2866 nameSuffix.append(Form("_Range[%f_%f]",rangeLo,rangeHi)) ;
2867
2868 } else if (pc.hasProcessed("RangeWithName")) {
2869
2870 for (const std::string& rangeNameToken : ROOT::Split(pc.getString("rangeName", "", false), ",")) {
2871 const char* thisRangeName = rangeNameToken.empty() ? nullptr : rangeNameToken.c_str();
2872 if (thisRangeName && !frame->getPlotVar()->hasRange(thisRangeName)) {
2873 coutE(Plotting) << "Range '" << rangeNameToken << "' not defined for variable '"
2874 << frame->getPlotVar()->GetName() << "'. Ignoring ..." << std::endl;
2875 continue;
2876 }
2877 rangeLim.push_back(frame->getPlotVar()->getRange(thisRangeName));
2878 }
2879 adjustNorm = pc.getInt("rangeWNAdjustNorm") ;
2880 hasCustomRange = true ;
2881
2882 coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") only plotting range '" << pc.getString("rangeName", "", false) << "'" ;
2883 if (!pc.hasProcessed("NormRange")) {
2884 ccoutI(Plotting) << ", curve is normalized to data in " << (adjustNorm?"given":"full") << " range" << endl ;
2885 } else {
2886 ccoutI(Plotting) << endl ;
2887 }
2888
2889 nameSuffix.append(Form("_Range[%s]",pc.getString("rangeName"))) ;
2890 }
2891 // Specification of a normalization range override those in a regular range
2892 if (pc.hasProcessed("NormRange")) {
2893 rangeLim.clear();
2894 for (const auto& rangeNameToken : ROOT::Split(pc.getString("normRangeName", "", false), ",")) {
2895 const char* thisRangeName = rangeNameToken.empty() ? nullptr : rangeNameToken.c_str();
2896 if (thisRangeName && !frame->getPlotVar()->hasRange(thisRangeName)) {
2897 coutE(Plotting) << "Range '" << rangeNameToken << "' not defined for variable '"
2898 << frame->getPlotVar()->GetName() << "'. Ignoring ..." << std::endl;
2899 continue;
2900 }
2901 rangeLim.push_back(frame->getPlotVar()->getRange(thisRangeName));
2902 }
2903 adjustNorm = true ;
2904 hasCustomRange = true ;
2905 coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") p.d.f. curve is normalized using explicit choice of ranges '" << pc.getString("normRangeName", "", false) << "'" << endl ;
2906
2907 nameSuffix.append(Form("_NormRange[%s]",pc.getString("rangeName"))) ;
2908
2909 }
2910
2911 if (hasCustomRange && adjustNorm) {
2912 // If overlapping ranges were given, remove them now
2913 const std::size_t oldSize = rangeLim.size();
2914 removeRangeOverlap(rangeLim);
2915
2916 if (oldSize != rangeLim.size() && !pc.hasProcessed("NormRange")) {
2917 // User gave overlapping ranges. This leads to double-counting events and integrals, and must
2918 // therefore be avoided. If a NormRange has been given, the overlap is alreay gone.
2919 // It's safe to plot even with overlap now.
2920 coutE(Plotting) << "Requested plot/integration ranges overlap. For correct plotting, new ranges "
2921 "will be defined." << std::endl;
2922 auto plotVar = dynamic_cast<RooRealVar*>(frame->getPlotVar());
2923 assert(plotVar);
2924 std::string rangesNoOverlap;
2925 for (auto it = rangeLim.begin(); it != rangeLim.end(); ++it) {
2926 std::stringstream rangeName;
2927 rangeName << "Remove_overlap_range_" << it - rangeLim.begin();
2928 plotVar->setRange(rangeName.str().c_str(), it->first, it->second);
2929 if (!rangesNoOverlap.empty())
2930 rangesNoOverlap += ",";
2931 rangesNoOverlap += rangeName.str();
2932 }
2933
2934 auto rangeArg = static_cast<RooCmdArg*>(cmdList.FindObject("RangeWithName"));
2935 if (rangeArg)
2936 rangeArg->setString(0, rangesNoOverlap.c_str());
2937 else {
2938 plotRange = std::make_unique<RooCmdArg>(RooFit::Range(rangesNoOverlap.c_str()));
2939 cmdList.Add(plotRange.get());
2940 }
2941 }
2942
2943 double rangeNevt(0) ;
2944 for (const auto& riter : rangeLim) {
2945 double nevt= frame->getFitRangeNEvt(riter.first, riter.second);
2946 rangeNevt += nevt ;
2947 }
2948
2949 scaleFactor *= rangeNevt/nExpected ;
2950
2951 } else {
2952 scaleFactor *= frame->getFitRangeNEvt()/nExpected ;
2953 }
2954 } else if (stype==RelativeExpected) {
2955 scaleFactor *= nExpected ;
2956 } else if (stype==NumEvent) {
2957 scaleFactor /= nExpected ;
2958 }
2959 scaleFactor *= frame->getFitRangeBinW() ;
2960 }
2961 frame->updateNormVars(*frame->getPlotVar()) ;
2962
2963 // Append overriding scale factor command at end of original command list
2964 RooCmdArg tmp = RooFit::Normalization(scaleFactor,Raw) ;
2965 tmp.setInt(1,1) ; // Flag this normalization command as created for internal use (so that VisualizeError can strip it)
2966 cmdList.Add(&tmp) ;
2967
2968 // Was a component selected requested
2969 if (haveCompSel) {
2970
2971 // Get complete set of tree branch nodes
2972 RooArgSet branchNodeSet ;
2973 branchNodeServerList(&branchNodeSet) ;
2974
2975 // Discard any non-RooAbsReal nodes
2976 for (const auto arg : branchNodeSet) {
2977 if (!dynamic_cast<RooAbsReal*>(arg)) {
2978 branchNodeSet.remove(*arg) ;
2979 }
2980 }
2981
2982 // Obtain direct selection
2983 std::unique_ptr<RooArgSet> dirSelNodes;
2984 if (compSet) {
2985 dirSelNodes.reset(static_cast<RooArgSet*>(branchNodeSet.selectCommon(*compSet)));
2986 } else {
2987 dirSelNodes.reset(static_cast<RooArgSet*>(branchNodeSet.selectByName(compSpec)));
2988 }
2989 if (dirSelNodes->getSize()>0) {
2990 coutI(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") directly selected PDF components: " << *dirSelNodes << endl ;
2991
2992 // Do indirect selection and activate both
2993 plotOnCompSelect(dirSelNodes.get());
2994 } else {
2995 if (compSet) {
2996 coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") ERROR: component selection set " << *compSet << " does not match any components of p.d.f." << endl ;
2997 } else {
2998 coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName() << ") ERROR: component selection expression '" << compSpec << "' does not select any components of p.d.f." << endl ;
2999 }
3000 return 0 ;
3001 }
3002 }
3003
3004
3005 RooCmdArg cnsuffix("CurveNameSuffix",0,0,0,0,nameSuffix.c_str(),0,0,0) ;
3006 cmdList.Add(&cnsuffix);
3007
3008 RooPlot* ret = RooAbsReal::plotOn(frame,cmdList) ;
3009
3010 // Restore selection status ;
3011 if (haveCompSel) plotOnCompSelect(0) ;
3012
3013 return ret ;
3014}
3015
3016
3017//_____________________________________________________________________________
3018/// Plot oneself on 'frame'. In addition to features detailed in RooAbsReal::plotOn(),
3019/// the scale factor for a PDF can be interpreted in three different ways. The interpretation
3020/// is controlled by ScaleType
3021/// ```
3022/// Relative - Scale factor is applied on top of PDF normalization scale factor
3023/// NumEvent - Scale factor is interpreted as a number of events. The surface area
3024/// under the PDF curve will match that of a histogram containing the specified
3025/// number of event
3026/// Raw - Scale factor is applied to the raw (projected) probability density.
3027/// Not too useful, option provided for completeness.
3028/// ```
3029// coverity[PASS_BY_VALUE]
3031{
3032
3033 // Sanity checks
3034 if (plotSanityChecks(frame)) return frame ;
3035
3036 // More sanity checks
3037 double nExpected(1) ;
3038 if (o.stype==RelativeExpected) {
3039 if (!canBeExtended()) {
3040 coutE(Plotting) << "RooAbsPdf::plotOn(" << GetName()
3041 << "): ERROR the 'Expected' scale option can only be used on extendable PDFs" << endl ;
3042 return frame ;
3043 }
3044 frame->updateNormVars(*frame->getPlotVar()) ;
3045 nExpected = expectedEvents(frame->getNormVars()) ;
3046 }
3047
3048 // Adjust normalization, if so requested
3049 if (o.stype != Raw) {
3050
3051 if (frame->getFitRangeNEvt() && o.stype==Relative) {
3052 // If non-default plotting range is specified, adjust number of events in fit range
3053 o.scaleFactor *= frame->getFitRangeNEvt()/nExpected ;
3054 } else if (o.stype==RelativeExpected) {
3055 o.scaleFactor *= nExpected ;
3056 } else if (o.stype==NumEvent) {
3057 o.scaleFactor /= nExpected ;
3058 }
3059 o.scaleFactor *= frame->getFitRangeBinW() ;
3060 }
3061 frame->updateNormVars(*frame->getPlotVar()) ;
3062
3063 return RooAbsReal::plotOn(frame,o) ;
3064}
3065
3066
3067
3068
3069////////////////////////////////////////////////////////////////////////////////
3070/// The following named arguments are supported
3071/// <table>
3072/// <tr><th> Type of CmdArg <th> Effect on parameter box
3073/// <tr><td> `Parameters(const RooArgSet& param)` <td> Only the specified subset of parameters will be shown. By default all non-constant parameters are shown.
3074/// <tr><td> `ShowConstants(bool flag)` <td> Also display constant parameters
3075/// <tr><td> `Format(const char* what,...)` <td> Parameter formatting options.
3076/// | Parameter | Format
3077/// | ---------------------- | --------------------------
3078/// | `const char* what` | Controls what is shown. "N" adds name (alternatively, "T" adds the title), "E" adds error, "A" shows asymmetric error, "U" shows unit, "H" hides the value
3079/// | `FixedPrecision(int n)`| Controls precision, set fixed number of digits
3080/// | `AutoPrecision(int n)` | Controls precision. Number of shown digits is calculated from error + n specified additional digits (1 is sensible default)
3081/// <tr><td> `Label(const chat* label)` <td> Add label to parameter box. Use `\n` for multi-line labels.
3082/// <tr><td> `Layout(double xmin, double xmax, double ymax)` <td> Specify relative position of left/right side of box and top of box.
3083/// Coordinates are given as position on the pad between 0 and 1.
3084/// The lower end of the box is calculated automatically from the number of lines in the box.
3085/// </table>
3086///
3087///
3088/// Example use:
3089/// ```
3090/// pdf.paramOn(frame, Label("fit result"), Format("NEU",AutoPrecision(1)) ) ;
3091/// ```
3092///
3093
3094RooPlot* RooAbsPdf::paramOn(RooPlot* frame, const RooCmdArg& arg1, const RooCmdArg& arg2,
3095 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
3096 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
3097{
3098 // Stuff all arguments in a list
3099 RooLinkedList cmdList;
3100 cmdList.Add(const_cast<RooCmdArg*>(&arg1)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg2)) ;
3101 cmdList.Add(const_cast<RooCmdArg*>(&arg3)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg4)) ;
3102 cmdList.Add(const_cast<RooCmdArg*>(&arg5)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg6)) ;
3103 cmdList.Add(const_cast<RooCmdArg*>(&arg7)) ; cmdList.Add(const_cast<RooCmdArg*>(&arg8)) ;
3104
3105 // Select the pdf-specific commands
3106 RooCmdConfig pc(Form("RooAbsPdf::paramOn(%s)",GetName())) ;
3107 pc.defineString("label","Label",0,"") ;
3108 pc.defineDouble("xmin","Layout",0,0.65) ;
3109 pc.defineDouble("xmax","Layout",1,0.9) ;
3110 pc.defineInt("ymaxi","Layout",0,Int_t(0.9*10000)) ;
3111 pc.defineInt("showc","ShowConstants",0,0) ;
3112 pc.defineSet("params","Parameters",0,0) ;
3113 pc.defineInt("dummy","FormatArgs",0,0) ;
3114
3115 // Process and check varargs
3116 pc.process(cmdList) ;
3117 if (!pc.ok(true)) {
3118 return frame ;
3119 }
3120
3121 auto formatCmd = static_cast<RooCmdArg const*>(cmdList.FindObject("FormatArgs")) ;
3122
3123 const char* label = pc.getString("label") ;
3124 double xmin = pc.getDouble("xmin") ;
3125 double xmax = pc.getDouble("xmax") ;
3126 double ymax = pc.getInt("ymaxi") / 10000. ;
3127 int showc = pc.getInt("showc") ;
3128
3129 // Decode command line arguments
3130 std::unique_ptr<RooArgSet> params{getParameters(frame->getNormVars())} ;
3131 if(RooArgSet* requestedParams = pc.getSet("params")) {
3132 params = std::unique_ptr<RooArgSet>{static_cast<RooArgSet*>(params->selectCommon(*requestedParams))};
3133 }
3134 paramOn(frame,*params,showc,label,xmin,xmax,ymax,formatCmd);
3135
3136 return frame ;
3137}
3138
3139
3140////////////////////////////////////////////////////////////////////////////////
3141/// Add a text box with the current parameter values and their errors to the frame.
3142/// Observables of this PDF appearing in the 'data' dataset will be omitted.
3143///
3144/// An optional label will be inserted if passed. Multi-line labels can be generated
3145/// by adding `\n` to the label string. Use 'sigDigits'
3146/// to modify the default number of significant digits printed. The 'xmin,xmax,ymax'
3147/// values specify the initial relative position of the text box in the plot frame.
3148
3149RooPlot* RooAbsPdf::paramOn(RooPlot* frame, const RooArgSet& params, bool showConstants, const char *label,
3150 double xmin, double xmax ,double ymax, const RooCmdArg* formatCmd)
3151{
3152
3153 // parse the options
3154 bool showLabel= (label != 0 && strlen(label) > 0);
3155
3156 // calculate the box's size, adjusting for constant parameters
3157
3158 double ymin(ymax), dy(0.06);
3159 for (const auto param : params) {
3160 auto var = static_cast<RooRealVar*>(param);
3161 if(showConstants || !var->isConstant()) ymin-= dy;
3162 }
3163
3164 std::string labelString = label;
3165 unsigned int numLines = std::count(labelString.begin(), labelString.end(), '\n') + 1;
3166 if (showLabel) ymin -= numLines * dy;
3167
3168 // create the box and set its options
3169 TPaveText *box= new TPaveText(xmin,ymax,xmax,ymin,"BRNDC");
3170 if(!box) return 0;
3171 box->SetName(Form("%s_paramBox",GetName())) ;
3172 box->SetFillColor(0);
3173 box->SetBorderSize(0);
3174 box->SetTextAlign(12);
3175 box->SetTextSize(0.04F);
3176 box->SetFillStyle(0);
3177
3178 for (const auto param : params) {
3179 auto var = static_cast<const RooRealVar*>(param);
3180 if(var->isConstant() && !showConstants) continue;
3181
3182 std::unique_ptr<TString> formatted{formatCmd ? var->format(*formatCmd) : var->format(2, "NELU")};
3183 box->AddText(formatted->Data());
3184 }
3185
3186 // add the optional label if specified
3187 if (showLabel) {
3188 for (const auto& line : ROOT::Split(label, "\n")) {
3189 box->AddText(line.c_str());
3190 }
3191 }
3192
3193 // Add box to frame
3194 frame->addObject(box) ;
3195
3196 return frame ;
3197}
3198
3199
3200
3201
3202////////////////////////////////////////////////////////////////////////////////
3203/// Return expected number of events from this p.d.f for use in extended
3204/// likelihood calculations. This default implementation returns zero
3205
3207{
3208 return 0 ;
3209}
3210
3211
3212
3213////////////////////////////////////////////////////////////////////////////////
3214/// Change global level of verbosity for p.d.f. evaluations
3215
3217{
3218 _verboseEval = stat ;
3219}
3220
3221
3222
3223////////////////////////////////////////////////////////////////////////////////
3224/// Return global level of verbosity for p.d.f. evaluations
3225
3227{
3228 return _verboseEval ;
3229}
3230
3231
3232
3233////////////////////////////////////////////////////////////////////////////////
3234/// Destructor of normalization cache element. If this element
3235/// provides the 'current' normalization stored in RooAbsPdf::_norm
3236/// zero _norm pointer here before object pointed to is deleted here
3237
3239{
3240 // Zero _norm pointer in RooAbsPdf if it is points to our cache payload
3241 if (_owner) {
3242 RooAbsPdf* pdfOwner = static_cast<RooAbsPdf*>(_owner) ;
3243 if (pdfOwner->_norm == _norm) {
3244 pdfOwner->_norm = 0 ;
3245 }
3246 }
3247
3248 delete _norm ;
3249}
3250
3251
3252
3253////////////////////////////////////////////////////////////////////////////////
3254/// Return a p.d.f that represent a projection of this p.d.f integrated over given observables
3255
3257{
3258 // Construct name for new object
3259 std::string name(GetName()) ;
3260 name.append("_Proj[") ;
3261 if (iset.getSize()>0) {
3262 bool first = true;
3263 for(auto const& arg : iset) {
3264 if (first) {
3265 first = false ;
3266 } else {
3267 name.append(",") ;
3268 }
3269 name.append(arg->GetName()) ;
3270 }
3271 }
3272 name.append("]") ;
3273
3274 // Return projected p.d.f.
3275 return new RooProjectedPdf(name.c_str(),name.c_str(),*this,iset) ;
3276}
3277
3278
3279
3280////////////////////////////////////////////////////////////////////////////////
3281/// Create a cumulative distribution function of this p.d.f in terms
3282/// of the observables listed in iset. If no nset argument is given
3283/// the c.d.f normalization is constructed over the integrated
3284/// observables, so that its maximum value is precisely 1. It is also
3285/// possible to choose a different normalization for
3286/// multi-dimensional p.d.f.s: eg. for a pdf f(x,y,z) one can
3287/// construct a partial cdf c(x,y) that only when integrated itself
3288/// over z results in a maximum value of 1. To construct such a cdf pass
3289/// z as argument to the optional nset argument
3290
3292{
3293 return createCdf(iset,RooFit::SupNormSet(nset)) ;
3294}
3295
3296
3297
3298////////////////////////////////////////////////////////////////////////////////
3299/// Create an object that represents the integral of the function over one or more observables listed in `iset`.
3300/// The actual integration calculation is only performed when the return object is evaluated. The name
3301/// of the integral object is automatically constructed from the name of the input function, the variables
3302/// it integrates and the range integrates over
3303///
3304/// The following named arguments are accepted
3305/// | Type of CmdArg | Effect on CDF
3306/// | ---------------------|-------------------
3307/// | SupNormSet(const RooArgSet&) | Observables over which should be normalized _in addition_ to the integration observables
3308/// | ScanNumCdf() | Apply scanning technique if cdf integral involves numeric integration [ default ]
3309/// | ScanAllCdf() | Always apply scanning technique
3310/// | ScanNoCdf() | Never apply scanning technique
3311/// | ScanParameters(Int_t nbins, Int_t intOrder) | Parameters for scanning technique of making CDF: number of sampled bins and order of interpolation applied on numeric cdf
3312
3314 const RooCmdArg& arg3, const RooCmdArg& arg4, const RooCmdArg& arg5,
3315 const RooCmdArg& arg6, const RooCmdArg& arg7, const RooCmdArg& arg8)
3316{
3317 // Define configuration for this method
3318 RooCmdConfig pc(Form("RooAbsReal::createCdf(%s)",GetName())) ;
3319 pc.defineSet("supNormSet","SupNormSet",0,0) ;
3320 pc.defineInt("numScanBins","ScanParameters",0,1000) ;
3321 pc.defineInt("intOrder","ScanParameters",1,2) ;
3322 pc.defineInt("doScanNum","ScanNumCdf",0,1) ;
3323 pc.defineInt("doScanAll","ScanAllCdf",0,0) ;
3324 pc.defineInt("doScanNon","ScanNoCdf",0,0) ;
3325 pc.defineMutex("ScanNumCdf","ScanAllCdf","ScanNoCdf") ;
3326
3327 // Process & check varargs
3328 pc.process(arg1,arg2,arg3,arg4,arg5,arg6,arg7,arg8) ;
3329 if (!pc.ok(true)) {
3330 return 0 ;
3331 }
3332
3333 // Extract values from named arguments
3334 const RooArgSet* snset = pc.getSet("supNormSet",0);
3335 RooArgSet nset ;
3336 if (snset) {
3337 nset.add(*snset) ;
3338 }
3339 Int_t numScanBins = pc.getInt("numScanBins") ;
3340 Int_t intOrder = pc.getInt("intOrder") ;
3341 Int_t doScanNum = pc.getInt("doScanNum") ;
3342 Int_t doScanAll = pc.getInt("doScanAll") ;
3343 Int_t doScanNon = pc.getInt("doScanNon") ;
3344
3345 // If scanning technique is not requested make integral-based cdf and return
3346 if (doScanNon) {
3347 return createIntRI(iset,nset) ;
3348 }
3349 if (doScanAll) {
3350 return createScanCdf(iset,nset,numScanBins,intOrder) ;
3351 }
3352 if (doScanNum) {
3353 std::unique_ptr<RooAbsReal> tmp{createIntegral(iset)} ;
3354 Int_t isNum= !static_cast<RooRealIntegral&>(*tmp).numIntRealVars().empty();
3355
3356 if (isNum) {
3357 coutI(NumIntegration) << "RooAbsPdf::createCdf(" << GetName() << ") integration over observable(s) " << iset << " involves numeric integration," << endl
3358 << " constructing cdf though numeric integration of sampled pdf in " << numScanBins << " bins and applying order "
3359 << intOrder << " interpolation on integrated histogram." << endl
3360 << " To override this choice of technique use argument ScanNone(), to change scan parameters use ScanParameters(nbins,order) argument" << endl ;
3361 }
3362
3363 return isNum ? createScanCdf(iset,nset,numScanBins,intOrder) : createIntRI(iset,nset) ;
3364 }
3365 return 0 ;
3366}
3367
3369{
3370 string name = string(GetName()) + "_NUMCDF_" + integralNameSuffix(iset,&nset).Data() ;
3371 RooRealVar* ivar = (RooRealVar*) iset.first() ;
3372 ivar->setBins(numScanBins,"numcdf") ;
3373 auto ret = std::make_unique<RooNumCdf>(name.c_str(),name.c_str(),*this,*ivar,"numcdf");
3374 ret->setInterpolationOrder(intOrder) ;
3375 return RooFit::Detail::owningPtr(std::move(ret));
3376}
3377
3378
3379
3380
3381////////////////////////////////////////////////////////////////////////////////
3382/// This helper function finds and collects all constraints terms of all component p.d.f.s
3383/// and returns a RooArgSet with all those terms.
3384
3385RooArgSet* RooAbsPdf::getAllConstraints(const RooArgSet& observables, RooArgSet& constrainedParams,
3386 bool stripDisconnected, bool removeConstraintsFromPdf) const
3387{
3388 RooArgSet* ret = new RooArgSet("AllConstraints") ;
3389
3390 std::unique_ptr<RooArgSet> comps(getComponents());
3391 for (const auto arg : *comps) {
3392 auto pdf = dynamic_cast<const RooAbsPdf*>(arg) ;
3393 if (pdf && !ret->find(pdf->GetName())) {
3394 std::unique_ptr<RooArgSet> compRet(
3395 pdf->getConstraints(observables,constrainedParams,stripDisconnected,removeConstraintsFromPdf));
3396 if (compRet) {
3397 ret->add(*compRet,false) ;
3398 }
3399 }
3400 }
3401
3402 return ret ;
3403}
3404
3405
3406////////////////////////////////////////////////////////////////////////////////
3407/// Returns the default numeric MC generator configuration for all RooAbsReals
3408
3410{
3412}
3413
3414
3415////////////////////////////////////////////////////////////////////////////////
3416/// Returns the specialized integrator configuration for _this_ RooAbsReal.
3417/// If this object has no specialized configuration, a null pointer is returned
3418
3420{
3421 return _specGeneratorConfig.get();
3422}
3423
3424
3425
3426////////////////////////////////////////////////////////////////////////////////
3427/// Returns the specialized integrator configuration for _this_ RooAbsReal.
3428/// If this object has no specialized configuration, a null pointer is returned,
3429/// unless createOnTheFly is true in which case a clone of the default integrator
3430/// configuration is created, installed as specialized configuration, and returned
3431
3433{
3434 if (!_specGeneratorConfig && createOnTheFly) {
3435 _specGeneratorConfig = std::make_unique<RooNumGenConfig>(*defaultGeneratorConfig()) ;
3436 }
3437 return _specGeneratorConfig.get();
3438}
3439
3440
3441
3442////////////////////////////////////////////////////////////////////////////////
3443/// Return the numeric MC generator configuration used for this object. If
3444/// a specialized configuration was associated with this object, that configuration
3445/// is returned, otherwise the default configuration for all RooAbsReals is returned
3446
3448{
3449 const RooNumGenConfig* config = specialGeneratorConfig() ;
3450 if (config) return config ;
3451 return defaultGeneratorConfig() ;
3452}
3453
3454
3455
3456////////////////////////////////////////////////////////////////////////////////
3457/// Set the given configuration as default numeric MC generator
3458/// configuration for this object
3459
3461{
3462 _specGeneratorConfig = std::make_unique<RooNumGenConfig>(config);
3463}
3464
3465
3466
3467////////////////////////////////////////////////////////////////////////////////
3468/// Remove the specialized numeric MC generator configuration associated
3469/// with this object
3470
3472{
3473 _specGeneratorConfig.reset();
3474}
3475
3476
3477
3478////////////////////////////////////////////////////////////////////////////////
3479
3481{
3482 delete _genContext ;
3483}
3484
3485
3486////////////////////////////////////////////////////////////////////////////////
3487
3488RooAbsPdf::GenSpec::GenSpec(RooAbsGenContext* context, const RooArgSet& whatVars, RooDataSet* protoData, Int_t nGen,
3489 bool extended, bool randProto, bool resampleProto, TString dsetName, bool init) :
3490 _genContext(context), _whatVars(whatVars), _protoData(protoData), _nGen(nGen), _extended(extended),
3491 _randProto(randProto), _resampleProto(resampleProto), _dsetName(dsetName), _init(init)
3492{
3493}
3494
3495
3496namespace {
3497
3498void sterilizeClientCaches(RooAbsArg & arg) {
3499 auto const& clients = arg.clients();
3500 for(std::size_t iClient = 0; iClient < clients.size(); ++iClient) {
3501
3502 const std::size_t oldClientsSize = clients.size();
3503 RooAbsArg* client = clients[iClient];
3504
3505 for(int iCache = 0; iCache < client->numCaches(); ++iCache) {
3506 if(auto cacheMgr = dynamic_cast<RooObjCacheManager*>(client->getCache(iCache))) {
3507 cacheMgr->sterilize();
3508 }
3509 }
3510
3511 // It can happen that the objects cached by the client are also clients of
3512 // the arg itself! In that case, the position of the client in the client
3513 // list might have changed, and we need to find the new index.
3514 if(clients.size() != oldClientsSize) {
3515 auto clientIter = std::find(clients.begin(), clients.end(), client);
3516 if(clientIter == clients.end()) {
3517 throw std::runtime_error("After a clients caches were cleared, the client was gone! This should not happen.");
3518 }
3519 iClient = std::distance(clients.begin(), clientIter);
3520 }
3521 }
3522}
3523
3524} // namespace
3525
3526
3527////////////////////////////////////////////////////////////////////////////////
3528
3529void RooAbsPdf::setNormRange(const char* rangeName)
3530{
3531 if (rangeName) {
3532 _normRange = rangeName ;
3533 } else {
3534 _normRange.Clear() ;
3535 }
3536
3537 // the stuff that the clients have cached may depend on the normalization range
3538 sterilizeClientCaches(*this);
3539
3540 if (_norm) {
3542 _norm = 0 ;
3543 }
3544}
3545
3546
3547////////////////////////////////////////////////////////////////////////////////
3548
3549void RooAbsPdf::setNormRangeOverride(const char* rangeName)
3550{
3551 if (rangeName) {
3552 _normRangeOverride = rangeName ;
3553 } else {
3555 }
3556
3557 // the stuff that the clients have cached may depend on the normalization range
3558 sterilizeClientCaches(*this);
3559
3560 if (_norm) {
3562 _norm = 0 ;
3563 }
3564}
3565
3566
3567////////////////////////////////////////////////////////////////////////////////
3568/// Hook function intercepting redirectServer calls. Discard current
3569/// normalization object if any server is redirected
3570
3571bool RooAbsPdf::redirectServersHook(const RooAbsCollection & newServerList, bool mustReplaceAll,
3572 bool nameChange, bool isRecursiveStep)
3573{
3574 // If servers are redirected, the cached normalization integrals and
3575 // normalization sets are most likely invalid.
3577
3578 // Object is own by _normCacheManager that will delete object as soon as cache
3579 // is sterilized by server redirect
3580 _norm = nullptr ;
3581
3582 // Similar to the situation with the normalization integral above: if a
3583 // server is redirected, the cached normalization set might not point to
3584 // the right observables anymore. We need to reset it.
3585 setActiveNormSet(nullptr);
3586 return RooAbsReal::redirectServersHook(newServerList, mustReplaceAll, nameChange, isRecursiveStep);
3587}
3588
3589
3590std::unique_ptr<RooAbsArg>
3592{
3593 if (normSet.empty() || selfNormalized()) {
3594 return RooAbsReal::compileForNormSet(normSet, ctx);
3595 }
3596 std::unique_ptr<RooAbsPdf> pdfClone(static_cast<RooAbsPdf *>(this->Clone()));
3597 ctx.compileServers(*pdfClone, normSet);
3598
3599 auto newArg = std::make_unique<RooNormalizedPdf>(*pdfClone, normSet);
3600
3601 // The direct servers are this pdf and the normalization integral, which
3602 // don't need to be compiled further.
3603 for (RooAbsArg *server : newArg->servers()) {
3604 server->setAttribute("_COMPILED");
3605 }
3606 newArg->setAttribute("_COMPILED");
3607 newArg->addOwnedComponents(std::move(pdfClone));
3608 return newArg;
3609}
header file containing the templated implementation of matrix inversion routines for use with ROOT's ...
std::unique_ptr< RooAbsReal > createConstraintTerm(std::string const &name, RooAbsPdf const &pdf, RooAbsData const &data, RooArgSet const *constrainedParameters, RooArgSet const *externalConstraints, RooArgSet const *globalObservables, const char *globalObservablesTag, bool takeGlobalObservablesFromData, bool removeConstraintsFromPdf)
Create the parameter constraint sum to add to the negative log-likelihood.
#define e(i)
Definition RSha256.hxx:103
#define coutI(a)
#define cxcoutI(a)
#define cxcoutD(a)
#define coutP(a)
#define oocoutW(o, a)
#define coutW(a)
#define oocoutI(o, a)
#define coutE(a)
#define ccoutI(a)
#define ccoutD(a)
int Int_t
Definition RtypesCore.h:45
#define ClassImp(name)
Definition Rtypes.h:377
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 input
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
float ymin
float xmax
float ymax
char * Form(const char *fmt,...)
Formats a string in a circular formatting buffer.
Definition TString.cxx:2468
class to compute the Cholesky decomposition of a matrix
bool Invert(M &m) const
place the inverse into m
RooAbsArg is the common abstract base class for objects that represent a value and a "shape" in RooFi...
Definition RooAbsArg.h:78
void clearValueAndShapeDirty() const
Definition RooAbsArg.h:599
void Print(Option_t *options=nullptr) const override
Print the object to the defaultPrintStream().
Definition RooAbsArg.h:321
bool dependsOn(const RooAbsCollection &serverList, const RooAbsArg *ignoreArg=nullptr, bool valueOnly=false) const
Test whether we depend on (ie, are served by) any object in the specified collection.
friend class RooDataSet
Definition RooAbsArg.h:669
void setOperMode(OperMode mode, bool recurseADirty=true)
Set the operation mode of this node.
void setStringAttribute(const Text_t *key, const Text_t *value)
Associate string 'value' to this object under key 'key'.
RooFit::OwningPtr< RooArgSet > getParameters(const RooAbsData *data, bool stripDisconnected=true) const
Create a list of leaf nodes in the arg tree starting with ourself as top node that don't match any of...
RooFit::OwningPtr< RooArgSet > getObservables(const RooArgSet &set, bool valueOnly=true) const
Given a set of possible observables, return the observables that this PDF depends on.
bool addOwnedComponents(const RooAbsCollection &comps)
Take ownership of the contents of 'comps'.
const Text_t * getStringAttribute(const Text_t *key) const
Get string attribute mapped under key 'key'.
virtual std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const
RooFit::OwningPtr< RooArgSet > getComponents() const
Create a RooArgSet with all components (branch nodes) of the expression tree headed by this object.
void removeStringAttribute(const Text_t *key)
Delete a string attribute with a given key.
bool getAttribute(const Text_t *name) const
Check if a named attribute is set. By default, all attributes are unset.
RooFit::OwningPtr< RooArgSet > getVariables(bool stripDisconnected=true) const
Return RooArgSet with all variables (tree leaf nodes of expresssion tree)
RooAbsCache * getCache(Int_t index) const
Return registered cache object by index.
const RefCountList_t & clients() const
List of all clients of this object.
Definition RooAbsArg.h:194
bool isValueDirty() const
Definition RooAbsArg.h:421
virtual void applyWeightSquared(bool flag)
Disables or enables the usage of squared weights.
void setAttribute(const Text_t *name, bool value=true)
Set (default) or clear a named boolean attribute of this object.
void setProxyNormSet(const RooArgSet *nset)
Forward a change in the cached normalization argset to all the registered proxies.
void branchNodeServerList(RooAbsCollection *list, const RooAbsArg *arg=nullptr, bool recurseNonDerived=false) const
Fill supplied list with all branch nodes of the arg tree starting with ourself as top node.
TObject * Clone(const char *newname=nullptr) const override
Make a clone of an object using the Streamer facility.
Definition RooAbsArg.h:90
RefCountList_t _serverList
Definition RooAbsArg.h:634
Int_t numCaches() const
Return number of registered caches.
RooAbsArg * findServer(const char *name) const
Return server of this with name name. Returns nullptr if not found.
Definition RooAbsArg.h:212
OperMode operMode() const
Query the operation mode of this node.
Definition RooAbsArg.h:484
RooAbsArg * _owner
! Pointer to owning RooAbsArg
void setInterpolationOrder(Int_t order)
Set interpolation order of RooHistFunct representing cache histogram.
RooAbsCategoryLValue is the common abstract base class for objects that represent a discrete value th...
RooAbsCollection is an abstract container object that can hold multiple RooAbsArg objects.
virtual bool remove(const RooAbsArg &var, bool silent=false, bool matchByNameOnly=false)
Remove the specified argument from our list.
Int_t getSize() const
Return the number of elements in the collection.
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
void assign(const RooAbsCollection &other) const
Sets the value, cache and constant attribute of any argument in our set that also appears in the othe...
RooAbsArg * first() const
RooAbsCollection * selectByName(const char *nameList, bool verbose=false) const
Create a subset of the current collection, consisting only of those elements with names matching the ...
bool selectCommon(const RooAbsCollection &refColl, RooAbsCollection &outColl) const
Create a subset of the current collection, consisting only of those elements that are contained as we...
std::string contentsString() const
Return comma separated list of contained object names as STL string.
RooAbsArg * find(const char *name) const
Find object with given name in list.
void Print(Option_t *options=nullptr) const override
This method must be overridden when a class wants to print itself.
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition RooAbsData.h:58
virtual const RooArgSet * get() const
Definition RooAbsData.h:102
virtual Int_t numEntries() const
Return number of entries in dataset, i.e., count unweighted entries.
RooAbsGenContext is the abstract base class for generator contexts of RooAbsPdf objects.
virtual void setExpectedData(bool)
virtual RooDataSet * generate(double nEvents=0, bool skipInit=false, bool extendedMode=false)
Generate the specified number of events with nEvents>0 and and return a dataset containing the genera...
bool isValid() const
virtual void setProtoDataOrder(Int_t *lut)
Set the traversal order of prototype data to that in the lookup tables passed as argument.
Normalization set with for above integral.
Definition RooAbsPdf.h:370
~CacheElem() override
Destructor of normalization cache element.
RooAbsReal * _norm
Definition RooAbsPdf.h:375
RooArgSet _whatVars
Definition RooAbsPdf.h:86
RooAbsGenContext * _genContext
Definition RooAbsPdf.h:85
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...
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...
int calcSumW2CorrectedCovariance(RooMinimizer &minimizer, RooAbsReal &nll) const
Apply correction to errors and covariance matrix.
double getNorm(const RooArgSet &nset) const
Get normalisation term needed to normalise the raw values returned by getVal().
Definition RooAbsPdf.h:248
std::unique_ptr< RooAbsArg > compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext &ctx) const override
RooObjCacheManager _normMgr
Definition RooAbsPdf.h:377
std::unique_ptr< RooNumGenConfig > _specGeneratorConfig
! MC generator configuration specific for this object
Definition RooAbsPdf.h:389
void getLogProbabilities(RooSpan< const double > pdfValues, double *output) const
double getValV(const RooArgSet *set=nullptr) const override
Return current value, normalized by integrating over the observables in nset.
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.
bool _selectComp
Component selection flag for RooAbsPdf::plotCompOn.
Definition RooAbsPdf.h:387
virtual void generateEvent(Int_t code)
Interface for generation of an event using the algorithm corresponding to the specified code.
RooFit::OwningPtr< RooAbsReal > createScanCdf(const RooArgSet &iset, const RooArgSet &nset, Int_t numScanBins, Int_t intOrder)
void logBatchComputationErrors(RooSpan< const double > &outputs, std::size_t begin) const
Scan through outputs and fix+log all nans and negative values.
void setGeneratorConfig()
Remove the specialized numeric MC generator configuration associated with this object.
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...
static int verboseEval()
Return global level of verbosity for p.d.f. evaluations.
RooFit::OwningPtr< 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.
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:347
virtual double expectedEvents(const RooArgSet *nset) const
Return expected number of events to be used in calculation of extended likelihood.
virtual RooAbsGenContext * binnedGenContext(const RooArgSet &vars, bool verbose=false) const
Return a binned generator context.
TString _normRange
Normalization range.
Definition RooAbsPdf.h:391
virtual bool isDirectGenSafe(const RooAbsArg &arg) const
Check if given observable can be safely generated using the pdfs internal generator mechanism (if tha...
virtual RooFit::OwningPtr< RooAbsReal > createNLL(RooAbsData &data, const RooLinkedList &cmdList={})
Construct representation of -log(L) of PDF with given dataset.
Int_t * randomizeProtoOrder(Int_t nProto, Int_t nGen, bool resample=false) const
Return lookup table with randomized order for nProto prototype events.
virtual RooFit::OwningPtr< RooFitResult > fitTo(RooAbsData &data, const RooLinkedList &cmdList={})
Fit PDF to given dataset.
void setNormRange(const char *rangeName)
~RooAbsPdf() override
Destructor.
RooArgSet const * _normSet
Normalization integral (owned by _normMgr)
Definition RooAbsPdf.h:368
RooFitResult * chi2FitTo(RooDataHist &data, const RooLinkedList &cmdList) override
Calls RooAbsPdf::createChi2(RooDataSet& data, const RooLinkedList& cmdList) and returns fit result.
RooNumGenConfig * specialGeneratorConfig() const
Returns the specialized integrator configuration for this RooAbsReal.
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.
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:260
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Print multi line detailed information of this RooAbsPdf.
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:384
bool canBeExtended() const
If true, PDF can provide extended likelihood term.
Definition RooAbsPdf.h:271
RooAbsReal * _norm
Definition RooAbsPdf.h:367
int calcAsymptoticCorrectedCovariance(RooMinimizer &minimizer, RooAbsData const &data)
Use the asymptotically correct approach to estimate errors in the presence of weights.
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...
Int_t _errorCount
Number of errors remaining to print.
Definition RooAbsPdf.h:383
@ CanBeExtended
Definition RooAbsPdf.h:265
@ MustBeExtended
Definition RooAbsPdf.h:265
@ CanNotBeExtended
Definition RooAbsPdf.h:265
double _rawValue
Definition RooAbsPdf.h:366
Int_t _negCount
Number of negative probablities remaining to print.
Definition RooAbsPdf.h:385
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.
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 ...
void setActiveNormSet(RooArgSet const *normSet) const
Setter for the _normSet member, which should never be set directly.
Definition RooAbsPdf.h:335
double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName=nullptr) const override
Analytical integral with normalization (see RooAbsReal::analyticalIntegralWN() for further informatio...
void setNormRangeOverride(const char *rangeName)
virtual RooDataSet * generateSimGlobal(const RooArgSet &whatVars, Int_t nEvents)
Special generator interface for generation of 'global observables' – for RooStats tools.
double normalizeWithNaNPacking(double rawVal, double normVal) const
virtual RooAbsGenContext * autoGenContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr, const RooArgSet *auxProto=nullptr, bool verbose=false, bool autoBinned=true, const char *binnedTag="") const
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....
const RooNumGenConfig * getGeneratorConfig() const
Return the numeric MC generator configuration used for this object.
virtual void initGenerator(Int_t code)
Interface for one-time initialization to setup the generator for the specified code.
virtual ExtendMode extendMode() const
Returns ability of PDF to provide extended likelihood terms.
Definition RooAbsPdf.h:269
RooAbsPdf()
Default constructor.
bool traceEvalPdf(double value) const
Check that passed value is positive and not 'not-a-number'.
static RooNumGenConfig * defaultGeneratorConfig()
Returns the default numeric MC generator configuration for all RooAbsReals.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursiveStep) override
The cache manager.
void printValue(std::ostream &os) const override
Print value of p.d.f, also print normalization integral that was last used, if any.
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.
static TString _normRangeOverride
Definition RooAbsPdf.h:392
static Int_t _verboseEval
Definition RooAbsPdf.h:362
double extendedTerm(double sumEntries, double expected, double sumEntriesW2=0.0, bool doOffset=false) const
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...
virtual RooAbsPdf * createProjection(const RooArgSet &iset)
Return a p.d.f that represent a projection of this p.d.f integrated over given observables.
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...
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
bool hasRange(const char *name) const override
Check if variable has a binning with given name.
std::pair< double, double > getRange(const char *name=nullptr) const
Get low and high bound of the variable.
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition RooAbsReal.h:61
RooDataHist * fillDataHist(RooDataHist *hist, const RooArgSet *nset, double scaleFactor, bool correctForBinVolume=false, bool showProgress=false) const
Fill a RooDataHist with values sampled from this function at the bin centers.
void plotOnCompSelect(RooArgSet *selNodes) const
Helper function for plotting of composite p.d.fs.
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition RooAbsReal.h:105
RooFit::OwningPtr< RooAbsReal > createIntegral(const RooArgSet &iset, const RooCmdArg &arg1, 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
Create an object that represents the integral of the function over one or more observables std::liste...
bool plotSanityChecks(RooPlot *frame) const
Utility function for plotOn(), perform general sanity check on frame to ensure safe plotting operatio...
RooDerivative * derivative(RooRealVar &obs, Int_t order=1, double eps=0.001)
Return function representing first, second or third order derivative of this function.
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg(), const RooCmdArg &arg10=RooCmdArg()) const
Plot (project) PDF on specified frame.
RooFitResult * chi2FitDriver(RooAbsReal &fcn, RooLinkedList &cmdList)
Internal driver function for chi2 fits.
void printMultiline(std::ostream &os, Int_t contents, bool verbose=false, TString indent="") const override
Structure printing.
bool redirectServersHook(const RooAbsCollection &newServerList, bool mustReplaceAll, bool nameChange, bool isRecursiveStep) override
A buffer for reading values from trees.
double _value
Cache for current value of object.
Definition RooAbsReal.h:487
virtual double analyticalIntegral(Int_t code, const char *rangeName=nullptr) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
static void setEvalErrorLoggingMode(ErrorLoggingMode m)
Set evaluation error logging mode.
TString integralNameSuffix(const RooArgSet &iset, const RooArgSet *nset=nullptr, const char *rangeName=nullptr, bool omitEmpty=false) const
Construct std::string with unique suffix name to give to integral object that encodes integrated obse...
virtual double evaluate() const =0
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
void logEvalError(const char *message, const char *serverValueString=nullptr) const
Log evaluation error message.
const RooNumIntConfig * getIntegratorConfig() const
Return the numeric integration configuration used for this object.
virtual bool isBinnedDistribution(const RooArgSet &) const
Tests if the distribution is binned. Unless overridden by derived classes, this always returns false.
Definition RooAbsReal.h:359
RooFit::OwningPtr< RooAbsReal > createIntRI(const RooArgSet &iset, const RooArgSet &nset={})
Utility function for createRunningIntegral.
virtual double offset() const
Definition RooAbsReal.h:384
RooAddition calculates the sum of a set of RooAbsReal terms, or when constructed with two sets,...
Definition RooAddition.h:27
RooArgList is a container object that can hold multiple RooAbsArg objects.
Definition RooArgList.h:22
RooAbsArg * at(Int_t idx) const
Return object at given index, or nullptr if index is out of range.
Definition RooArgList.h:110
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition RooArgSet.h:55
RooBinnedGenContext is an efficient implementation of the generator context specific for binned pdfs.
Int_t setObj(const RooArgSet *nset, T *obj, const TNamed *isetRangeName=nullptr)
Setter function without integration set.
T * getObj(const RooArgSet *nset, Int_t *sterileIndex=nullptr, const TNamed *isetRangeName=nullptr)
Getter function without integration set.
RooCachedReal is an implementation of RooAbsCachedReal that can cache any external RooAbsReal input f...
void setCacheSource(bool flag)
RooChi2Var implements a simple calculation from a binned dataset and a PDF.
Definition RooChi2Var.h:25
RooCmdArg is a named container for two doubles, two integers two object points and three string point...
Definition RooCmdArg.h:26
void setInt(Int_t idx, Int_t value)
Definition RooCmdArg.h:72
Int_t getInt(Int_t idx) const
Definition RooCmdArg.h:86
void setString(Int_t idx, const char *value)
Definition RooCmdArg.h:78
Class RooCmdConfig is a configurable parser for RooCmdArg named arguments.
bool defineObject(const char *name, const char *argName, Int_t setNum, const TObject *obj=nullptr, bool isArray=false)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
double getDouble(const char *name, double defaultValue=0.0)
Return double property registered with name 'name'.
RooArgSet * getSet(const char *name, RooArgSet *set=nullptr)
Return RooArgSet property registered with name 'name'.
void defineMutex(const char *head, Args_t &&... tail)
Define arguments where any pair is mutually exclusive.
bool process(const RooCmdArg &arg)
Process given RooCmdArg.
bool hasProcessed(const char *cmdName) const
Return true if RooCmdArg with name 'cmdName' has been processed.
Int_t getInt(const char *name, Int_t defaultValue=0)
Return integer property registered with name 'name'.
static void stripCmdList(RooLinkedList &cmdList, const char *cmdsToPurge)
Utility function that strips command names listed (comma separated) in cmdsToPurge from cmdList.
bool defineInt(const char *name, const char *argName, Int_t intNum, Int_t defValue=0)
Define integer property name 'name' mapped to integer in slot 'intNum' in RooCmdArg with name argName...
bool defineString(const char *name, const char *argName, Int_t stringNum, const char *defValue="", bool appendMode=false)
Define double property name 'name' mapped to double in slot 'stringNum' in RooCmdArg with name argNam...
bool ok(bool verbose) const
Return true of parsing was successful.
TObject * getObject(const char *name, TObject *obj=nullptr)
Return TObject property registered with name 'name'.
const char * getString(const char *name, const char *defaultValue="", bool convEmptyToNull=false)
Return string property registered with name 'name'.
void allowUndefined(bool flag=true)
If flag is true the processing of unrecognized RooCmdArgs is not considered an error.
bool defineSet(const char *name, const char *argName, Int_t setNum, const RooArgSet *set=nullptr)
Define TObject property name 'name' mapped to object in slot 'setNum' in RooCmdArg with name argName ...
RooLinkedList filterCmdList(RooLinkedList &cmdInList, const char *cmdNameList, bool removeFromInList=true) const
Utility function to filter commands listed in cmdNameList from cmdInList.
bool defineDouble(const char *name, const char *argName, Int_t doubleNum, double defValue=0.0)
Define double property name 'name' mapped to double in slot 'doubleNum' in RooCmdArg with name argNam...
The RooDataHist is a container class to hold N-dimensional binned data.
Definition RooDataHist.h:39
double weight(std::size_t i) const
Return weight of i-th bin.
void set(std::size_t binNumber, double weight, double wgtErr)
Set bin content of bin that was last loaded with get(std::size_t).
const RooArgSet * get() const override
Get bin centre of current bin.
Definition RooDataHist.h:76
double sumEntries() const override
Sum the weights of all bins.
RooDataSet is a container class to hold unbinned data.
Definition RooDataSet.h:57
RooFitResult is a container class to hold the input and output of a PDF fit to a dataset.
void compileServers(RooAbsArg &arg, RooArgSet const &normSet)
Class RooGenContext implement a universal generator context for all RooAbsPdf classes that do not hav...
Switches the message service to a different level while the instance is alive.
Definition RooHelpers.h:42
RooLinkedList is an collection class for internal use, storing a collection of RooAbsArg pointers in ...
virtual void Add(TObject *arg)
TObject * FindObject(const char *name) const override
Return pointer to obejct with given name.
RooMinimizer is a wrapper class around ROOT::Fit:Fitter that provides a seamless interface between th...
int hesse()
Execute HESSE.
RooFitResult * save(const char *name=nullptr, const char *title=nullptr)
Save and return a RooFitResult snapshot of current minimizer status.
void applyCovarianceMatrix(TMatrixDSym const &V)
Apply results of given external covariance matrix.
static const char * str(const TNamed *ptr)
Return C++ string corresponding to given TNamed pointer.
Definition RooNameReg.h:37
RooNumGenConfig holds the configuration parameters of the various numeric integrators used by RooReal...
static RooNumGenConfig & defaultConfig()
Return reference to instance of default numeric integrator configuration object.
Class RooObjCacheManager is an implementation of class RooCacheManager<RooAbsCacheElement> and specia...
void sterilize() override
Clear the cache payload but retain slot mapping w.r.t to normalization and integration sets.
A RooPlot is a plot frame and a container for graphics objects within that frame.
Definition RooPlot.h:43
void addObject(TObject *obj, Option_t *drawOptions="", bool invisible=false)
Add a generic object to this plot.
Definition RooPlot.cxx:383
double getFitRangeNEvt() const
Return the number of events in the fit range.
Definition RooPlot.h:139
const RooArgSet * getNormVars() const
Definition RooPlot.h:146
RooAbsRealLValue * getPlotVar() const
Definition RooPlot.h:137
void updateNormVars(const RooArgSet &vars)
Install the given set of observables are reference normalization variables for this frame.
Definition RooPlot.cxx:368
double getFitRangeBinW() const
Return the bin width that is being used to normalise the PDF.
Definition RooPlot.h:142
virtual void printStream(std::ostream &os, Int_t contents, StyleOption style, TString indent="") const
Print description of object on ostream, printing contents set by contents integer,...
Class RooProjectedPdf is a RooAbsPdf implementation that represent a projection of a given input p....
static UInt_t integer(UInt_t max, TRandom *generator=randomGenerator())
Return an integer uniformly distributed from [0,n-1].
Definition RooRandom.cxx:99
static TRandom * randomGenerator()
Return a pointer to a singleton random-number generator implementation.
Definition RooRandom.cxx:51
RooRealIntegral performs hybrid numerical/analytical integrals of RooAbsReal objects.
const RooArgSet & numIntRealVars() const
RooRealVar represents a variable that can be changed from the outside.
Definition RooRealVar.h:40
static TClass * Class()
void setRange(const char *name, double min, double max)
Set a fit or plotting range.
void setBins(Int_t nBins, const char *name=nullptr)
Create a uniform binning under name 'name' for this variable.
Container_t::const_iterator begin() const
Iterator over contained objects.
Container_t::const_iterator end() const
End of contained objects.
std::size_t size() const
Number of contained objects (neglecting the ref count).
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
constexpr std::size_t size() const noexcept
Definition RooSpan.h:119
RooXYChi2Var implements a simple chi^2 calculation from an unbinned dataset with values x,...
Int_t GetNrows() const
TMatrixTSym< Element > & Similarity(const TMatrixT< Element > &n)
Calculate B * (*this) * B^T , final matrix will be (nrowsb x nrowsb) This is a similarity transform w...
The TNamed class is the base class for all named ROOT classes.
Definition TNamed.h:29
const char * GetName() const override
Returns name of object.
Definition TNamed.h:47
const char * GetTitle() const override
Returns title of object.
Definition TNamed.h:48
Mother of all ROOT objects.
Definition TObject.h:41
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition TObject.cxx:207
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition TObject.cxx:525
A Pave (see TPave) with text, lines or/and boxes inside.
Definition TPaveText.h:21
virtual Int_t Poisson(Double_t mean)
Generates a random integer N according to a Poisson law.
Definition TRandom.cxx:402
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition TRandom.cxx:672
virtual UInt_t Integer(UInt_t imax)
Returns a random integer uniformly distributed on the interval [ 0, imax-1 ].
Definition TRandom.cxx:360
Basic string class.
Definition TString.h:139
Ssiz_t Length() const
Definition TString.h:421
void Clear()
Clear string without changing its capacity.
Definition TString.cxx:1222
const char * Data() const
Definition TString.h:380
TLine * line
void box(Int_t pat, Double_t x1, Double_t y1, Double_t x2, Double_t y2)
Definition fillpatterns.C:1
RooCmdArg SupNormSet(const RooArgSet &nset)
RooCmdArg WeightVar(const char *name="weight", bool reinterpretAsWeight=false)
RooCmdArg Hesse(bool flag=true)
RooCmdArg ModularL(bool flag=false)
RooCmdArg PrintLevel(Int_t code)
RooCmdArg BatchMode(std::string const &batchMode="cpu")
RooCmdArg NormRange(const char *rangeNameList)
RooCmdArg Range(const char *rangeName, bool adjustNorm=true)
RooCmdArg Normalization(double scaleFactor)
std::vector< std::string > Split(std::string_view str, std::string_view delims, bool skipEmpty=false)
Splits a string at each character in delims.
std::unique_ptr< RooAbsReal > createNLL(RooAbsPdf &pdf, RooAbsData &data, std::unique_ptr< RooAbsReal > &&constraints, std::string const &rangeName, RooArgSet const &projDeps, bool isExtended, double integrateOverBinsPrecision, RooFit::OffsetMode offset)
OwningPtr< T > owningPtr(std::unique_ptr< T > &&ptr)
Internal helper to turn a std::unique_ptr<T> into an OwningPtr.
Definition Config.h:50
std::unique_ptr< T > compileForNormSet(T const &arg, RooArgSet const &normSet)
std::string & defaultBatchMode()
Get a handle on the default BatchMode option that is used when creating likelihoods.
OffsetMode
For setting the offset mode with the Offset() command argument to RooAbsPdf::fitTo()
BatchModeOption
For setting the batch mode flag with the BatchMode() command argument to RooAbsPdf::fitTo()
T * OwningPtr
An alias for raw pointers for indicating that the return type of a RooFit function is an owning point...
Definition Config.h:43
std::unique_ptr< T > cloneTreeWithSameParameters(T const &arg, RooArgSet const *observables=nullptr)
Clone RooAbsArg object and reattach to original parameters.
Definition RooHelpers.h:149
RooArgSet selectFromArgSet(RooArgSet const &, std::string const &names)
Construct a RooArgSet of objects in a RooArgSet whose names match to those in the names string.
std::string getColonSeparatedNameString(RooArgSet const &argSet)
Create a string with all sorted names of RooArgSet elements separated by colons.
Bool_t IsNaN(Double_t x)
Definition TMath.h:892
Double_t QuietNaN()
Returns a quiet NaN as defined by IEEE 754.
Definition TMath.h:902
Definition first.py:1
Configuration struct for RooAbsPdf::minimizeNLL with all the default.
Definition RooAbsPdf.h:171
const RooArgSet * minosSet
Definition RooAbsPdf.h:193
std::string rangeName
Stores the configuration parameters for RooAbsTestStatistic.
Config argument to RooMinimizer ctor.
__roodevice__ static __roohost__ double packFloatIntoNaN(float payload)
Pack float into mantissa of a NaN.
TMarker m
Definition textangle.C:8
TLine l
Definition textangle.C:4
static void output()