Logo ROOT  
Reference Guide
RooRealL.cxx
Go to the documentation of this file.
1/*
2 * Project: RooFit
3 * Authors:
4 * PB, Patrick Bos, Netherlands eScience Center, p.bos@esciencecenter.nl
5 *
6 * Copyright (c) 2021, CERN
7 *
8 * Redistribution and use in source and binary forms,
9 * with or without modification, are permitted according to the terms
10 * listed in LICENSE (http://roofit.sourceforge.net/license.txt)
11 */
12
15#include <RooRealVar.h>
16
17namespace RooFit {
18namespace TestStatistics {
19
20/** \class RooRealL
21 * \ingroup Roofitcore
22 *
23 * \brief RooAbsReal that wraps RooAbsL likelihoods for use in RooFit outside of the RooMinimizer context
24 *
25 * This class provides a simple wrapper to evaluate RooAbsL derived likelihood objects like a regular RooFit real value.
26 * Whereas the RooAbsL objects are meant to be used within the context of minimization, RooRealL can be used in any
27 * RooFit context, like plotting. The value can be accessed through getVal(), like with other RooFit real variables.
28 **/
29
30RooRealL::RooRealL(const char *name, const char *title, std::shared_ptr<RooAbsL> likelihood)
31 : RooAbsReal(name, title), likelihood_(std::move(likelihood)),
32 vars_proxy_("varsProxy", "proxy set of parameters", this)
33{
34 vars_obs_.add(*likelihood_->getParameters());
35 vars_proxy_.add(*likelihood_->getParameters());
36}
37
38RooRealL::RooRealL(const RooRealL &other, const char *name)
39 : RooAbsReal(other, name), likelihood_(other.likelihood_), vars_proxy_("varsProxy", this, other.vars_proxy_)
40{
41 vars_obs_.add(other.vars_obs_) ;
42}
43
44double RooRealL::evaluate() const
45{
46 // Transfer values from proxy variables to internal variables of likelihood
47 if (!vars_proxy_.empty()) {
48 for (auto i = 0u; i < vars_obs_.size(); ++i) {
49 auto harg = vars_obs_[i];
50 const auto parg = vars_proxy_[i];
51
52 if (harg != parg) {
53 ((RooAbsRealLValue*)harg)->setVal(((RooAbsReal*)parg)->getVal());
54 }
55 }
56 }
57 // Evaluate as straight FUNC
58 std::size_t last_component = likelihood_->getNComponents();
59
60 auto ret_kahan = likelihood_->evaluatePartition({0, 1}, 0, last_component);
61
62 const double norm = globalNormalization();
63 double ret = ret_kahan.Sum() / norm;
64 eval_carry = ret_kahan.Carry() / norm;
65
66 return ret;
67}
68
69} // namespace TestStatistics
70} // namespace RooFit
char name[80]
Definition: TGX11.cxx:110
bool empty() const
virtual bool add(const RooAbsArg &var, bool silent=false)
Add the specified argument to list.
Storage_t::size_type size() const
RooAbsRealLValue is the common abstract base class for objects that represent a real value that may a...
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:64
double getVal(const RooArgSet *normalisationSet=nullptr) const
Evaluate object.
Definition: RooAbsReal.h:94
bool add(const RooAbsArg &var, bool valueServer, bool shapeServer, bool silent)
Overloaded RooCollection_t::add() method insert object into set and registers object as server to own...
RooAbsReal that wraps RooAbsL likelihoods for use in RooFit outside of the RooMinimizer context.
Definition: RooRealL.h:28
std::shared_ptr< RooAbsL > likelihood_
Definition: RooRealL.h:48
double evaluate() const override
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
Definition: RooRealL.cxx:44
RooSetProxy vars_proxy_
sets up client-server connections
Definition: RooRealL.h:50
double globalNormalization() const
Definition: RooRealL.h:35
RooArgSet vars_obs_
list of observables
Definition: RooRealL.h:51
RooRealL(const char *name, const char *title, std::shared_ptr< RooAbsL > likelihood)
Definition: RooRealL.cxx:30
The namespace RooFit contains mostly switches that change the behaviour of functions of PDFs (or othe...
Definition: Common.h:18