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