19#include "HistogramAccumulatorMerge.hh"
25#include "BDSBH4DBase.hh"
38 const std::string& resultHistNameIn,
39 const std::string& resultHistTitleIn):
56 unsigned long oldEntries;
57 unsigned long newEntries;
61 oldEntries = (
unsigned long)
static_cast<BDSBH4DBase*
>(mean)->GetEntries_BDSBH4D();
62 newEntries = (
unsigned long)
static_cast<BDSBH4DBase*
>(newValue)->GetEntries_BDSBH4D();
66 oldEntries = (
unsigned long)mean->GetEntries();
67 newEntries = (
unsigned long)newValue->GetEntries();
70 unsigned long newTotalEntries = oldEntries + newEntries;
71 const double nD = (double)newEntries;
72 const double factor = nD * (nD - 1);
78 TH1D* h1 =
dynamic_cast<TH1D*
>(mean);
79 TH1D* h1e =
dynamic_cast<TH1D*
>(variance);
80 TH1D* ht =
dynamic_cast<TH1D*
>(newValue);
81 for (
int j = 0; j <= h1->GetNbinsX() + 1; ++j)
83 var = std::pow(ht->GetBinError(j), 2) * factor;
85 h1e->GetBinContent(j),
88 oldEntries, newEntries,
90 h1->SetBinContent(j, newMean);
91 h1e->SetBinContent(j, newVari);
97 TH2D* h1 =
dynamic_cast<TH2D*
>(mean);
98 TH2D* h1e =
dynamic_cast<TH2D*
>(variance);
99 TH2D* ht =
dynamic_cast<TH2D*
>(newValue);
100 for (
int j = 0; j <= h1->GetNbinsX() + 1; ++j)
102 for (
int k = 0; k <= h1->GetNbinsY() + 1; ++k)
104 var = std::pow(ht->GetBinError(j,k), 2) * factor;
106 h1e->GetBinContent(j,k),
107 ht->GetBinContent(j,k),
109 oldEntries, newEntries,
111 h1->SetBinContent(j, k, newMean);
112 h1e->SetBinContent(j, k, newVari);
119 TH3D* h1 =
dynamic_cast<TH3D*
>(mean);
120 TH3D* h1e =
dynamic_cast<TH3D*
>(variance);
121 TH3D* ht =
dynamic_cast<TH3D*
>(newValue);
122 for (
int j = 0; j <= h1->GetNbinsX() + 1; ++j)
124 for (
int k = 0; k <= h1->GetNbinsY() + 1; ++k)
126 for (
int l = 0; l <= h1->GetNbinsZ() + 1; ++l)
128 var = std::pow(ht->GetBinError(j,k,l), 2) * factor;
130 h1e->GetBinContent(j,k,l),
131 ht->GetBinContent(j,k,l),
133 oldEntries, newEntries,
135 h1->SetBinContent(j, k, l, newMean);
136 h1e->SetBinContent(j, k, l, newVari);
148 for (
int j = -1; j <= h1->GetNbinsX(); ++j)
150 for (
int k = -1; k <= h1->GetNbinsY(); ++k)
152 for (
int l = -1; l <= h1->GetNbinsZ(); ++l)
154 for (
int e = -1; e <= h1->GetNbinsE(); ++e)
156 var = std::pow(ht->AtError(j,k,l,e), 2) * factor;
161 oldEntries, newEntries,
163 h1->Set_BDSBH4D(j, k, l, e, newMean);
164 h1e->Set_BDSBH4D(j, k, l, e, newVari);
177 dynamic_cast<BDSBH4DBase*
>(mean)->SetEntries_BDSBH4D(newTotalEntries);
178 dynamic_cast<BDSBH4DBase*
>(variance)->SetEntries_BDSBH4D(newTotalEntries);
182 mean->SetEntries(newTotalEntries);
183 variance->SetEntries(newTotalEntries);
192 unsigned long nEntriesAccumulated,
193 unsigned long nEntriesToAccumulate,
195 double& newVari)
const
197 double dMean = x - oldMean;
198 double dMean2 = std::pow(dMean, 2);
199 double nA = (double)nEntriesAccumulated;
200 double nB = (double)nEntriesToAccumulate;
203 newMean = oldMean + nB * (dMean / nT);
205 double q = (nA * nB) * (dMean2 / nT);
206 newVari = oldVari + xVari + q;
Base class for the 4D histogram classes.
Class to combine histograms of mean with error on mean.
virtual void AccumulateSingleValue(double oldMean, double oldVari, double x, double xVari, unsigned long nEntriesAccumulated, unsigned long nEntriesToAccumulate, double &newMean, double &newVari) const
virtual void Accumulate(TH1 *newValue)
HistogramAccumulatorMerge()
Default constructor only for ROOT reflexivity - not intended for use.
Class to accumulate and merge histograms in different ways.
int nDimensions
Number of dimensions.