ntupleStats.hh

Go to the documentation of this file.

#ifndef NPSTAT_NTUPLESTATS_HH_
#define NPSTAT_NTUPLESTATS_HH_
 
/*!
// \file ntupleStats.hh
//
// \brief Convenience functions for calculating ntuple means and covariances
//
// Author: I. Volobouev
//
// October 2022
*/
 
#include "npstat/stat/AbsNtuple.hh"
#include "npstat/stat/MultivariateSumsqAccumulator.hh"
#include "npstat/stat/MultivariateWeightedSumsqAccumulator.hh"
 
namespace npstat {
    template<typename Numeric>
    inline void ntupleMean(const AbsNtuple<Numeric>& nt,
                           double* mean, const unsigned long lengthMean)
    {
        const unsigned long nCols = nt.nColumns();
        assert(lengthMean >= nCols);
        assert(mean);
 
        MultivariateSumAccumulator<> sums;
        nt.cycleOverRows(sums);
        for (unsigned long i=0; i<nCols; ++i)
            mean[i] = sums.mean(i);
    }
 
    template<typename Numeric, class WeightCalc>
    inline void ntupleWeightedMean(const AbsNtuple<Numeric>& nt,
                                   const WeightCalc& wcalc,
                                   double* mean, const unsigned long lengthMean)
    {
        const unsigned long nCols = nt.nColumns();
        assert(lengthMean >= nCols);
        assert(mean);
 
        MultivariateWeightedSumAccumulator<> sums;
        nt.weightedCycleOverRows(sums, wcalc);
        for (unsigned long i=0; i<nCols; ++i)
            mean[i] = sums.mean(i);
    }
 
    template<typename Numeric>
    inline Matrix<double> ntupleCovariance(const AbsNtuple<Numeric>& nt)
    {
        MultivariateSumAccumulator<> sums;
        nt.cycleOverRows(sums);
        MultivariateSumsqAccumulator<> sumsqs(sums);
        nt.cycleOverRows(sumsqs);
        return sumsqs.covMat();
    }
 
    template<typename Numeric, class WeightCalc>
    inline Matrix<double> ntupleWeightedCovariance(
        const AbsNtuple<Numeric>& nt, const WeightCalc& wcalc)
    {
        MultivariateWeightedSumAccumulator<> sums;
        nt.weightedCycleOverRows(sums, wcalc);
        MultivariateWeightedSumsqAccumulator<> sumsqs(sums);
        nt.weightedCycleOverRows(sumsqs, wcalc);
        return sumsqs.covMat();
    }
 
#ifdef SWIG
    template<typename Numeric>
    inline std::vector<double> ntupleMean2(const AbsNtuple<Numeric>& nt)
    {
        const unsigned long nCols = nt.nColumns();
        std::vector<double> tmp(nCols);
        ntupleMean(nt, &tmp[0], nCols);
        return tmp;
    }
 
    template<typename Numeric, class WeightCalc>
    inline std::vector<double> ntupleWeightedMean2(const AbsNtuple<Numeric>& nt,
                                                   const WeightCalc& wcalc)
    {
        const unsigned long nCols = nt.nColumns();
        std::vector<double> tmp(nCols);
        ntupleWeightedMean(nt, wcalc, &tmp[0], nCols);
        return tmp;
    }
#endif // SWIG
}
 
#endif // NPSTAT_NTUPLESTATS_HH_