19#include "HistogramAccumulator.hh"
25#include "BDSBH4DBase.hh"
43 const std::string& resultHistNameIn,
44 const std::string& resultHistTitleIn):
45 HistogramAccumulator(baseHistogram, (int)baseHistogram->GetDimension(), resultHistNameIn, resultHistTitleIn)
50 const std::string& resultHistNameIn,
51 const std::string& resultHistTitleIn):
52 nDimensions(nDimensionsIn),
55 resultHistName(resultHistNameIn),
56 resultHistTitle(resultHistTitleIn),
67 mean =
dynamic_cast<TH1D*
>(baseHistogram->Clone(meanName.c_str()));
68 variance =
dynamic_cast<TH1D*
>(baseHistogram->Clone(variName.c_str()));
69 result =
dynamic_cast<TH1D*
>(baseHistogram->Clone(
resultHistName.c_str()));
74 mean =
dynamic_cast<TH2D*
>(baseHistogram->Clone(meanName.c_str()));
75 variance =
dynamic_cast<TH2D*
>(baseHistogram->Clone(variName.c_str()));
76 result =
dynamic_cast<TH2D*
>(baseHistogram->Clone(
resultHistName.c_str()));
81 mean =
dynamic_cast<TH3D*
>(baseHistogram->Clone(meanName.c_str()));
82 variance =
dynamic_cast<TH3D*
>(baseHistogram->Clone(variName.c_str()));
83 result =
dynamic_cast<TH3D*
>(baseHistogram->Clone(
resultHistName.c_str()));
89 mean =
dynamic_cast<BDSBH4DBase*
>(baseHistogram)->Clone(meanName.c_str());
90 variance =
dynamic_cast<BDSBH4DBase*
>(baseHistogram)->Clone(variName.c_str());
96 {
throw std::domain_error(
"Invalid number of dimensions");
break;}
98 if (mean && variance && result)
108 mean->SetTitle(meanName.c_str());
109 variance->SetTitle(variName.c_str());
116 static_cast<BDSBH4DBase*
>(variance)->Reset_BDSBH4D();
117 static_cast<BDSBH4DBase*
>(result)->Reset_BDSBH4D();
120 static_cast<BDSBH4DBase*
>(mean)->SetTitle(meanName.c_str());
121 static_cast<BDSBH4DBase*
>(variance)->SetTitle(variName.c_str());
140 const double error = 0;
141 const unsigned long nEntriesToAccumulate = 1;
148 TH1D* h1 =
dynamic_cast<TH1D*
>(mean);
149 TH1D* h1e =
dynamic_cast<TH1D*
>(variance);
150 TH1D* ht =
dynamic_cast<TH1D*
>(newValue);
151 for (
int j = 0; j <= h1->GetNbinsX() + 1; ++j)
154 h1e->GetBinContent(j),
155 ht->GetBinContent(j),
156 error,
n, nEntriesToAccumulate,
158 h1->SetBinContent(j, newMean);
159 h1e->SetBinContent(j, newVari);
165 TH2D* h1 =
dynamic_cast<TH2D*
>(mean);
166 TH2D* h1e =
dynamic_cast<TH2D*
>(variance);
167 TH2D* ht =
dynamic_cast<TH2D*
>(newValue);
168 for (
int j = 0; j <= h1->GetNbinsX() + 1; ++j)
170 for (
int k = 0; k <= h1->GetNbinsY() + 1; ++k)
173 h1e->GetBinContent(j,k),
174 ht->GetBinContent(j,k),
175 error,
n, nEntriesToAccumulate,
177 h1->SetBinContent(j, k, newMean);
178 h1e->SetBinContent(j, k, newVari);
185 TH3D* h1 =
dynamic_cast<TH3D*
>(mean);
186 TH3D* h1e =
dynamic_cast<TH3D*
>(variance);
187 TH3D* ht =
dynamic_cast<TH3D*
>(newValue);
188 for (
int j = 0; j <= h1->GetNbinsX() + 1; ++j)
190 for (
int k = 0; k <= h1->GetNbinsY() + 1; ++k)
192 for (
int l = 0; l <= h1->GetNbinsZ() + 1; ++l)
195 h1e->GetBinContent(j,k,l),
196 ht->GetBinContent(j,k,l),
197 error,
n, nEntriesToAccumulate,
199 h1->SetBinContent(j, k, l, newMean);
200 h1e->SetBinContent(j, k, l, newVari);
212 for (
int j = -1; j <= h1->GetNbinsX(); ++j)
214 for (
int k = -1; k <= h1->GetNbinsY(); ++k)
216 for (
int l = -1; l <= h1->GetNbinsZ(); ++l)
218 for (
int e = -1; e <= h1->GetNbinsE(); ++e)
223 error,
n, nEntriesToAccumulate,
225 h1->Set_BDSBH4D(j,k,l,e, newMean);
226 h1e->Set_BDSBH4D(j,k,l,e, newVari);
244 const double nD = (double)
n;
245 const double factor = std::sqrt(1./(nD * (nD - 1)));
256 for (
int j = 0; j <= result->GetNbinsX() + 1; ++j)
258 mn = mean->GetBinContent(j);
259 var = variance->GetBinContent(j);
260 err =
n > 1 ? factor*std::sqrt(var) : 0;
261 result->SetBinContent(j, mn);
262 result->SetBinError(j, err);
268 for (
int j = 0; j <= result->GetNbinsX() + 1; ++j)
270 for (
int k = 0; k <= result->GetNbinsY() + 1; ++k)
272 mn = mean->GetBinContent(j,k);
273 var = variance->GetBinContent(j, k);
274 err =
n > 1 ? factor*std::sqrt(var) : 0;
275 result->SetBinContent(j, k, mn);
276 result->SetBinError(j, k, err);
283 for (
int j = 0; j <= result->GetNbinsX() + 1; ++j)
285 for (
int k = 0; k <= result->GetNbinsY() + 1; ++k)
287 for (
int l = 0; l <= result->GetNbinsZ() + 1; ++l)
289 mn = mean->GetBinContent(j,k,l);
290 var = variance->GetBinContent(j, k, l);
291 err =
n > 1 ? factor*std::sqrt(var) : 0;
292 result->SetBinContent(j,k,l, mn);
293 result->SetBinError(j,k,l, err);
302 auto histCast =
dynamic_cast<BDSBH4DBase*
>(result);
304 auto varCast =
dynamic_cast<BDSBH4DBase*
>(variance);
306 int nBinsX = histCast->GetNbinsX();
307 int nBinsY = histCast->GetNbinsY();
308 int nBinsZ = histCast->GetNbinsZ();
309 int nBinsE = histCast->GetNbinsE();
310 for (
int j = -1; j <= nBinsX; ++j)
312 for (
int k = -1; k <= nBinsY; ++k)
314 for (
int l = -1; l <= nBinsZ; ++l)
316 for (
int e = -1; e <= nBinsE; ++e)
318 mn = mnCast->At(j, k, l, e);
319 var = varCast->At(j, k, l, e);
320 err =
n > 1 ? factor*std::sqrt(var) : 0;
321 resCast->Set_BDSBH4D(j, k, l, e, mn);
322 resCast->SetError_BDSBH4D(j, k, l, e, err);
334 {
dynamic_cast<BDSBH4DBase*
>(result)->SetEntries_BDSBH4D((
double)
n);}
336 {result->SetEntries((
double)
n);}
345 unsigned long nEntriesAccumulated,
348 double& newVari)
const
350 newMean = oldMean + ((x - oldMean) / (
double)nEntriesAccumulated);
351 newVari = oldVari + ((x - oldMean) * (x - newMean));
Base class for the 4D histogram classes.
Class to accumulate and merge histograms in different ways.
const std::string resultHistName
Name for resultant histogram.
int nDimensions
Number of dimensions.
virtual ~HistogramAccumulator()
virtual void AccumulateSingleValue(double oldMean, double oldVari, double x, double xVari, unsigned long nEntriesAccumulated, unsigned long nEntriesToAccumulate, double &newMean, double &newVari) const
const std::string resultHistTitle
Title for resultant histogram.
virtual TH1 * Terminate()
HistogramAccumulator()
Default constructor only for ROOT reflexivity - not intended for use.
virtual void Accumulate(TH1 *newValue)