[RF] Move updateCoefficients() from RooAddPdf/Model to RooAddHelpers

The same function was duplicated in RooAddPdf and RooAddModel, and it
can therefore be moved into the new `RooAddHelpers` file.

Author: guitargeek

roofit/roofitcore/inc/RooAbsReal.h

@@ -545,6 +545,7 @@ class RooAbsReal : public RooAbsArg {
 
   friend class RooRealSumPdf ;
   friend class RooRealSumFunc;
+  friend class RooAddHelpers ;
   friend class RooAddPdf ;
   friend class RooAddModel ;
   void selectComp(bool flag) {

roofit/roofitcore/inc/RooAddModel.h

@@ -33,7 +33,6 @@ class RooAddModel : public RooResolutionModel {
   RooAddModel(const RooAddModel& other, const char* name=nullptr) ;
   TObject* clone(const char* newname) const override { return new RooAddModel(*this,newname) ; }
   RooResolutionModel* convolution(RooFormulaVar* basis, RooAbsArg* owner) const override ;
-  ~RooAddModel() override ;
 
   double evaluate() const override ;
   bool checkObservables(const RooArgSet* nset) const override ;
@@ -93,10 +92,10 @@ class RooAddModel : public RooResolutionModel {
   void selectNormalizationRange(const char* rangeName=nullptr, bool force=false) override ;
 
   mutable RooSetProxy _refCoefNorm ;   ///<! Reference observable set for coefficient interpretation
-  mutable TNamed* _refCoefRangeName ;  ///<! Reference range name for coefficient interpretation
+  mutable TNamed* _refCoefRangeName = nullptr;  ///<! Reference range name for coefficient interpretation
 
-  bool _projectCoefs ;         ///< If true coefficients need to be projected for use in evaluate()
-  mutable double* _coefCache ; ///<! Transiet cache with transformed values of coefficients
+  bool _projectCoefs = false;  ///< If true coefficients need to be projected for use in evaluate()
+  mutable std::vector<double> _coefCache; ///<! Transiet cache with transformed values of coefficients
 
 
   mutable RooObjCacheManager _projCacheMgr ;  ///<! Manager of cache with coefficient projections and transformations
@@ -114,14 +113,14 @@ class RooAddModel : public RooResolutionModel {
 
   mutable RooObjCacheManager _intCacheMgr ; ///<! Manager of cache with integrals
 
-  mutable RooAICRegistry _codeReg ;  ///<! Registry of component analytical integration codes
+  mutable RooAICRegistry _codeReg = 10; ///<! Registry of component analytical integration codes
 
   RooListProxy _pdfList ;   ///<  List of component PDFs
   RooListProxy _coefList ;  ///<  List of coefficients
   mutable RooArgList* _snormList{nullptr};  ///<!  List of supplemental normalization factors
 
-  bool _haveLastCoef ;    ///<  Flag indicating if last PDFs coefficient was supplied in the ctor
-  bool _allExtendable ;   ///<  Flag indicating if all PDF components are extendable
+  bool _haveLastCoef = false;    ///<  Flag indicating if last PDFs coefficient was supplied in the ctor
+  bool _allExtendable = false;   ///<  Flag indicating if all PDF components are extendable
 
   mutable Int_t _coefErrCount ; ///<! Coefficient error counter
 

roofit/roofitcore/inc/RooAddPdf.h

@@ -101,7 +101,7 @@ class RooAddPdf : public RooAbsPdf {
   mutable TNamed* _refCoefRangeName = nullptr ;  ///< Reference range name for coefficient interpreation
 
   bool _projectCoefs = false;     ///< If true coefficients need to be projected for use in evaluate()
-  std::vector<double> _coefCache; ///<! Transient cache with transformed values of coefficients
+  mutable std::vector<double> _coefCache; ///<! Transient cache with transformed values of coefficients
 
 
   mutable RooObjCacheManager _projCacheMgr ;  //! Manager of cache with coefficient projections and transformations

roofit/roofitcore/src/RooAddHelpers.cxx

@@ -12,6 +12,7 @@
 
 #include <RooAbsPdf.h>
 #include <RooArgSet.h>
+#include <RooNaNPacker.h>
 #include <RooRealConstant.h>
 #include <RooRealIntegral.h>
 #include <RooRealVar.h>
@@ -239,3 +240,127 @@ RooArgList AddCacheElem::containedArgs(Action)
 
    return allNodes;
 }
+
+void RooAddHelpers::updateCoefficients(RooAbsPdf const &addPdf, RooArgList const &pdfList, RooArgList const &coefList,
+                                       AddCacheElem &cache, const RooArgSet *nset, bool projectCoefs,
+                                       RooArgSet const &refCoefNorm, bool allExtendable, std::vector<double> &coefCache,
+                                       int &coefErrCount)
+{
+   bool haveLastCoef = pdfList.size() == coefList.size();
+
+   coefCache.resize(haveLastCoef ? coefList.size() : pdfList.size(), 0.);
+
+   // Straight coefficients
+   if (allExtendable) {
+
+      // coef[i] = expectedEvents[i] / SUM(expectedEvents)
+      double coefSum(0);
+      std::size_t i = 0;
+      for (auto arg : pdfList) {
+         auto pdf = static_cast<RooAbsPdf *>(arg);
+         coefCache[i] = pdf->expectedEvents(!refCoefNorm.empty() ? &refCoefNorm : nset);
+         coefSum += coefCache[i];
+         i++;
+      }
+
+      if (coefSum == 0.) {
+         oocoutW(&addPdf, Eval) << addPdf.ClassName() << "::updateCoefCache(" << addPdf.GetName()
+                                << ") WARNING: total number of expected events is 0" << std::endl;
+      } else {
+         for (std::size_t j = 0; j < pdfList.size(); j++) {
+            coefCache[j] /= coefSum;
+         }
+      }
+
+   } else {
+      if (haveLastCoef) {
+
+         // coef[i] = coef[i] / SUM(coef)
+         double coefSum(0);
+         std::size_t i = 0;
+         for (auto coefArg : coefList) {
+            auto coef = static_cast<RooAbsReal *>(coefArg);
+            coefCache[i] = coef->getVal(nset);
+            coefSum += coefCache[i++];
+         }
+         if (coefSum == 0.) {
+            oocoutW(&addPdf, Eval) << addPdf.ClassName() << "::updateCoefCache(" << addPdf.GetName()
+                                   << ") WARNING: sum of coefficients is zero 0" << std::endl;
+         } else {
+            const double invCoefSum = 1. / coefSum;
+            for (std::size_t j = 0; j < coefList.size(); j++) {
+               coefCache[j] *= invCoefSum;
+            }
+         }
+      } else {
+
+         // coef[i] = coef[i] ; coef[n] = 1-SUM(coef[0...n-1])
+         double lastCoef(1);
+         std::size_t i = 0;
+         for (auto coefArg : coefList) {
+            auto coef = static_cast<RooAbsReal *>(coefArg);
+            coefCache[i] = coef->getVal(nset);
+            lastCoef -= coefCache[i++];
+         }
+         coefCache[coefList.size()] = lastCoef;
+
+         // Treat coefficient degeneration
+         const float coefDegen = lastCoef < 0. ? -lastCoef : (lastCoef > 1. ? lastCoef - 1. : 0.);
+         if (coefDegen > 1.E-5) {
+            coefCache[coefList.size()] = RooNaNPacker::packFloatIntoNaN(100.f * coefDegen);
+
+            std::stringstream msg;
+            if (coefErrCount-- > 0) {
+               msg << "RooAddPdf::updateCoefCache(" << addPdf.GetName()
+                   << " WARNING: sum of PDF coefficients not in range [0-1], value=" << 1 - lastCoef;
+               if (coefErrCount == 0) {
+                  msg << " (no more will be printed)";
+               }
+               oocoutW(&addPdf, Eval) << msg.str() << std::endl;
+            }
+         }
+      }
+   }
+
+   // Stop here if not projection is required or needed
+   if ((!projectCoefs && !addPdf.normRange()) || cache._projList.empty()) {
+      return;
+   }
+
+   // Adjust coefficients for given projection
+   double coefSum(0);
+   {
+      RooAbsReal::GlobalSelectComponentRAII compRAII(true);
+
+      for (std::size_t i = 0; i < pdfList.size(); i++) {
+
+         RooAbsReal *pp = ((RooAbsReal *)cache._projList.at(i));
+         RooAbsReal *sn = ((RooAbsReal *)cache._suppProjList.at(i));
+         RooAbsReal *r1 = ((RooAbsReal *)cache._refRangeProjList.at(i));
+         RooAbsReal *r2 = ((RooAbsReal *)cache._rangeProjList.at(i));
+
+         double proj = pp->getVal() / sn->getVal() * (r2->getVal() / r1->getVal());
+
+         coefCache[i] *= proj;
+         coefSum += coefCache[i];
+      }
+   }
+
+   if ((RooMsgService::_debugCount > 0) &&
+       RooMsgService::instance().isActive(&addPdf, RooFit::Caching, RooFit::DEBUG)) {
+      for (std::size_t i = 0; i < pdfList.size(); ++i) {
+         ooccoutD(&addPdf, Caching) << " ALEX:   POST-SYNC coef[" << i << "] = " << coefCache[i]
+                                    << " ( _coefCache[i]/coefSum = " << coefCache[i] * coefSum << "/" << coefSum
+                                    << " ) " << std::endl;
+      }
+   }
+
+   if (coefSum == 0.) {
+      oocoutE(&addPdf, Eval) << addPdf.ClassName() << "::updateCoefCache(" << addPdf.GetName()
+                             << ") sum of coefficients is zero." << std::endl;
+   }
+
+   for (std::size_t i = 0; i < pdfList.size(); i++) {
+      coefCache[i] /= coefSum;
+   }
+}

roofit/roofitcore/src/RooAddHelpers.h

@@ -34,4 +34,12 @@ class AddCacheElem : public RooAbsCacheElement {
    RooArgList _rangeProjList;    ///< Range integrals to be multiplied with coefficients (target range)
 };
 
+class RooAddHelpers {
+public:
+   static void updateCoefficients(RooAbsPdf const &addPdf, RooArgList const &pdfList, RooArgList const &coefList,
+                                  AddCacheElem &cache, const RooArgSet *nset, bool projectCoefs,
+                                  RooArgSet const &refCoefNorm, bool allExtendable, std::vector<double> &coefCache,
+                                  int &coefErrCount);
+};
+
 #endif