Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-09-16 09:31:12

0001 // -*- C++ -*-
0002 //
0003 // This file is part of YODA -- Yet more Objects for Data Analysis
0004 // Copyright (C) 2008-2025 The YODA collaboration (see AUTHORS for details)
0005 //
0006 #ifndef YODA_Counter_h
0007 #define YODA_Counter_h
0008 
0009 #include "YODA/AnalysisObject.h"
0010 #include "YODA/Fillable.h"
0011 #include "YODA/Dbn.h"
0012 #include "YODA/Point.h"
0013 #include "YODA/Estimate0D.h"
0014 #include "YODA/Scatter.h"
0015 #include "YODA/Exceptions.h"
0016 #include <vector>
0017 #include <string>
0018 #include <map>
0019 #include <tuple>
0020 #include <memory>
0021 
0022 namespace YODA {
0023 
0024 
0025   /// A weighted counter.
0026   class Counter : public AnalysisObject, public Fillable {
0027   public:
0028 
0029 
0030     using FillType = std::tuple<>;
0031     using Ptr = std::shared_ptr<Counter>;
0032     using AnalysisObject::operator =;
0033 
0034     /// @name Constructors
0035     /// @{
0036 
0037     /// Default constructor
0038     Counter(const std::string& path="", const std::string& title="")
0039       : AnalysisObject("Counter", path, title) { }
0040 
0041 
0042     /// @brief Constructor accepting an explicit Dbn0D.
0043     ///
0044     /// Intended both for internal persistency and user use.
0045     Counter(const Dbn0D& dbn,
0046             const std::string& path="", const std::string& title="")
0047             : AnalysisObject("Counter", path, title), _dbn(dbn) { }
0048 
0049 
0050     /// @brief Constructor accepting an explicit rvalue Dbn0D.
0051     ///
0052     /// Intended both for internal persistency and user use.
0053     Counter(Dbn0D&& dbn,
0054             const std::string& path="", const std::string& title="")
0055             : AnalysisObject("Counter", path, title), _dbn(std::move(dbn)) { }
0056 
0057 
0058     /// @brief Constructor accepting a double (treated as the weight of a single fill).
0059     ///
0060     /// Intended for user convenience only, so Counter can be treated as a number.
0061     Counter(double w,
0062             const std::string& path="", const std::string& title="")
0063             : AnalysisObject("Counter", path, title) { _dbn.fill(w); }
0064 
0065 
0066     /// Copy constructor with optional new path
0067     /// @todo Don't copy the path?
0068     Counter(const Counter& c, const std::string& path="")
0069             : AnalysisObject("Counter", (path.size() == 0) ? c.path() : path, c, c.title()),
0070               _dbn(c._dbn) { }
0071 
0072     /// Move constructor with optional new path
0073     Counter(Counter&& c, const std::string& path="")
0074             : AnalysisObject("Counter", (path.size() == 0) ? c.path() : path, c, c.title()),
0075               _dbn(std::move(c._dbn)) { }
0076 
0077 
0078     /// Assignment operator
0079     Counter& operator = (const Counter& c) noexcept {
0080       if (this != &c) {
0081         AnalysisObject::operator = (c);
0082         _dbn = c._dbn;
0083       }
0084       return *this;
0085     }
0086 
0087     /// Move operator
0088     Counter& operator = (Counter&& c) noexcept {
0089       if (this != &c) {
0090         AnalysisObject::operator = (c);
0091         _dbn = std::move(c._dbn);
0092       }
0093       return *this;
0094     }
0095 
0096     /// Make a copy on the stack
0097     Counter clone() const {
0098       return Counter(*this);
0099     }
0100 
0101     /// Make a copy on the heap, via 'new'
0102     Counter* newclone() const {
0103       return new Counter(*this);
0104     }
0105 
0106     /// @}
0107 
0108 
0109     /// @name I/O
0110     /// @{
0111 
0112     // @brief Render information about this AO
0113     void _renderYODA(std::ostream& os, const int width = 13) const noexcept {
0114       os << std::setw(width) << std::left << "# sumW" << "\t"
0115          << std::setw(width) << std::left << "sumW2"  << "\t"
0116          << "numEntries\n";
0117       os << std::setw(width) << std::left << sumW() << "\t";
0118       os << std::setw(width) << std::left << sumW2() << "\t";
0119       os << std::setw(width) << std::left << numEntries() << "\n";
0120     }
0121 
0122     // @brief Render scatter-like information about this AO
0123     void _renderFLAT(std::ostream& os, const int width) const noexcept {
0124       Scatter1D tmp = this->mkScatter();
0125       tmp._renderYODA(os, width);
0126     }
0127 
0128     /// @}
0129 
0130 
0131     /// @name Dimensions
0132     /// @{
0133 
0134     /// @brief Total dimension of this data object
0135     size_t dim() const noexcept { return 1; }
0136 
0137     /// @brief Fill dimension of this data object
0138     size_t fillDim() const noexcept { return 0; }
0139 
0140     /// @}
0141 
0142 
0143     /// @name Modifiers
0144     /// @{
0145 
0146     /// Fill histo by value and weight
0147     virtual int fill(double weight=1.0, double fraction=1.0) {
0148       _dbn.fill(weight, fraction);
0149       return 0;
0150     }
0151 
0152     virtual int fill(FillType&&, double weight=1.0, double fraction=1.0) {
0153       return fill(weight, fraction);
0154     }
0155 
0156     /// @brief Reset the histogram.
0157     ///
0158     /// Keep the binning but set all bin contents and related quantities to zero
0159     virtual void reset() {
0160       _dbn.reset();
0161     }
0162 
0163 
0164     /// Rescale as if all fill weights had been different by factor @a scalefactor.
0165     void scaleW(double scalefactor) {
0166       setAnnotation("ScaledBy", annotation<double>("ScaledBy", 1.0) * scalefactor);
0167       _dbn.scaleW(scalefactor);
0168     }
0169 
0170     /// @}
0171 
0172 
0173     /// @name Data access
0174     /// @{
0175 
0176     /// Get the number of fills
0177     double numEntries(bool=false) const { return _dbn.numEntries(); }
0178 
0179     /// Get the effective number of fills
0180     double effNumEntries(bool=false) const { return _dbn.effNumEntries(); }
0181 
0182     /// Get the sum of weights
0183     double sumW(bool=false) const { return _dbn.sumW(); }
0184 
0185     /// Get the sum of squared weights
0186     double sumW2(bool=false) const { return _dbn.sumW2(); }
0187 
0188     /// Get the value
0189     double val(bool=false) const { return sumW(); }
0190 
0191     /// Get the uncertainty on the value
0192     /// @todo Implement on Dbn0D and feed through to this and Dbn1D, 2D, etc.
0193     double err() const { return sqrt(sumW2()); }
0194 
0195     /// Get the relative uncertainty on the value
0196     /// @todo Implement on Dbn0D and feed through to this and Dbn1D, 2D, etc.
0197     // double err() const { return _dbn.err(); }
0198     double relErr() const {
0199       /// @todo Throw excp if sumW2 is 0?
0200       return sumW2() != 0 ? err()/sumW() : 0;
0201     }
0202 
0203     /// @}
0204 
0205 
0206     /// @name Internal state access and modification (mainly for persistency use)
0207     /// @{
0208 
0209     /// Get the internal distribution object
0210     const Dbn0D& dbn() const {
0211       return _dbn;
0212     }
0213 
0214     /// Set the internal distribution object: CAREFUL!
0215     void setDbn(const Dbn0D& dbn) {
0216       _dbn = dbn;
0217     }
0218 
0219     /// Set the internal distribution object: CAREFUL!
0220     void setDbn(Dbn0D&& dbn) {
0221       _dbn = std::move(dbn);
0222     }
0223 
0224     /// Set the internal distribution object: CAREFUL!
0225     void set(const double numEntries, const double sumW, const double sumW2) {
0226       _dbn.set(numEntries, sumW, sumW2);
0227     }
0228 
0229     // /// Set the whole object state
0230     // void setState(const Dbn0D& dbn, const AnalysisObject::Annotations& anns=AnalysisObject::Annotations()) {
0231     //   setDbn(dbn);
0232     //   setAnnotations(anns);
0233     // }
0234 
0235     /// @}
0236 
0237 
0238     /// @name Adding and subtracting counters
0239     /// @{
0240 
0241     /// Add another counter to this
0242     Counter& operator += (const Counter& toAdd) {
0243       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0244       _dbn += toAdd._dbn;
0245       return *this;
0246     }
0247     //
0248     Counter& operator += (Counter&& toAdd) {
0249       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0250       _dbn += std::move(toAdd._dbn);
0251       return *this;
0252     }
0253 
0254     /// Subtract another counter from this
0255     Counter& operator -= (const Counter& toSubtract) {
0256       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0257       _dbn -= toSubtract._dbn;
0258       return *this;
0259     }
0260     //
0261     Counter& operator -= (Counter&& toSubtract) {
0262       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0263       _dbn -= std::move(toSubtract._dbn);
0264       return *this;
0265     }
0266 
0267     /// Increment as if by a fill of weight = 1
0268     /// @note This is post-increment only, i.e. cn++ not ++cn
0269     Counter& operator ++ () {
0270       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0271       *this += 1;
0272       return *this;
0273     }
0274 
0275     /// Increment as if by a fill of weight = -1
0276     /// @note This is post-decrement only, i.e. cn-- not --cn
0277     Counter& operator -- () {
0278       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0279       *this -= 1;
0280       return *this;
0281     }
0282 
0283     /// Scale by a double (syntactic sugar for @c scaleW(s))
0284     Counter& operator *= (double s) {
0285       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0286       scaleW(s);
0287       return *this;
0288     }
0289 
0290     /// Inverse-scale by a double (syntactic sugar for @c scaleW(1/s))
0291     Counter& operator /= (double s) {
0292       if (AO::hasAnnotation("ScaledBy")) AO::rmAnnotation("ScaledBy");
0293       scaleW(1/s);
0294       return *this;
0295     }
0296 
0297     /// @}
0298 
0299     /// @name Type reduction
0300     /// @{
0301 
0302     /// @brief Make an inert estimate with the value and error of the counter
0303     inline Estimate0D mkEstimate(const std::string& path = "", const std::string& source = "") const {
0304       Estimate0D rtn;
0305       for (const std::string& a : annotations()) {
0306         if (a != "Type")  rtn.setAnnotation(a, annotation(a));
0307       }
0308       rtn.setAnnotation("Path", path);
0309 
0310       rtn.setVal(val());
0311       if (numEntries()) { // only set uncertainty for filled Dbns
0312         rtn.setErr(err(), source);
0313       }
0314       return rtn;
0315     }
0316 
0317     /// @brief Make an inert estimate with the value and error of the counter
0318     ///
0319     /// @note The uncertainty on the effective entries is currently set to zero,
0320     /// as we don't have a formula for it! So the source arg (for uniformity)
0321     /// is currently ignored.
0322     inline Estimate0D mkEffNumEntries(const std::string& path = "", const std::string& /* source */ = "") const {
0323       Estimate0D rtn;
0324       for (const std::string& a : annotations()) {
0325         if (a != "Type")  rtn.setAnnotation(a, annotation(a));
0326       }
0327       rtn.setAnnotation("Path", path);
0328 
0329       rtn.setVal(effNumEntries());
0330       /// @todo What's the uncertainty on Neff?
0331       // if (numEntries()) { // only set uncertainty for filled Dbns
0332       //   rtn.setErr(err(), source);
0333       // }
0334 
0335       return rtn;
0336     }
0337 
0338     /// @brief Make an inert scatter with the value and error of the counter
0339     inline Scatter1D mkScatter(const std::string& path = "") const {
0340       Scatter1D rtn;
0341       for (const std::string& a : annotations()) {
0342         if (a != "Type")  rtn.setAnnotation(a, annotation(a));
0343       }
0344       rtn.setAnnotation("Path", path);
0345 
0346       rtn.addPoint(Point1D(val(), err()));
0347       return rtn;
0348     }
0349 
0350     /// @brief Return an inert version of the analysis object (e.g. scatter, estimate)
0351     AnalysisObject* mkInert(const std::string& path = "",
0352                             const std::string& source = "") const noexcept {
0353       return mkEstimate(path, source).newclone();
0354     }
0355 
0356     /// @}
0357 
0358     /// @name MPI (de-)serialisation
0359     //@{
0360 
0361     size_t lengthContent(bool = false) const noexcept {
0362       return _dbn._lengthContent();
0363     }
0364 
0365     std::vector<double> serializeContent(bool = false) const noexcept {
0366       return _dbn._serializeContent();
0367     }
0368 
0369     void deserializeContent(const std::vector<double>& data) {
0370       _dbn._deserializeContent(data);
0371     }
0372 
0373     // @}
0374 
0375   private:
0376 
0377     /// @name Data
0378     /// @{
0379 
0380     /// Contained 0D distribution
0381     Dbn0D _dbn;
0382 
0383     /// @}
0384 
0385   };
0386 
0387 
0388   /// @name Combining counters: global operators
0389   /// @{
0390 
0391   /// Add two counters
0392   inline Counter operator + (Counter first, const Counter& second) {
0393     first += second;
0394     return first;
0395   }
0396 
0397   /// Add two counters
0398   inline Counter operator + (Counter first, Counter&& second) {
0399     first += std::move(second);
0400     return first;
0401   }
0402 
0403   /// Subtract two counters
0404   inline Counter operator - (Counter first, const Counter& second) {
0405     first -= second;
0406     return first;
0407   }
0408 
0409   /// Subtract two counters
0410   inline Counter operator - (Counter first, Counter&& second) {
0411     first -= std::move(second);
0412     return first;
0413   }
0414 
0415   /// Divide two counters, with an uncorrelated error treatment
0416   Estimate0D divide(const Counter& numer, const Counter& denom);
0417 
0418   /// Divide two Counter objects
0419   inline Estimate0D operator / (const Counter& numer, const Counter& denom) {
0420     return divide(numer, denom);
0421   }
0422 
0423   /// Divide two Counter objects
0424   inline Estimate0D operator / (Counter&& numer, const Counter& denom) {
0425     return divide(std::move(numer), denom);
0426   }
0427 
0428   /// Divide two Counter objects
0429   inline Estimate0D operator / (const Counter& numer, Counter&& denom) {
0430     return divide(numer, std::move(denom));
0431   }
0432 
0433   /// Divide two Counter objects
0434   inline Estimate0D operator / (Counter&& numer, Counter&& denom) {
0435     return divide(std::move(numer), std::move(denom));
0436   }
0437 
0438   /// @todo Add divide functions/operators on pointers
0439 
0440   /// @brief Calculate an efficiency ratio of two counters
0441   ///
0442   /// Note that an efficiency is not the same thing as a standard division of two
0443   /// histograms: the errors must be treated as correlated
0444   Estimate0D efficiency(const Counter& accepted, const Counter& total);
0445 
0446   /// @}
0447 
0448 
0449 }
0450 
0451 #endif