root-project · guitargeek · Jan 24, 2023 · Jan 18, 2023 · Jan 19, 2023 · Jan 12, 2023
@@ -192,6 +192,13 @@ RNTuple is still experimental and is scheduled to become production grade in 202
 
 ## Math Libraries
 
+### KahanSum updates
+
+The `ROOT::Math::KahanSum` class was slightly modified:
+- The behavior of `operator-=` and `operator+=` on a `KahanSum` were not symmetric, leading to slight bit-wise inaccuracies. In fits, where such operations are done a lot of times (e.g. through the offsetting mechanism in RooFit which subtracts a constant `KahanSum` term after each likelihood evaluation), this can add up to significant numerical divergence. An improved algorithm was implemented, based on an algorithm for combining Kahan sums and carry terms (Tian et al. 2012). (PR #11940)
+- The auto-conversion to type `T` and implicit type `T` constructor in `KahanSum` made it hard to debug `KahanSum`, because it is easy to overlook implicit conversions in code, especially in lines where the type of the return value is `auto`. These auto-conversions were removed. Where necessary, they should be replaced with an explicit construction or explicit conversion to double via `Sum()`. (PR #11941)
+- Binary addition and subtraction operators were added, as well as a unary negation operator. (PR #11940)
+- Comparison operators `operator==` and `operator!=` were added.
 
 ## RooFit Libraries
 

@@ -113,6 +113,9 @@ class TGFileBrowser : public TGMainFrame, public TBrowserImp {
    void        ToggleSort();
    void        Update();
 
+   TGListTree     *GetListTree() const { return fListTree; }
+   TGListTreeItem *GetRootDir() const { return fRootDir; }
+
    ClassDefOverride(TGFileBrowser, 0) // File browser.
 };
 

@@ -146,6 +146,7 @@ class TRootBrowser : public TGMainFrame, public TBrowserImp {
    void              HandleMenu(Int_t id);
    void              RecursiveReparent(TGPopupMenu *popup);
    void              RemoveTab(Int_t pos, Int_t subpos);
+   TBrowserImp      *GetActBrowser() const { return fActBrowser; }
    void              SetActBrowser(TBrowserImp *b) { fActBrowser = b; }
    void              ShowMenu(TGCompositeFrame *menu);
    void              StartEmbedding(Int_t pos = kRight, Int_t subpos = -1) override;

@@ -264,6 +264,15 @@ namespace ROOT {
        /// Note that while Tian et al. add the carry in the first step, we subtract
        /// the carry, in accordance with the Add(Indexed) implementation(s) above.
        /// This is purely an implementation choice that has no impact on performance.
+       ///
+       /// \note Take care when using += (and -=) to add other KahanSums into a zero-initialized
+       ///       KahanSum. The operator behaves correctly in this case, but the result may be slightly
+       ///       off if you expect 0 + x to yield exactly x (where 0 is the zero-initialized KahanSum
+       ///       and x another KahanSum). In particular, x's carry term may get lost. This doesn't
+       ///       just happen with zero-initialized KahanSums; see the SubtractWithABitTooSmallCarry
+       ///       test case in the testKahan unittest for other examples. This behavior is internally
+       ///       consistent: the carry also gets lost if you switch the operands and it also happens with
+       ///       other KahanSum operators.
        template<typename U, unsigned int M>
        KahanSum<T, N>& operator+=(const KahanSum<U, M>& other) {
          U corrected_arg_sum = other.Sum() - (fCarry[0] + other.Carry());
@@ -292,6 +301,16 @@ namespace ROOT {
           return {-this->fSum[0], -this->fCarry[0]};
        }
 
+       template<typename U, unsigned int M>
+       bool operator ==(KahanSum<U, M> const& other) const {
+          return (this->Sum() == other.Sum()) && (this->Carry() == other.Carry());
+       }
+
+       template<typename U, unsigned int M>
+       bool operator !=(KahanSum<U, M> const& other) const {
+          return !(*this == other);
+       }
+
      private:
        T fSum[N];
        T fCarry[N];

@@ -16,6 +16,7 @@ endif()
 ROOT_STANDARD_LIBRARY_PACKAGE(RooFitCore
   HEADERS
     RooFit/Detail/DataMap.h
+    RooFit/Detail/NormalizationHelpers.h
     RooFit/Floats.h
     Roo1DTable.h
     RooAbsAnaConvPdf.h

@@ -24,6 +24,7 @@
 #include "RooAbsCache.h"
 #include "RooNameReg.h"
 #include "RooLinkedListIter.h"
+#include <RooFit/Detail/NormalizationHelpers.h>
 #include <RooStringView.h>
 
 #include <deque>
@@ -266,7 +267,9 @@ class RooAbsArg : public TNamed, public RooPrintable {
 
   void addServer(RooAbsArg& server, bool valueProp=true, bool shapeProp=false, std::size_t refCount = 1);
   void addServerList(RooAbsCollection& serverList, bool valueProp=true, bool shapeProp=false) ;
-  void replaceServer(RooAbsArg& oldServer, RooAbsArg& newServer, bool valueProp, bool shapeProp) ;
+  void
+  R__SUGGEST_ALTERNATIVE("This interface is unsafe! Use RooAbsArg::redirectServers()")
+  replaceServer(RooAbsArg& oldServer, RooAbsArg& newServer, bool valueProp, bool shapeProp) ;
   void changeServer(RooAbsArg& server, bool valueProp, bool shapeProp) ;
   void removeServer(RooAbsArg& server, bool force=false) ;
   RooAbsArg *findNewServer(const RooAbsCollection &newSet, bool nameChange) const;
@@ -578,7 +581,7 @@ class RooAbsArg : public TNamed, public RooPrintable {
 
   virtual void applyWeightSquared(bool flag);
 
-  virtual std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& normSet, RooAbsArg const& server) const;
+  virtual std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const;
 
   virtual bool isCategory() const { return false; }
 

@@ -57,6 +57,7 @@ class RooAbsCachedPdf : public RooAbsPdf {
   Int_t getAnalyticalIntegralWN(RooArgSet& allVars, RooArgSet& analVars, const RooArgSet* normSet, const char* rangeName=nullptr) const override ;
   double analyticalIntegralWN(Int_t code, const RooArgSet* normSet, const char* rangeName=nullptr) const override ;
 
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
 
   class PdfCacheElem : public RooAbsCacheElement {
   public:

@@ -321,6 +321,8 @@ class RooAbsPdf : public RooAbsReal {
   virtual RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=nullptr, const RooArgSet* auxProto=nullptr,
                   bool verbose=false, bool autoBinned=true, const char* binnedTag="") const ;
 
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
+
 private:
 
   RooDataSet *generate(RooAbsGenContext& context, const RooArgSet& whatVars, const RooDataSet* prototype,

@@ -106,12 +106,6 @@ class RooBinSamplingPdf : public RooAbsPdf {
   const RooAbsPdf& pdf() const { return _pdf.arg(); }
   const RooAbsReal& observable() const { return _observable.arg(); }
 
-  std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& /*normSet*/,
-                         RooAbsArg const& /*server*/) const override {
-    // servers are evaluated unnormalized
-    return std::make_unique<RooArgSet>();
-  }
-
 protected:
   double evaluate() const override;
   RooSpan<double> evaluateSpan(RooBatchCompute::RunContext& evalData, const RooArgSet* normSet) const override;

@@ -43,7 +43,7 @@ class RooConstraintSum : public RooAbsReal {
 
   void computeBatch(cudaStream_t*, double* output, size_t size, RooFit::Detail::DataMap const&) const override;
 
-  std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& normSet, RooAbsArg const& server) const override;
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
 
 protected:
 

@@ -0,0 +1,63 @@
+/*
+ * Project: RooFit
+ * Authors:
+ *   Jonas Rembser, CERN 2022
+ *
+ * Copyright (c) 2022, CERN
+ *
+ * Redistribution and use in source and binary forms,
+ * with or without modification, are permitted according to the terms
+ * listed in LICENSE (http://roofit.sourceforge.net/license.txt)
+ */
+
+#ifndef RooFit_Detail_NormalizationHelpers_h
+#define RooFit_Detail_NormalizationHelpers_h
+
+#include <memory>
+#include <unordered_map>
+
+class RooAbsArg;
+class RooArgSet;
+
+class TNamed;
+
+namespace RooFit {
+
+namespace Detail {
+
+class CompileContext {
+public:
+   CompileContext(RooArgSet const &topLevelNormSet);
+
+   ~CompileContext();
+
+   template <class T>
+   T *compile(T &arg, RooAbsArg &owner, RooArgSet const &normSet)
+   {
+      return static_cast<T *>(compileImpl(arg, owner, normSet));
+   }
+
+   void compileServers(RooAbsArg &arg, RooArgSet const &normSet);
+
+private:
+   RooAbsArg *compileImpl(RooAbsArg &arg, RooAbsArg &owner, RooArgSet const &normSet);
+   void add(RooAbsArg &arg);
+   RooAbsArg *find(RooAbsArg &arg) const;
+
+   RooArgSet const &_topLevelNormSet;
+   std::unordered_map<TNamed const *, RooAbsArg *> _clonedArgsSet;
+};
+
+template <class T>
+std::unique_ptr<T> compileForNormSet(T const &arg, RooArgSet const &normSet)
+{
+   RooFit::Detail::CompileContext ctx{normSet};
+   std::unique_ptr<RooAbsArg> head = arg.compileForNormSet(normSet, ctx);
+   return std::unique_ptr<T>{static_cast<T *>(head.release())};
+}
+
+} // namespace Detail
+
+} // namespace RooFit
+
+#endif
@@ -98,10 +98,12 @@ class RooProdPdf : public RooAbsPdf {
 
   void writeCacheToStream(std::ostream& os, RooArgSet const* nset) const;
 
-  std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& normSet, RooAbsArg const& server) const override;
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
 
 private:
 
+  std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& normSet, RooAbsArg const& server) const;
+
   double evaluate() const override ;
 
   std::unique_ptr<RooAbsReal> makeCondPdfRatioCorr(RooAbsReal& term, const RooArgSet& termNset, const RooArgSet& termImpSet, const char* normRange, const char* refRange) const ;

@@ -41,6 +41,7 @@ class RooProjectedPdf : public RooAbsPdf {
 
   void printMetaArgs(std::ostream& os) const override ;
 
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
 
 protected:
 

@@ -81,8 +81,6 @@ class RooRealIntegral : public RooAbsReal {
   void setAllowComponentSelection(bool allow);
   bool getAllowComponentSelection() const;
 
-  std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& normSet, RooAbsArg const& server) const override;
-
 protected:
 
   mutable bool _valid = false;

@@ -58,9 +58,7 @@ class RooRealSumFunc : public RooAbsReal {
    CacheMode canNodeBeCached() const override { return RooAbsArg::NotAdvised; };
    void setCacheAndTrackHints(RooArgSet &) override;
 
-   std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& /*normSet*/, RooAbsArg const& /*server*/) const override {
-     return std::make_unique<RooArgSet>();
-   }
+   std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
 
 protected:
    mutable RooObjCacheManager _normIntMgr; //! The integration cache manager

@@ -68,9 +68,7 @@ class RooRealSumPdf : public RooAbsPdf {
   CacheMode canNodeBeCached() const override { return RooAbsArg::NotAdvised ; } ;
   void setCacheAndTrackHints(RooArgSet&) override ;
 
-  std::unique_ptr<RooArgSet> fillNormSetForServer(RooArgSet const& /*normSet*/, RooAbsArg const& /*server*/) const override {
-     return std::make_unique<RooArgSet>();
-  }
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
 
 protected:
 

@@ -57,6 +57,9 @@ class RooSTLRefCountList {
     ///Add an object or increase refCount if it is already present. Only compares
     ///pointers to check for existing objects
     void Add(T * obj, std::size_t initialCount = 1) {
+      // Nothing to add because `refCount` would be zero.
+      if(initialCount == 0) return;
+
       auto foundItem = findByPointer(obj);
 
       if (foundItem != _storage.end()) {
@@ -191,12 +194,18 @@ class RooSTLRefCountList {
     ///Decrease ref count of given object. Shrink list if ref count reaches 0.
     ///\param obj Decrease ref count of given object. Compare by pointer.
     ///\param force If true, remove irrespective of ref count.
-    void Remove(const T * obj, bool force = false) {
+    ///Returns by how much the `refCount` for the element to be removed was
+    ///decreased (zero if nothing was removed). If `force == false`, it can
+    ///only be zero or one, if `force == true`, it can be the full `refCount`
+    ///for that element.
+    int Remove(const T * obj, bool force = false) {
       auto item = findByPointer(obj);
 
       if (item != _storage.end()) {
         const std::size_t pos = item - _storage.begin();
 
+        const UInt_t origRefCount = _refCount[pos];
+
         if (force || --_refCount[pos] == 0) {
           //gcc4.x doesn't know how to erase at the position of a const_iterator
           //Therefore, erase at begin + pos instead of 'item'
@@ -211,8 +220,31 @@ class RooSTLRefCountList {
             // _storage at the beginning of Remove().
             _orderedStorage.erase(std::find(_orderedStorage.begin(), _orderedStorage.end(), obj));
           }
+          return origRefCount;
         }
+        return 1;
       }
+
+      return 0;
+    }
+
+
+    ///Replace an element with a new value, keeping the same `refCount`. Will
+    ///return the `refCount` for that element if the replacement succeeded,
+    ///otherwise returns zero in case the `oldObj` could not be found in the
+    ///collection.
+    int Replace(const T * oldObj, T * newObj) {
+      auto item = findByPointer(oldObj);
+
+      if (item != _storage.end()) {
+        const std::size_t pos = item - _storage.begin();
+        _storage[pos] = newObj;
+        // The content has changed, so the ordered-by-name storage was invalidated.
+        _orderedStorage.clear();
+        return _refCount[pos];
+      }
+
+      return 0;
     }
 
 

@@ -89,12 +89,22 @@ class RooSimultaneous : public RooAbsPdf {
                                  std::map<std::string, double> const& precisions,
                                  bool useCategoryNames=false);
 
+  RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=nullptr, const RooArgSet* auxProto=nullptr,
+                  bool verbose=false, bool autoBinned=true, const char* binnedTag="") const override ;
+  RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr,
+                                  const RooArgSet* auxProto=nullptr, bool verbose= false) const override ;
+
+  std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;
+
 protected:
 
   void initialize(RooAbsCategoryLValue& inIndexCat, std::map<std::string,RooAbsPdf*> pdfMap) ;
 
   void selectNormalization(const RooArgSet* depSet=nullptr, bool force=false) override ;
   void selectNormalizationRange(const char* rangeName=nullptr, bool force=false) override ;
+
+  RooArgSet const& flattenedCatList() const;
+
   mutable RooSetProxy _plotCoefNormSet ;
   const TNamed* _plotCoefNormRange = nullptr;
 
@@ -109,14 +119,12 @@ class RooSimultaneous : public RooAbsPdf {
 
   friend class RooSimGenContext ;
   friend class RooSimSplitGenContext ;
-  RooAbsGenContext* autoGenContext(const RooArgSet &vars, const RooDataSet* prototype=nullptr, const RooArgSet* auxProto=nullptr,
-                  bool verbose=false, bool autoBinned=true, const char* binnedTag="") const override ;
-  RooAbsGenContext* genContext(const RooArgSet &vars, const RooDataSet *prototype=nullptr,
-                                  const RooArgSet* auxProto=nullptr, bool verbose= false) const override ;
 
   RooCategoryProxy _indexCat ; ///< Index category
   TList    _pdfProxyList ;     ///< List of PDF proxies (named after applicable category state)
   Int_t    _numPdf = 0;        ///< Number of registered PDFs
+private:
+  mutable std::unique_ptr<RooArgSet> _indexCatSet ; ///<! Index category wrapped in a RooArgSet if needed internally
 
   ClassDefOverride(RooSimultaneous,3)  // Simultaneous operator p.d.f, functions like C++  'switch()' on input p.d.fs operating on index category5A
 };
Original file line number	Diff line number	Diff line change
Expand Up		@@ -41,6 +41,7 @@ class RooProjectedPdf : public RooAbsPdf {

		void printMetaArgs(std::ostream& os) const override ;

		std::unique_ptr<RooAbsArg> compileForNormSet(RooArgSet const &normSet, RooFit::Detail::CompileContext & ctx) const override;

		protected:

Expand Down