From 6980f6f5e650ce83443f527d4f3c20d66b5e1f55 Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Tue, 19 Nov 2024 17:35:24 +0100 Subject: [PATCH] Add public getters to useful information in combine classes This will help for example in the serialization to JSON or the C++ code generation for automatic differentiation. --- interface/AsymPow.h | 6 +++-- interface/ProcessNormalization.h | 8 ++++++- interface/VerticalInterpHistPdf.h | 40 ++++++++++++++++++++++--------- src/AsymPow.cc | 7 ------ src/ProcessNormalization.cc | 14 ++++------- 5 files changed, 44 insertions(+), 31 deletions(-) diff --git a/interface/AsymPow.h b/interface/AsymPow.h index 8346bdadba6..6cda05c5522 100644 --- a/interface/AsymPow.h +++ b/interface/AsymPow.h @@ -26,9 +26,11 @@ class AsymPow : public RooAbsReal { AsymPow(const char *name, const char *title, RooAbsReal &kappaLow, RooAbsReal &kappaHigh, RooAbsReal &theta) ; AsymPow(const AsymPow &other, const char *newname=0) ; - ~AsymPow() override ; + TObject * clone(const char *newname) const override { return new AsymPow(*this, newname); } - TObject * clone(const char *newname) const override ; + RooAbsReal const &kappaLow() const { return kappaLow_.arg(); } + RooAbsReal const &kappaHigh() const { return kappaHigh_.arg(); } + RooAbsReal const &theta() const { return theta_.arg(); } protected: Double_t evaluate() const override; diff --git a/interface/ProcessNormalization.h b/interface/ProcessNormalization.h index 24916cc3792..edf2ec3f15c 100644 --- a/interface/ProcessNormalization.h +++ b/interface/ProcessNormalization.h @@ -19,7 +19,6 @@ class ProcessNormalization : public RooAbsReal { ProcessNormalization(const char *name, const char *title, double nominal=1) ; ProcessNormalization(const char *name, const char *title, RooAbsReal &nominal) ; ProcessNormalization(const ProcessNormalization &other, const char *newname = 0) ; - ~ProcessNormalization() override ; TObject * clone(const char *newname) const override { return new ProcessNormalization(*this, newname); } @@ -29,6 +28,13 @@ class ProcessNormalization : public RooAbsReal { void addOtherFactor(RooAbsReal &factor) ; void dump() const ; + double nominalValue() const { return nominalValue_; } + std::vector const &logKappa() const { return logKappa_; } + RooArgList const &thetaList() const { return thetaList_; } + std::vector> const &logAsymmKappa() const { return logAsymmKappa_; } + RooArgList const &asymmThetaList() const { return asymmThetaList_; } + RooArgList const &otherFactorList() const { return otherFactorList_; } + protected: Double_t evaluate() const override; diff --git a/interface/VerticalInterpHistPdf.h b/interface/VerticalInterpHistPdf.h index f1172b6f6d0..c90fdc051ff 100644 --- a/interface/VerticalInterpHistPdf.h +++ b/interface/VerticalInterpHistPdf.h @@ -34,7 +34,12 @@ class VerticalInterpHistPdf : public RooAbsPdf { const RooArgList& coefList() const { return _coefList ; } const RooRealProxy &x() const { return _x; } + + double smoothRegion() const { return _smoothRegion; } + int smoothAlgo() const { return _smoothAlgo; } + friend class FastVerticalInterpHistPdf; + protected: RooRealProxy _x; RooListProxy _funcList ; // List of component FUNCs @@ -84,6 +89,11 @@ class FastVerticalInterpHistPdfBase : public RooAbsPdf { /// Must be public, for serialization struct Morph { FastTemplate sum; FastTemplate diff; }; + std::vector const &morphs() const { return _morphs; } + + double smoothRegion() const { return _smoothRegion; } + int smoothAlgo() const { return _smoothAlgo; } + friend class FastVerticalInterpHistPdf2Base; protected: RooRealProxy _x; @@ -127,21 +137,22 @@ class FastVerticalInterpHistPdf : public FastVerticalInterpHistPdfBase { FastVerticalInterpHistPdf() : FastVerticalInterpHistPdfBase() {} FastVerticalInterpHistPdf(const char *name, const char *title, const RooRealVar &x, const RooArgList& funcList, const RooArgList& coefList, Double_t smoothRegion=1., Int_t smoothAlgo=1) : FastVerticalInterpHistPdfBase(name,title,RooArgSet(x),funcList,coefList,smoothRegion,smoothAlgo), - _x("x","Independent variable",this,const_cast(x)), - _cache(), _cacheNominal(), _cacheNominalLog() {} + _x("x","Independent variable",this,const_cast(x)) {} FastVerticalInterpHistPdf(const FastVerticalInterpHistPdf& other, const char* name=0) : FastVerticalInterpHistPdfBase(other, name), _x("x",this,other._x), _cache(other._cache), _cacheNominal(other._cacheNominal), _cacheNominalLog(other._cacheNominalLog) {} TObject* clone(const char* newname) const override { return new FastVerticalInterpHistPdf(*this,newname) ; } - ~FastVerticalInterpHistPdf() override {} Double_t evaluate() const override ; Bool_t hasCache() const { return _cache.size() > 0; } Bool_t isCacheReady() const { return _cache.size() > 0 && _init; } + FastHisto const &cacheNominal() const { return _cacheNominal; } + friend class FastVerticalInterpHistPdfV; friend class FastVerticalInterpHistPdf2; + protected: RooRealProxy _x; @@ -186,14 +197,12 @@ class FastVerticalInterpHistPdf2D : public FastVerticalInterpHistPdfBase { FastVerticalInterpHistPdfBase(name,title,RooArgSet(x,y),funcList,coefList,smoothRegion,smoothAlgo), _x("x","Independent variable",this,const_cast(x)), _y("y","Independent variable",this,const_cast(y)), - _conditional(conditional), - _cache(), _cacheNominal(), _cacheNominalLog() {} + _conditional(conditional) {} FastVerticalInterpHistPdf2D(const FastVerticalInterpHistPdf2D& other, const char* name=0) : FastVerticalInterpHistPdfBase(other, name), _x("x",this,other._x), _y("y",this,other._y), _conditional(other._conditional), _cache(other._cache), _cacheNominal(other._cacheNominal), _cacheNominalLog(other._cacheNominalLog) {} TObject* clone(const char* newname) const override { return new FastVerticalInterpHistPdf2D(*this,newname) ; } - ~FastVerticalInterpHistPdf2D() override {} const RooRealVar & x() const { return dynamic_cast(_x.arg()); } const RooRealVar & y() const { return dynamic_cast(_y.arg()); } @@ -204,6 +213,8 @@ class FastVerticalInterpHistPdf2D : public FastVerticalInterpHistPdfBase { Bool_t hasCache() const { return _cache.size() > 0; } Bool_t isCacheReady() const { return _cache.size() > 0 && _init; } + FastHisto2D const &cacheNominal() const { return _cacheNominal; } + friend class FastVerticalInterpHistPdf2D2; protected: RooRealProxy _x, _y; @@ -245,8 +256,15 @@ class FastVerticalInterpHistPdf2Base : public RooAbsPdf { virtual void setActiveBins(unsigned int bins) {} bool cacheIsGood() const { return _sentry.good() && _initBase; } + + double smoothRegion() const { return _smoothRegion; } + int smoothAlgo() const { return _smoothAlgo; } + /// Must be public, for serialization typedef FastVerticalInterpHistPdfBase::Morph Morph; + + std::vector const &morphs() const { return _morphs; } + protected: RooListProxy _coefList ; // List of coefficients Double_t _smoothRegion; @@ -302,7 +320,6 @@ class FastVerticalInterpHistPdf2 : public FastVerticalInterpHistPdf2Base { _cache(other._cache), _cacheNominal(other._cacheNominal), _cacheNominalLog(other._cacheNominalLog) {} explicit FastVerticalInterpHistPdf2(const FastVerticalInterpHistPdf& other, const char* name=0) ; TObject* clone(const char* newname) const override { return new FastVerticalInterpHistPdf2(*this,newname) ; } - ~FastVerticalInterpHistPdf2() override {} const RooRealVar & x() const { return dynamic_cast(_x.arg()); } @@ -311,6 +328,8 @@ class FastVerticalInterpHistPdf2 : public FastVerticalInterpHistPdf2Base { FastHisto const& cache() const { return _cache; } + FastHisto const &cacheNominal() const { return _cacheNominal; } + friend class FastVerticalInterpHistPdf2V; protected: RooRealProxy _x; @@ -358,7 +377,6 @@ class FastVerticalInterpHistPdf2D2 : public FastVerticalInterpHistPdf2Base { _cache(other._cache), _cacheNominal(other._cacheNominal), _cacheNominalLog(other._cacheNominalLog) {} explicit FastVerticalInterpHistPdf2D2(const FastVerticalInterpHistPdf2D& other, const char* name=0) ; TObject* clone(const char* newname) const override { return new FastVerticalInterpHistPdf2D2(*this,newname) ; } - ~FastVerticalInterpHistPdf2D2() override {} const RooRealVar & x() const { return dynamic_cast(_x.arg()); } const RooRealVar & y() const { return dynamic_cast(_y.arg()); } @@ -369,6 +387,8 @@ class FastVerticalInterpHistPdf2D2 : public FastVerticalInterpHistPdf2Base { Double_t evaluate() const override ; + FastHisto2D const &cacheNominal() const { return _cacheNominal; } + protected: RooRealProxy _x, _y; bool _conditional; @@ -400,14 +420,12 @@ class FastVerticalInterpHistPdf3D : public FastVerticalInterpHistPdfBase { _x("x","Independent variable",this,const_cast(x)), _y("y","Independent variable",this,const_cast(y)), _z("z","Independent variable",this,const_cast(z)), - _conditional(conditional), - _cache(), _cacheNominal(), _cacheNominalLog() {} + _conditional(conditional) {} FastVerticalInterpHistPdf3D(const FastVerticalInterpHistPdf3D& other, const char* name=0) : FastVerticalInterpHistPdfBase(other, name), _x("x",this,other._x), _y("y",this,other._y), _z("z",this,other._z), _conditional(other._conditional), _cache(other._cache), _cacheNominal(other._cacheNominal), _cacheNominalLog(other._cacheNominalLog) {} TObject* clone(const char* newname) const override { return new FastVerticalInterpHistPdf3D(*this,newname) ; } - ~FastVerticalInterpHistPdf3D() override {} const RooRealVar & x() const { return dynamic_cast(_x.arg()); } const RooRealVar & y() const { return dynamic_cast(_y.arg()); } diff --git a/src/AsymPow.cc b/src/AsymPow.cc index fb7ba951c13..1a3ed897e20 100644 --- a/src/AsymPow.cc +++ b/src/AsymPow.cc @@ -20,13 +20,6 @@ AsymPow::AsymPow(const AsymPow &other, const char *newname) : theta_("theta",this,other.theta_) { } -AsymPow::~AsymPow() {} - -TObject *AsymPow::clone(const char *newname) const -{ - return new AsymPow(*this,newname); -} - Double_t AsymPow::evaluate() const { return RooFit::Detail::MathFuncs::asymPow(theta_, kappaLow_, kappaHigh_); } diff --git a/src/ProcessNormalization.cc b/src/ProcessNormalization.cc index 5baa229b3a1..b1fbb936388 100644 --- a/src/ProcessNormalization.cc +++ b/src/ProcessNormalization.cc @@ -15,14 +15,10 @@ ProcessNormalization::ProcessNormalization(const char *name, const char *title, { } -ProcessNormalization::ProcessNormalization(const char *name, const char *title, RooAbsReal &nominal) : - RooAbsReal(name,title), - nominalValue_(1.0), - thetaList_("thetaList", "List of nuisances for symmetric kappas", this), - asymmThetaList_("asymmThetaList", "List of nuisances for asymmetric kappas", this), - otherFactorList_("otherFactorList", "Other multiplicative terms", this) -{ - otherFactorList_.add(nominal); +ProcessNormalization::ProcessNormalization(const char *name, const char *title, RooAbsReal &nominal) + : ProcessNormalization{name, title, 1.0} +{ + otherFactorList_.add(nominal); } ProcessNormalization::ProcessNormalization(const ProcessNormalization &other, const char *newname) : @@ -36,8 +32,6 @@ ProcessNormalization::ProcessNormalization(const ProcessNormalization &other, co { } -ProcessNormalization::~ProcessNormalization() {} - void ProcessNormalization::addLogNormal(double kappa, RooAbsReal &theta) { if (kappa != 0.0 && kappa != 1.0) { logKappa_.push_back(std::log(kappa));