From bcc3b3df378999612ad14eadcd38f9db676fb5ed Mon Sep 17 00:00:00 2001 From: Rahul Balasubramanian Date: Wed, 16 Feb 2022 22:04:57 +0100 Subject: [PATCH 1/6] [RF] Cleanup RooLagrangianMorphFunc --- roofit/roofit/inc/RooLagrangianMorphFunc.h | 45 +++--- roofit/roofit/src/RooLagrangianMorphFunc.cxx | 144 ++++++------------- 2 files changed, 65 insertions(+), 124 deletions(-) diff --git a/roofit/roofit/inc/RooLagrangianMorphFunc.h b/roofit/roofit/inc/RooLagrangianMorphFunc.h index 5622720e2b45f..3ce2ca4f11943 100644 --- a/roofit/roofit/inc/RooLagrangianMorphFunc.h +++ b/roofit/roofit/inc/RooLagrangianMorphFunc.h @@ -96,7 +96,7 @@ class RooLagrangianMorphFunc : public RooAbsReal { RooArgList prodCouplings; RooArgList folders; std::vector vertices; - std::vector nonInterfering; + std::vector> nonInterfering; bool allowNegativeYields = true; }; @@ -108,24 +108,21 @@ class RooLagrangianMorphFunc : public RooAbsReal { virtual ~RooLagrangianMorphFunc(); - virtual std::list * - binBoundaries(RooAbsRealLValue & /*obs*/, Double_t /*xlo*/, Double_t /*xhi*/) const override; - virtual std::list * - plotSamplingHint(RooAbsRealLValue & /*obs*/, Double_t /*xlo*/, Double_t /*xhi*/) const override; - virtual Bool_t isBinnedDistribution(const RooArgSet &obs) const override; - virtual Double_t evaluate() const override; - virtual TObject *clone(const char *newname) const override; - virtual Double_t getValV(const RooArgSet *set = 0) const override; - - virtual Bool_t checkObservables(const RooArgSet *nset) const override; - virtual Bool_t forceAnalyticalInt(const RooAbsArg &arg) const override; - virtual Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, - const char *rangeName = 0) const override; - virtual Double_t - analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName = 0) const override; - virtual void printMetaArgs(std::ostream &os) const override; - virtual RooAbsArg::CacheMode canNodeBeCached() const override; - virtual void setCacheAndTrackHints(RooArgSet &) override; + std::list *binBoundaries(RooAbsRealLValue & /*obs*/, double /*xlo*/, double /*xhi*/) const override; + std::list *plotSamplingHint(RooAbsRealLValue & /*obs*/, double /*xlo*/, double /*xhi*/) const override; + bool isBinnedDistribution(const RooArgSet &obs) const override; + double evaluate() const override; + TObject *clone(const char *newname) const override; + double getValV(const RooArgSet *set = 0) const override; + + bool checkObservables(const RooArgSet *nset) const override; + bool forceAnalyticalInt(const RooAbsArg &arg) const override; + Int_t getAnalyticalIntegralWN(RooArgSet &allVars, RooArgSet &numVars, const RooArgSet *normSet, + const char *rangeName = 0) const override; + double analyticalIntegralWN(Int_t code, const RooArgSet *normSet, const char *rangeName = 0) const override; + void printMetaArgs(std::ostream &os) const override; + RooAbsArg::CacheMode canNodeBeCached() const override; + void setCacheAndTrackHints(RooArgSet &) override; void insert(RooWorkspace *ws); @@ -178,11 +175,13 @@ class RooLagrangianMorphFunc : public RooAbsReal { TH1 *createTH1(const std::string &name); TH1 *createTH1(const std::string &name, bool correlateErrors); -protected: +private: class CacheElem; void init(); void setup(bool ownParams = true); bool _ownParameters = false; + void disableInterference(const std::vector &nonInterfering); + void disableInterferences(const std::vector> &nonInterfering); mutable RooObjCacheManager _cacheMgr; //! The cache manager @@ -206,7 +205,6 @@ class RooLagrangianMorphFunc : public RooAbsReal { int countSamples(std::vector &vertices); int countSamples(int nprod, int ndec, int nboth); - TPair *makeCrosssectionContainer(double xs, double unc); std::map createWeightStrings(const ParamMap &inputs, const std::vector> &vertices); std::map @@ -252,7 +250,7 @@ class RooLagrangianMorphFunc : public RooAbsReal { static std::unique_ptr makeRatio(const char *name, const char *title, RooArgList &nr, RooArgList &dr); -protected: +private: double _scale = 1; std::map _sampleMap; RooListProxy _physics; @@ -262,7 +260,8 @@ class RooLagrangianMorphFunc : public RooAbsReal { RooListProxy _flags; Config _config; std::vector> _diagrams; - mutable const RooArgSet *_curNormSet; //! + mutable const RooArgSet *_curNormSet = nullptr; //! + std::vector _nonInterfering; ClassDefOverride(RooLagrangianMorphFunc, 1) diff --git a/roofit/roofit/src/RooLagrangianMorphFunc.cxx b/roofit/roofit/src/RooLagrangianMorphFunc.cxx index 7f3bb9c2afc72..3530cddb9f38c 100644 --- a/roofit/roofit/src/RooLagrangianMorphFunc.cxx +++ b/roofit/roofit/src/RooLagrangianMorphFunc.cxx @@ -638,43 +638,6 @@ inline void extractCouplings(const T1 &inCouplings, T2 &outCouplings) } } -/////////////////////////////////////////////////////////////////////////////// -/// find and, if necessary, create a parameter from a list - -template -inline RooAbsArg &get(T &operators, const char *name, double defaultval = 0) -{ - RooAbsArg *kappa = operators.find(name); - if (kappa) - return *kappa; - RooRealVar *newKappa = new RooRealVar(name, name, defaultval); - newKappa->setConstant(false); - operators.add(*newKappa); - return *newKappa; -} - -/////////////////////////////////////////////////////////////////////////////// -/// find and, if necessary, create a parameter from a list - -template -inline RooAbsArg &get(T &operators, const std::string &name, double defaultval = 0) -{ - return get(operators, name.c_str(), defaultval); -} - -//////////////////////////////////////////////////////////////////////////////// -/// create a new coupling and add it to the set - -template -inline void addCoupling(T &set, const TString &name, const TString &formula, const RooArgList &components, bool isNP) -{ - if (!set.find(name)) { - RooFormulaVar *c = new RooFormulaVar(name, formula, components); - c->setAttribute("NewPhysics", isNP); - set.add(*c); - } -} - //////////////////////////////////////////////////////////////////////////////// /// set parameter values first set all values to defaultVal (if value not /// present in param_card then it should be 0) @@ -1396,9 +1359,8 @@ class RooLagrangianMorphFunc::CacheElem : public RooAbsCacheElement { extractCouplings(*vertex, this->_couplings); } } - extractOperators(this->_couplings, operators); - this->_formulas = - ::createFormulas(funcname, inputParameters, inputFlags, diagrams, this->_couplings, flags, nonInterfering); + extractOperators(_couplings, operators); + _formulas = ::createFormulas(funcname, inputParameters, inputFlags, diagrams, _couplings, flags, nonInterfering); } ////////////////////////////////////////////////////////////////////////////// @@ -1771,18 +1733,15 @@ void RooLagrangianMorphFunc::collectInputs(TDirectory *file) RooRealVar *observable = this->setupObservable(obsName.c_str(), mode, obj); if (classname.find("TH1") != std::string::npos) { - collectHistograms(this->GetName(), file, this->_sampleMap, this->_physics, *observable, obsName, - this->_config.paramCards); + collectHistograms(this->GetName(), file, _sampleMap, _physics, *observable, obsName, _config.paramCards); } else if (classname.find("RooHistFunc") != std::string::npos || classname.find("RooParamHistFunc") != std::string::npos || classname.find("PiecewiseInterpolation") != std::string::npos) { collectRooAbsReal(this->GetName(), file, this->_sampleMap, this->_physics, obsName, this->_config.paramCards); } else if (classname.find("TParameter") != std::string::npos) { - collectCrosssections(this->GetName(), file, this->_sampleMap, this->_physics, obsName, - this->_config.paramCards); + collectCrosssections(this->GetName(), file, _sampleMap, _physics, obsName, _config.paramCards); } else if (classname.find("TParameter") != std::string::npos) { - collectCrosssections(this->GetName(), file, this->_sampleMap, this->_physics, obsName, - this->_config.paramCards); + collectCrosssections(this->GetName(), file, _sampleMap, _physics, obsName, _config.paramCards); } else if (classname.find("TPair") != std::string::npos) { collectCrosssectionsTPair(this->GetName(), file, this->_sampleMap, this->_physics, obsName, folderNames[0], this->_config.paramCards); @@ -1868,6 +1827,7 @@ RooLagrangianMorphFunc::RooLagrangianMorphFunc(const char *name, const char *tit _binWidths("binWidths", "set of binWidth objects", this, kTRUE, kFALSE), _config(config), _curNormSet(0) { this->init(); + this->disableInterferences(_config.nonInterfering); this->setup(false); TRACE_CREATE @@ -1899,42 +1859,9 @@ RooLagrangianMorphFunc::RooLagrangianMorphFunc(const char *name, const char *tit void RooLagrangianMorphFunc::setup(bool own) { - this->_ownParameters = own; - auto diagrams = this->_diagrams; - - if (diagrams.size() > 0) { - RooArgList operators; - for (auto const &v : diagrams) { - for (const RooArgList *t : v) { - extractOperators(*t, operators); - } - } - - if (own) { - this->_operators.addOwned(operators); - } else { - this->_operators.add(operators); - } - - for (size_t j = 0; j < diagrams.size(); ++j) { - std::vector diagram; - for (size_t i = 0; i < diagrams[j].size(); ++i) { - std::stringstream name; - name << "!vertex" << i; - std::stringstream title; - title << "set of couplings in the vertex " << i; - diagram.push_back(new RooListProxy(name.str().c_str(), title.str().c_str(), this, kTRUE, kFALSE)); - if (own) { - diagram[i]->addOwned(*diagrams[j][i]); - } else { - diagram[i]->add(*diagrams[j][i]); - } - } - this->_diagrams.push_back(diagram); - } - } + _ownParameters = own; - else if (this->_config.couplings.size() > 0) { + if (_config.couplings.size() > 0) { RooArgList operators; std::vector vertices; extractOperators(this->_config.couplings, operators); @@ -1973,7 +1900,36 @@ void RooLagrangianMorphFunc::setup(bool own) } //////////////////////////////////////////////////////////////////////////////// -/// (-?-) +/// disable interference between terms + +void RooLagrangianMorphFunc::disableInterference(const std::vector &nonInterfering) +{ + // disable interference between the listed operators + std::stringstream name; + name << "noInteference"; + for (auto c : nonInterfering) { + name << c; + } + RooListProxy *p = new RooListProxy(name.str().c_str(), name.str().c_str(), this, kTRUE, kFALSE); + this->_nonInterfering.push_back(p); + for (auto c : nonInterfering) { + p->addOwned(*(new RooStringVar(c, c, c))); + } +} + +//////////////////////////////////////////////////////////////////////////////// +/// disable interference between terms + +void RooLagrangianMorphFunc::disableInterferences(const std::vector> &nonInterfering) +{ + // disable interferences between the listed groups of operators + for (size_t i = 0; i < nonInterfering.size(); ++i) { + this->disableInterference(nonInterfering[i]); + } +} + +//////////////////////////////////////////////////////////////////////////////// +/// initialise inputs required for the morphing function void RooLagrangianMorphFunc::init() { @@ -2104,9 +2060,6 @@ int RooLagrangianMorphFunc::countSamples(int nprod, int ndec, int nboth) //////////////////////////////////////////////////////////////////////////////// /// calculate the number of samples needed to morph a certain physics process -/// usage: -/// countSamples ( { RooLagrangianMorphFunc::makeHCggfCouplings(), -/// RooLagrangianMorphFunc::makeHCHZZCouplings() } ) int RooLagrangianMorphFunc::countSamples(std::vector &vertices) { @@ -2122,15 +2075,6 @@ int RooLagrangianMorphFunc::countSamples(std::vector &vertices) return morphfuncpattern.size(); } -//////////////////////////////////////////////////////////////////////////////// -/// create TPair containers of the type expected by the RooLagrangianMorph - -TPair *RooLagrangianMorphFunc::makeCrosssectionContainer(double xs, double unc) -{ - TPair *v = new TPair(new TParameter("xsection", xs), new TParameter("uncertainty", unc)); - return v; -} - //////////////////////////////////////////////////////////////////////////////// /// create only the weight formulas. static function for external usage. @@ -2618,24 +2562,22 @@ RooRealVar *RooLagrangianMorphFunc::getBinWidth() const //////////////////////////////////////////////////////////////////////////////// /// retrieve a histogram output of the current morphing settings -TH1 *RooLagrangianMorphFunc::createTH1(const std::string &name - /*RooFitResult *r*/) +TH1 *RooLagrangianMorphFunc::createTH1(const std::string &name) { - return this->createTH1(name, false /*r*/); + return this->createTH1(name, false); } //////////////////////////////////////////////////////////////////////////////// /// retrieve a histogram output of the current morphing settings -TH1 *RooLagrangianMorphFunc::createTH1(const std::string &name, bool correlateErrors - /* RooFitResult *r*/) +TH1 *RooLagrangianMorphFunc::createTH1(const std::string &name, bool correlateErrors) { auto mf = std::make_unique(*(this->getFunc())); RooRealVar *observable = this->getObservable(); const int nbins = observable->getBins(); - TH1 *hist = new TH1F(name.c_str(), name.c_str(), nbins, observable->getBinning().array()); + auto hist = std::make_unique(name.c_str(), name.c_str(), nbins, observable->getBinning().array()); RooArgSet *args = mf->getComponents(); for (int i = 0; i < nbins; ++i) { @@ -2663,7 +2605,7 @@ TH1 *RooLagrangianMorphFunc::createTH1(const std::string &name, bool correlateEr hist->SetBinContent(i + 1, val); hist->SetBinError(i + 1, correlateErrors ? unc : sqrt(unc2)); } - return hist; + return hist.release(); } //////////////////////////////////////////////////////////////////////////////// From 393d26a9061c91d911dd1c4c86068ae87ab10a7c Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Fri, 15 Jul 2022 13:55:54 +0200 Subject: [PATCH 2/6] [RF] Fixup to parameter index calculation in ParamHistFunc Last year, with commit 3657e7c, the parameter index calculation was changed to be on the fly instead of using a look-up map, which is much faster. However, the implemented formula was not correct for two or three dimensions, which is fixed by this commit. To make sure that the index computation is correct this time, the new code was tested in this code snippet with various inputs: ```C++ void runTest(int nx = 42, int ny = 42, int nz = 42) { const int nxy = nx * ny; const int nyz = ny * nz; for (int i = 0; i < nx; ++i) { for (int j = 0; j < ny; ++j) { for (int k = 0; k < nz; ++k) { const int index = k + j * nz + i * ny * nz; const int gammaIndex = i + j * nx + k * nx * ny; const int i2 = index / nyz; const int tmp = index % nyz; const int j2 = tmp / nz; const int k2 = tmp % nz; const int gammaIndex2 = i2 + j2 * nx + k2 * nxy; if (gammaIndex2 != gammaIndex) { std::cout << "The unraveled indices were not correct!" << std::endl; return; } } } } } ``` Needs to be backported to the 6.26 branch to get into the 6.26.06 patch release. This commit the following problem reported on the forum: https://root-forum.cern.ch/t/cpycppyy-segfault-on-mac-m1/50822 --- roofit/histfactory/src/ParamHistFunc.cxx | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/roofit/histfactory/src/ParamHistFunc.cxx b/roofit/histfactory/src/ParamHistFunc.cxx index 9106e6860cab6..812deb83db55e 100644 --- a/roofit/histfactory/src/ParamHistFunc.cxx +++ b/roofit/histfactory/src/ParamHistFunc.cxx @@ -210,9 +210,13 @@ RooAbsReal& ParamHistFunc::getParameter( Int_t index ) const { _numBinsPerDim = getNumBinsPerDim(_dataVars); } - int i = index / n.yz ; - int j = (index % n.y) / n.z; - int k = index % (n.yz); + // Unravel the index to 3D coordinates. We can't use the index directly, + // because in the parameter set the dimensions are ordered in reverse order + // compared to the RooDataHist (for historical reasons). + const int i = index / n.yz; + const int tmp = index % n.yz; + const int j = tmp / n.z; + const int k = tmp % n.z; return static_cast(_paramSet[i + j * n.x + k * n.xy]); } From b5237a2f2671e235bf07993e57b7ac448269c765 Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Thu, 21 Jul 2022 19:26:52 +0200 Subject: [PATCH 3/6] [RF] Skip out-of-range events in RDataFrame to RooDataSet Helper The RDataFrameHelper should be consistent with creating a RooDataSet from a TTree, meaning out-of-range events should be skipped. This is implemented in this commit, borrowing the logic from `RooTreeDataStore::loadValues()`. A unit test is also implemented. The previous logic of just taking just all values to fill the dataset was very dangerous, because these values then clipped to the RooRealVar limits and biased the number of events observed at the boundaries. Closes #11017. --- .../RDataFrameHelpers/inc/RooAbsDataHelper.h | 38 +++++++++++++- .../test/testActionHelpers.cxx | 51 +++++++++++++++++++ 2 files changed, 87 insertions(+), 2 deletions(-) diff --git a/roofit/RDataFrameHelpers/inc/RooAbsDataHelper.h b/roofit/RDataFrameHelpers/inc/RooAbsDataHelper.h index 78288c651bf94..b16374f57be51 100644 --- a/roofit/RDataFrameHelpers/inc/RooAbsDataHelper.h +++ b/roofit/RDataFrameHelpers/inc/RooAbsDataHelper.h @@ -22,6 +22,7 @@ #include #include #include +#include #include #include @@ -71,6 +72,7 @@ class RooAbsDataHelper : public ROOT::Detail::RDF::RActionImpl _dataset; std::mutex _mutex_dataset; + std::size_t _numInvalid = 0; std::vector> _events; // One vector of values per data-processing slot const std::size_t _eventSize; // Number of variables in dataset @@ -146,6 +148,11 @@ class RooAbsDataHelper : public ROOT::Detail::RDF::RActionImpl0) { + const auto prefix = std::string(_dataset->ClassName()) + "Helper::Finalize(" + _dataset->GetName() + ") "; + oocoutW(static_cast(nullptr), DataHandling) << prefix << "Ignored " << _numInvalid << " out-of-range events\n"; + } } @@ -163,10 +170,37 @@ class RooAbsDataHelper : public ROOT::Detail::RDF::RActionImplget(); for (std::size_t i = 0; i < events.size(); i += eventSize) { + + // Creating a RooDataSet from an RDataFrame should be consistent with the + // creation from a TTree. The construction from a TTree discards entries + // outside the variable definition range, so we have to do that too (see + // also RooTreeDataStore::loadValues). + + bool allOK = true; for (std::size_t j=0; j < eventSize; ++j) { - static_cast(argSet[j])->setVal(events[i+j]); + auto * destArg = static_cast(argSet[j]); + double sourceVal = events[i+j]; + + if (!destArg->inRange(sourceVal, nullptr)) { + _numInvalid++ ; + allOK = false; + const auto prefix = std::string(_dataset->ClassName()) + "Helper::FillDataSet(" + _dataset->GetName() + ") "; + if (_numInvalid < 5) { + // Unlike in the TreeVectorStore case, we don't log the event + // number here because we don't know it anyway, because of + // RDataFrame slots and multithreading. + oocoutI(static_cast(nullptr), DataHandling) << prefix << "Skipping event because " << destArg->GetName() + << " cannot accommodate the value " << sourceVal << "\n"; + } else if (_numInvalid == 5) { + oocoutI(static_cast(nullptr), DataHandling) << prefix << "Skipping ...\n"; + } + break ; + } + destArg->setVal(sourceVal); + } + if(allOK) { + _dataset->add(argSet); } - _dataset->add(argSet); } } }; diff --git a/roofit/RDataFrameHelpers/test/testActionHelpers.cxx b/roofit/RDataFrameHelpers/test/testActionHelpers.cxx index fc4db76b78e5d..62153b5de583d 100644 --- a/roofit/RDataFrameHelpers/test/testActionHelpers.cxx +++ b/roofit/RDataFrameHelpers/test/testActionHelpers.cxx @@ -4,8 +4,12 @@ /// \author Stephan Hageboeck (CERN) #include +#include +#include #include +#include +#include #include #include "gtest/gtest.h" @@ -76,3 +80,50 @@ TEST(RooAbsDataHelper, MTConstruction) EXPECT_NEAR(rooDataHist->moment(y, 2.), 0.25, 1.E-2); // Variance is affected in a binned distribution } + +/// This test verifies that out-of-range events are correctly skipped, +/// consistent with the construction of a RooDataSet from a TTree. +TEST(RooAbsDataHelper, SkipEventsOutOfRange) { + + RooMsgService::instance().getStream(0).removeTopic(RooFit::DataHandling); + RooMsgService::instance().getStream(1).removeTopic(RooFit::DataHandling); + + std::size_t nEvents = 100; + const char * filename = "testRooAbsDataHelper_SkipEventsOutOfRange_tree.root"; + const char * treename = "tree"; + + { + // Create the ROOT file with the dataset + ROOT::RDataFrame rdf(nEvents); + auto rdf_x = rdf.Define("x", [](){ return gRandom->Gaus(0.0, 1.0); }); + rdf_x.Snapshot(treename, filename); + // We can't reuse the same RDataFrame now, because when we rerun the event + // loop it would generate new random values. So this scope ends here and we + // open a new RDF from the file later. + } + + // Open dataset with RDataFrame and TTree + std::unique_ptr file{TFile::Open(filename, "READ")}; + auto tree = file->Get(treename); + ROOT::RDataFrame rdf(treename, filename); + + // Create RooFit variable + RooRealVar x{"x", "x", 0.0, -2.0, 2.0}; + + // Create a RooDataset from the TTree, and one from the RDataFrame + RooDataSet dataSetTree{"dataSetTree", "dataSetTree", tree, x}; + auto dataSetRDF = rdf.Book(RooDataSetHelper("dataSetRDF", "dataSetRDF", RooArgSet(x)), {"x"}); + + // Check if in the creation of the datasets, the entries outside the + // variable range were sucessfully discarded. + double nPassing = *rdf.Filter("x >= -2 && x <= 2.0").Count(); + + EXPECT_EQ(dataSetRDF->numEntries(), nPassing); + EXPECT_EQ(dataSetTree.numEntries(), nPassing); + + file.reset(); // Close file + gSystem->Unlink(filename); // delete temporary file + + RooMsgService::instance().getStream(0).addTopic(RooFit::DataHandling); + RooMsgService::instance().getStream(1).addTopic(RooFit::DataHandling); +} From c48e7c483e4c87fd1ccd44a7aec6c5b472b06a7d Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Fri, 22 Jul 2022 13:34:40 +0200 Subject: [PATCH 4/6] [RF] Remove unnecessary data members from RooLagrangianMorphFunc This commit removes the `_curNormSet` and `_ownParameters` member variables from the RooLagrangianMorphFunc. The `_ownParameters` was simply unused, and the `_curNormSet` was redundant because there is already `RooAbsReal::_lastNSet` that also points to the current normalization set. Now, the `RooLagrangianMorphFunc::getValV()` can also be removed because the sole purpose was to set `_curNormSet`, and also the private `RooLagrangianMorphFunc::getCache()` interface can be changed to take no normSet parameter (as it was unused anyway). The `RooLagrangianMorphFunc::_cacheMgr` declaration is also moved to the bottom of the file together with the other member variables. --- roofit/roofit/inc/RooLagrangianMorphFunc.h | 12 ++-- roofit/roofit/src/RooLagrangianMorphFunc.cxx | 61 ++++++++------------ 2 files changed, 28 insertions(+), 45 deletions(-) diff --git a/roofit/roofit/inc/RooLagrangianMorphFunc.h b/roofit/roofit/inc/RooLagrangianMorphFunc.h index 3ce2ca4f11943..3548ead4ba7fe 100644 --- a/roofit/roofit/inc/RooLagrangianMorphFunc.h +++ b/roofit/roofit/inc/RooLagrangianMorphFunc.h @@ -113,7 +113,6 @@ class RooLagrangianMorphFunc : public RooAbsReal { bool isBinnedDistribution(const RooArgSet &obs) const override; double evaluate() const override; TObject *clone(const char *newname) const override; - double getValV(const RooArgSet *set = 0) const override; bool checkObservables(const RooArgSet *nset) const override; bool forceAnalyticalInt(const RooAbsArg &arg) const override; @@ -179,16 +178,13 @@ class RooLagrangianMorphFunc : public RooAbsReal { class CacheElem; void init(); void setup(bool ownParams = true); - bool _ownParameters = false; void disableInterference(const std::vector &nonInterfering); void disableInterferences(const std::vector> &nonInterfering); - mutable RooObjCacheManager _cacheMgr; //! The cache manager - void addFolders(const RooArgList &folders); bool hasCache() const; - RooLagrangianMorphFunc::CacheElem *getCache(const RooArgSet *nset) const; + RooLagrangianMorphFunc::CacheElem *getCache() const; void updateSampleWeights(); RooRealVar *setupObservable(const char *obsname, TClass *mode, TObject *inputExample); @@ -251,7 +247,9 @@ class RooLagrangianMorphFunc : public RooAbsReal { static std::unique_ptr makeRatio(const char *name, const char *title, RooArgList &nr, RooArgList &dr); private: - double _scale = 1; + + mutable RooObjCacheManager _cacheMgr; //! The cache manager + double _scale = 1.0; std::map _sampleMap; RooListProxy _physics; RooSetProxy _operators; @@ -260,8 +258,6 @@ class RooLagrangianMorphFunc : public RooAbsReal { RooListProxy _flags; Config _config; std::vector> _diagrams; - mutable const RooArgSet *_curNormSet = nullptr; //! - std::vector _nonInterfering; ClassDefOverride(RooLagrangianMorphFunc, 1) diff --git a/roofit/roofit/src/RooLagrangianMorphFunc.cxx b/roofit/roofit/src/RooLagrangianMorphFunc.cxx index 3530cddb9f38c..df1b320977cd1 100644 --- a/roofit/roofit/src/RooLagrangianMorphFunc.cxx +++ b/roofit/roofit/src/RooLagrangianMorphFunc.cxx @@ -1669,7 +1669,7 @@ inline void RooLagrangianMorphFunc::updateSampleWeights() { //#ifdef USE_MULTIPRECISION_LC int sampleidx = 0; - auto cache = this->getCache(_curNormSet); + auto cache = this->getCache(); const size_t n(size(cache->_inverse)); for (auto sampleit : this->_config.paramCards) { const std::string sample(sampleit.first); @@ -1824,7 +1824,7 @@ RooLagrangianMorphFunc::RooLagrangianMorphFunc(const char *name, const char *tit : RooAbsReal(name, title), _cacheMgr(this, 10, kTRUE, kTRUE), _operators("operators", "set of operators", this, kTRUE, kFALSE), _observables("observables", "morphing observables", this, kTRUE, kFALSE), - _binWidths("binWidths", "set of binWidth objects", this, kTRUE, kFALSE), _config(config), _curNormSet(0) + _binWidths("binWidths", "set of binWidth objects", this, kTRUE, kFALSE), _config(config) { this->init(); this->disableInterferences(_config.nonInterfering); @@ -1841,7 +1841,7 @@ RooLagrangianMorphFunc::RooLagrangianMorphFunc(const char *name, const char *tit : RooAbsReal(name, title), _cacheMgr(this, 10, kTRUE, kTRUE), _operators("operators", "set of operators", this, kTRUE, kFALSE), _observables("observables", "morphing observables", this, kTRUE, kFALSE), - _binWidths("binWidths", "set of binWidth objects", this, kTRUE, kFALSE), _curNormSet(0) + _binWidths("binWidths", "set of binWidth objects", this, kTRUE, kFALSE) { this->_config.fileName = filename; this->_config.observableName = observableName; @@ -1859,8 +1859,6 @@ RooLagrangianMorphFunc::RooLagrangianMorphFunc(const char *name, const char *tit void RooLagrangianMorphFunc::setup(bool own) { - _ownParameters = own; - if (_config.couplings.size() > 0) { RooArgList operators; std::vector vertices; @@ -1975,7 +1973,7 @@ RooLagrangianMorphFunc::RooLagrangianMorphFunc(const RooLagrangianMorphFunc &oth _operators(other._operators.GetName(), this, other._operators), _observables(other._observables.GetName(), this, other._observables), _binWidths(other._binWidths.GetName(), this, other._binWidths), _flags(other._flags.GetName(), this, other._flags), - _config(other._config), _curNormSet(0) + _config(other._config) { for (size_t j = 0; j < other._diagrams.size(); ++j) { std::vector diagram; @@ -2209,7 +2207,7 @@ std::vector RooLagrangianMorphFunc::getSamples() const RooAbsReal *RooLagrangianMorphFunc::getSampleWeight(const char *name) { - auto cache = this->getCache(_curNormSet); + auto cache = this->getCache(); auto wname = std::string("w_") + name + "_" + this->GetName(); return dynamic_cast(cache->_weights.find(wname.c_str())); } @@ -2227,8 +2225,8 @@ void RooLagrangianMorphFunc::printWeights() const void RooLagrangianMorphFunc::printSampleWeights() const { - auto *cache = this->getCache(this->_curNormSet); - for (const auto &sample : this->_sampleMap) { + auto *cache = this->getCache(); + for (const auto &sample : _sampleMap) { auto weightName = std::string("w_") + sample.first + "_" + this->GetName(); auto weight = static_cast(cache->_weights.find(weightName.c_str())); if (!weight) @@ -2261,7 +2259,7 @@ void RooLagrangianMorphFunc::randomizeParameters(double z) bool RooLagrangianMorphFunc::updateCoefficients() { - auto cache = this->getCache(_curNormSet); + auto cache = this->getCache(); std::string filename = this->_config.fileName; TDirectory *file = openFile(filename.c_str()); @@ -2329,7 +2327,7 @@ bool RooLagrangianMorphFunc::useCoefficients(const char *filename) cache = RooLagrangianMorphFunc::CacheElem::createCache(this, readMatrixFromFileT(filename)); if (!cache) coutE(Caching) << "unable to create cache!" << std::endl; - this->_cacheMgr.setObj(0, 0, cache, 0); + _cacheMgr.setObj(nullptr, nullptr, cache, nullptr); return true; } @@ -2338,7 +2336,7 @@ bool RooLagrangianMorphFunc::useCoefficients(const char *filename) bool RooLagrangianMorphFunc::writeCoefficients(const char *filename) { - auto cache = this->getCache(_curNormSet); + auto cache = this->getCache(); if (!cache) return false; writeMatrixToFileT(cache->_inverse, filename); @@ -2348,7 +2346,7 @@ bool RooLagrangianMorphFunc::writeCoefficients(const char *filename) //////////////////////////////////////////////////////////////////////////////// /// retrieve the cache object -typename RooLagrangianMorphFunc::CacheElem *RooLagrangianMorphFunc::getCache(const RooArgSet * /*nset*/) const +typename RooLagrangianMorphFunc::CacheElem *RooLagrangianMorphFunc::getCache() const { auto cache = static_cast(_cacheMgr.getObj(0, (RooArgSet *)0)); if (!cache) { @@ -2356,7 +2354,7 @@ typename RooLagrangianMorphFunc::CacheElem *RooLagrangianMorphFunc::getCache(con cxcoutP(Caching) << "current storage has size " << this->_sampleMap.size() << std::endl; cache = RooLagrangianMorphFunc::CacheElem::createCache(this); if (cache) - this->_cacheMgr.setObj(0, 0, cache, 0); + _cacheMgr.setObj(nullptr, nullptr, cache, nullptr); else coutE(Caching) << "unable to create cache!" << std::endl; } @@ -2368,7 +2366,7 @@ typename RooLagrangianMorphFunc::CacheElem *RooLagrangianMorphFunc::getCache(con bool RooLagrangianMorphFunc::hasCache() const { - return (bool)(_cacheMgr.getObj(0, (RooArgSet *)0)); + return (bool)(_cacheMgr.getObj(nullptr, static_cast(nullptr))); } //////////////////////////////////////////////////////////////////////////////// @@ -2487,7 +2485,7 @@ double RooLagrangianMorphFunc::getParameterValue(const char *name) const if (param) { return param->getVal(); } - return 0; + return 0.0; } //////////////////////////////////////////////////////////////////////////////// @@ -2690,7 +2688,7 @@ int RooLagrangianMorphFunc::nParameters() const int RooLagrangianMorphFunc::nPolynomials() const { // return the number of samples in this morphing function - auto cache = getCache(_curNormSet); + auto cache = getCache(); return cache->_formulas.size(); } @@ -2735,7 +2733,7 @@ const RooArgSet *RooLagrangianMorphFunc::getParameterSet() const const RooArgList *RooLagrangianMorphFunc::getCouplingSet() const { - auto cache = getCache(_curNormSet); + auto cache = getCache(); return &(cache->_couplings); } @@ -2777,7 +2775,7 @@ void RooLagrangianMorphFunc::setParameters(const ParamSet ¶ms) std::unique_ptr RooLagrangianMorphFunc::createPdf() const { - auto cache = getCache(_curNormSet); + auto cache = getCache(); auto func = std::make_unique(*(cache->_sumFunc)); // create a wrapper on the roorealsumfunc @@ -2789,7 +2787,7 @@ std::unique_ptr RooLagrangianMorphFunc::createPdf() const RooRealSumFunc *RooLagrangianMorphFunc::getFunc() const { - auto cache = getCache(_curNormSet); + auto cache = getCache(); return cache->_sumFunc.get(); } @@ -2835,7 +2833,7 @@ Double_t RooLagrangianMorphFunc::expectedEvents(const RooArgSet &nset) const double RooLagrangianMorphFunc::expectedUncertainty() const { RooRealVar *observable = this->getObservable(); - auto cache = this->getCache(_curNormSet); + auto cache = this->getCache(); double unc2 = 0; for (const auto &sample : this->_sampleMap) { RooAbsArg *phys = this->_physics.at(sample.second); @@ -2923,23 +2921,12 @@ std::list *RooLagrangianMorphFunc::plotSamplingHint(RooAbsRealLValue & //////////////////////////////////////////////////////////////////////////////// /// call getVal on the internal function -Double_t RooLagrangianMorphFunc::getValV(const RooArgSet *set) const -{ - // cout << "XX RooLagrangianMorphFunc::getValV(" << this << ") set = " << set - // << std::endl ; - this->_curNormSet = set; - return RooAbsReal::getValV(set); -} - -//////////////////////////////////////////////////////////////////////////////// -/// call getVal on the internal function - -Double_t RooLagrangianMorphFunc::evaluate() const +double RooLagrangianMorphFunc::evaluate() const { // call getVal on the internal function RooRealSumFunc *pdf = this->getFunc(); if (pdf) - return this->_scale * pdf->getVal(_curNormSet); + return _scale * pdf->getVal(_lastNSet); else std::cerr << "unable to acquire in-built function!" << std::endl; return 0.; @@ -3015,7 +3002,7 @@ void RooLagrangianMorphFunc::setCacheAndTrackHints(RooArgSet &arg) TMatrixD RooLagrangianMorphFunc::getMatrix() const { - auto cache = getCache(_curNormSet); + auto cache = getCache(); if (!cache) coutE(Caching) << "unable to retrieve cache!" << std::endl; return makeRootMatrix(cache->_matrix); @@ -3026,7 +3013,7 @@ TMatrixD RooLagrangianMorphFunc::getMatrix() const TMatrixD RooLagrangianMorphFunc::getInvertedMatrix() const { - auto cache = getCache(_curNormSet); + auto cache = getCache(); if (!cache) coutE(Caching) << "unable to retrieve cache!" << std::endl; return makeRootMatrix(cache->_inverse); @@ -3039,7 +3026,7 @@ TMatrixD RooLagrangianMorphFunc::getInvertedMatrix() const double RooLagrangianMorphFunc::getCondition() const { - auto cache = getCache(_curNormSet); + auto cache = getCache(); if (!cache) coutE(Caching) << "unable to retrieve cache!" << std::endl; return cache->_condition; From 81d3934b7d43a3a3d57c35bb3b214d2f70ac7b88 Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Fri, 22 Jul 2022 13:38:25 +0200 Subject: [PATCH 5/6] [RF] Create dictionaries for RooLagrangianMorphFunc::Config This is to enable IO for the RooLagrangianMorphFunc. --- roofit/roofit/inc/LinkDef1.h | 1 + 1 file changed, 1 insertion(+) diff --git a/roofit/roofit/inc/LinkDef1.h b/roofit/roofit/inc/LinkDef1.h index bc59baba0df95..f91889d7e3b36 100644 --- a/roofit/roofit/inc/LinkDef1.h +++ b/roofit/roofit/inc/LinkDef1.h @@ -61,6 +61,7 @@ #pragma link C++ class RooJeffreysPrior+ ; #pragma link C++ class RooJohnson+; #pragma link C++ class RooLagrangianMorphFunc+; +#pragma link C++ class RooLagrangianMorphFunc::Config+; #pragma link C++ class RooFunctorBinding+ ; #pragma link C++ class RooFunctor1DBinding+ ; #pragma link C++ class RooFunctorPdfBinding+ ; From e7ac26b6d0012fcd8fc866511171c70553d9e4c4 Mon Sep 17 00:00:00 2001 From: Jonas Rembser Date: Fri, 22 Jul 2022 13:45:36 +0200 Subject: [PATCH 6/6] [RF] Correctly set plot titles in the RooLagrangianMorphFunc tutorials The title for the x-axis of the first plot should be `"p_{T}^{V}"` and not `"c_{Hq^{(3)}}"`. Also, the title for the observable is changed to `"p_{T}^{V}"` to get the right axis label for the x-axes of the plots. --- tutorials/roofit/rf711_lagrangianmorph.C | 5 +---- tutorials/roofit/rf711_lagrangianmorph.py | 5 +---- 2 files changed, 2 insertions(+), 8 deletions(-) diff --git a/tutorials/roofit/rf711_lagrangianmorph.C b/tutorials/roofit/rf711_lagrangianmorph.C index a78d915ed0eba..5247e2ffc74e0 100644 --- a/tutorials/roofit/rf711_lagrangianmorph.C +++ b/tutorials/roofit/rf711_lagrangianmorph.C @@ -38,7 +38,7 @@ void rf711_lagrangianmorph() std::string observablename = "pTV"; // Setup observable that is morphed - RooRealVar obsvar(observablename.c_str(), "observable of pTV", 10, 600); + RooRealVar obsvar(observablename.c_str(), "p_{T}^{V}", 10, 600); // Setup two couplings that enters the morphing function // kSM -> SM coupling set to constant (1) @@ -134,7 +134,6 @@ void rf711_lagrangianmorph() gPad->SetRightMargin(0.05); frame0->Draw(); - frame0->GetXaxis()->SetTitle("c_{Hq^{(3)}}"); TLegend *leg1 = new TLegend(0.55, 0.65, 0.94, 0.87); leg1->SetTextSize(0.04); leg1->SetFillColor(kWhite); @@ -151,7 +150,6 @@ void rf711_lagrangianmorph() gPad->SetRightMargin(0.05); frame1->Draw(); - frame1->GetXaxis()->SetTitle("p_{T}^{V}"); TLegend *leg2 = new TLegend(0.60, 0.65, 0.94, 0.87); leg2->SetTextSize(0.04); @@ -175,7 +173,6 @@ void rf711_lagrangianmorph() gPad->SetLogz(); hh_data->GetYaxis()->SetTitle("c_{Hq^{(3)}}"); hh_data->GetYaxis()->SetRangeUser(0, 0.5); - hh_data->GetXaxis()->SetTitle("p_{T}^{V}"); hh_data->GetZaxis()->SetTitleOffset(1.8); hh_data->Draw("COLZ"); c1->SaveAs("rf711_lagrangianmorph.png"); diff --git a/tutorials/roofit/rf711_lagrangianmorph.py b/tutorials/roofit/rf711_lagrangianmorph.py index 8183a8e5a3bf5..2320ac3dabf1b 100644 --- a/tutorials/roofit/rf711_lagrangianmorph.py +++ b/tutorials/roofit/rf711_lagrangianmorph.py @@ -23,7 +23,7 @@ observablename = "pTV" # Setup observable that is to be morphed -obsvar = ROOT.RooRealVar(observablename, "observable of pTV", 10, 600) +obsvar = ROOT.RooRealVar(observablename, "p_{T}^{V}", 10, 600) # Setup two couplings that enters the morphing function # kSM -> SM coupling set to constant (1) @@ -113,7 +113,6 @@ ROOT.gPad.SetRightMargin(0.05) frame0.Draw() -frame0.GetXaxis().SetTitle("c_{Hq^{(3)}}") leg1 = ROOT.TLegend(0.55, 0.65, 0.94, 0.87) leg1.SetTextSize(0.04) leg1.SetFillColor(ROOT.kWhite) @@ -130,7 +129,6 @@ ROOT.gPad.SetRightMargin(0.05) frame1.Draw() -frame1.GetXaxis().SetTitle("p_{T}^{V}") leg2 = ROOT.TLegend(0.62, 0.65, 0.94, 0.87) leg2.SetTextSize(0.04) @@ -155,7 +153,6 @@ ROOT.gPad.SetLogz() hh_data.GetYaxis().SetTitle("c_{Hq^{(3)}}") hh_data.GetYaxis().SetRangeUser(0, 0.5) -hh_data.GetXaxis().SetTitle("p_{T}^{V}") hh_data.GetZaxis().SetTitleOffset(1.8) hh_data.Draw("COLZ") c1.SaveAs("rf711_lagrangianmorph.png")