12 #ifndef ROOT_TProfileHelper
13 #define ROOT_TProfileHelper
84 Int_t nz = p->GetNbinsZ();
86 if ( nx != p1->GetNbinsX() || nx != p2->GetNbinsX() ||
87 ny != p1->GetNbinsY() || ny != p2->GetNbinsY() ||
88 nz != p1->GetNbinsZ() || nz != p2->GetNbinsZ() ) {
89 Error(
"TProfileHelper::Add",
"Attempt to add profiles with different number of bins");
96 p->fEntries = ac1*p1->GetEntries() + ac2*p2->GetEntries();
99 for (i=0;i<
TH1::kNstat;i++) {s0[i] = s1[i] = s2[i] = 0;}
104 if (i == 1) s0[i] = c1*c1*s1[i] + c2*c2*s2[i];
105 else s0[i] = ac1*s1[i] + ac2*s2[i];
116 if (p->fBinSumw2.fN == 0 && (p1->fBinSumw2.fN != 0 || p2->fBinSumw2.fN != 0) ) p->Sumw2();
118 if (ew1 == 0) ew1 = en1;
119 if (ew2 == 0) ew2 = en2;
120 for (bin =0;bin< p->fN;bin++) {
121 p->fArray[bin] = c1*cu1[bin] + c2*cu2[bin];
122 p->fSumw2.fArray[bin] = ac1*er1[bin] + ac2*er2[bin];
123 p->fBinEntries.fArray[bin] = ac1*en1[bin] + ac2*en2[bin];
124 if (p->fBinSumw2.fN ) p->fBinSumw2.fArray[bin] = ac1*ac1*ew1[bin] + ac2*ac2*ew2[bin];
129 template <
typename T>
136 p->fBinEntries.Set(p->fNcells);
137 p->fSumw2.Set(p->fNcells);
142 template <
typename T>
153 if (p->fBuffer) p->BufferEmpty();
155 if (bin < 0 || bin >= p->fNcells)
return 0;
156 double sumOfWeights = p->fBinEntries.fArray[bin];
157 if ( p->fBinSumw2.fN == 0 || p->fBinSumw2.fN != p->fNcells) {
162 double sumOfWeightsSquare = p->fBinSumw2.fArray[bin];
163 return ( sumOfWeightsSquare > 0 ? sumOfWeights * sumOfWeights / sumOfWeightsSquare : 0 );
166 template <
typename T>
204 allHaveLimits = allHaveLimits && hasLimits;
212 if (firstNonEmptyHist ) {
215 if (!p->SameLimitsAndNBins(p->fXaxis, *(h->GetXaxis())) )
216 p->fXaxis.Set(h->GetXaxis()->GetNbins(), h->GetXaxis()->GetXmin(),h->GetXaxis()->GetXmax());
217 if (!p->SameLimitsAndNBins(p->fYaxis, *(h->GetYaxis())) )
218 p->fYaxis.Set(h->GetYaxis()->GetNbins(), h->GetYaxis()->GetXmin(),h->GetYaxis()->GetXmax());
219 if (!p->SameLimitsAndNBins(p->fZaxis, *(h->GetZaxis())) )
220 p->fZaxis.Set(h->GetZaxis()->GetNbins(), h->GetZaxis()->GetXmin(),h->GetZaxis()->GetXmax());
222 firstNonEmptyHist =
kFALSE;
228 if (!initialLimitsFound) {
229 initialLimitsFound =
kTRUE;
230 newXAxis.
Set(h->GetXaxis()->GetNbins(), h->GetXaxis()->GetXmin(),
231 h->GetXaxis()->GetXmax());
232 if ( p->GetDimension() >= 2 )
233 newYAxis.
Set(h->GetYaxis()->GetNbins(), h->GetYaxis()->GetXmin(),
234 h->GetYaxis()->GetXmax());
235 if ( p->GetDimension() >= 3 )
236 newZAxis.
Set(h->GetZaxis()->GetNbins(), h->GetZaxis()->GetXmin(),
237 h->GetZaxis()->GetXmax());
241 if (!p->SameLimitsAndNBins(newXAxis, *(h->GetXaxis())) ||
242 !p->SameLimitsAndNBins(newYAxis, *(h->GetYaxis())) ||
243 !p->SameLimitsAndNBins(newZAxis, *(h->GetZaxis())) ) {
250 if (!p->RecomputeAxisLimits(newXAxis, *(h->GetXaxis()))) {
251 Error(
"TProfileHelper::Merge",
"Cannot merge profiles %d dim - limits are inconsistent:\n "
252 "first: (%d, %f, %f), second: (%d, %f, %f)",p->GetDimension(),
254 h->GetXaxis()->GetNbins(), h->GetXaxis()->GetXmin(),
255 h->GetXaxis()->GetXmax());
258 if (p->GetDimension() >= 2 && !p->RecomputeAxisLimits(newYAxis, *(h->GetYaxis()))) {
259 Error(
"TProfileHelper::Merge",
"Cannot merge profiles %d dim - limits are inconsistent:\n "
260 "first: (%d, %f, %f), second: (%d, %f, %f)",p->GetDimension(),
262 h->GetYaxis()->GetNbins(), h->GetYaxis()->GetXmin(),
263 h->GetYaxis()->GetXmax());
266 if (p->GetDimension() >= 3 && !p->RecomputeAxisLimits(newZAxis, *(h->GetZaxis()))) {
267 Error(
"TProfileHelper::Merge",
"Cannot merge profiles %d dim - limits are inconsistent:\n "
268 "first: (%d, %f, %f), second: (%d, %f, %f)",p->GetDimension(),
270 h->GetZaxis()->GetNbins(), h->GetZaxis()->GetXmin(),
271 h->GetZaxis()->GetXmax());
277 }
while ( ( h = dynamic_cast<T*> (
next() ) ) !=
NULL );
278 if (!h && (*next) ) {
279 Error(
"TProfileHelper::Merge",
"Attempt to merge object of class: %s to a %s",
280 (*next)->ClassName(),p->ClassName());
291 if (!allSameLimits) {
296 hclone = (
T*)p->IsA()->New();
298 hclone->SetDirectory(0);
307 if (!allSameLimits && initialLimitsFound) {
315 if (!allHaveLimits) {
317 while ( (h = dynamic_cast<T*> (
next()) ) ) {
318 if (h->GetXaxis()->GetXmin() >= h->GetXaxis()->GetXmax() && h->fBuffer) {
322 for (
Int_t i = 0; i < nbentries; i++)
323 if ( p->GetDimension() == 3 ) {
324 v[0] = h->fBuffer[5*i + 2];
325 v[1] = h->fBuffer[5*i + 3];
326 v[2] = h->fBuffer[5*i + 4];
327 v[3] = h->fBuffer[5*i + 5];
328 v[4] = h->fBuffer[5*i + 1];
330 }
else if ( p->GetDimension() == 2 ) {
331 v[0] = h->fBuffer[4*i + 2];
332 v[1] = h->fBuffer[4*i + 3];
333 v[2] = h->fBuffer[4*i + 4];
334 v[3] = h->fBuffer[4*i + 1];
338 else if ( p->GetDimension() == 1 ) {
339 v[0] = h->fBuffer[3*i + 2];
340 v[1] = h->fBuffer[3*i + 3];
341 v[2] = h->fBuffer[3*i + 1];
347 if (!initialLimitsFound) {
352 return (
Int_t) p->GetEntries();
360 p->GetStats(totstats);
362 Bool_t canExtend = p->CanExtendAllAxes();
365 while ( (h=static_cast<T*>(
next())) ) {
368 if (h->GetXaxis()->GetXmin() < h->GetXaxis()->GetXmax()) {
372 totstats[i] += stats[i];
373 nentries += h->GetEntries();
375 for (
Int_t hbin = 0; hbin < h->fN; ++hbin ) {
377 if (!allSameLimits) {
382 if ( h->GetW()[hbin] != 0 && (h->IsBinUnderflow(hbin) || h->IsBinOverflow(hbin)) ) {
384 Error(
"TProfileHelper::Merge",
"Cannot merge profiles - they have"
385 " different limits and undeflows/overflows are present."
386 " The initial profile is now broken!");
389 Int_t hbinx, hbiny, hbinz;
390 h->GetBinXYZ(hbin, hbinx, hbiny, hbinz);
392 pbin = p->GetBin( p->fXaxis.FindBin( h->GetXaxis()->GetBinCenter(hbinx) ),
393 p->fYaxis.FindBin( h->GetYaxis()->GetBinCenter(hbiny) ),
394 p->fZaxis.FindBin( h->GetZaxis()->GetBinCenter(hbinz) ) );
398 p->fArray[pbin] += h->GetW()[hbin];
399 p->fSumw2.fArray[pbin] += h->GetW2()[hbin];
400 p->fBinEntries.fArray[pbin] += h->GetB()[hbin];
401 if (p->fBinSumw2.fN) {
402 if ( h->GetB2() ) p->fBinSumw2.fArray[pbin] += h->GetB2()[hbin];
403 else p->fBinSumw2.fArray[pbin] += h->GetB()[hbin];
411 p->PutStats(totstats);
412 p->SetEntries(nentries);
420 template <
typename T>
435 if (axis->
GetNbins() <= 0)
return 0;
438 if (!p->FindNewAxisLimits(axis, x, xmin, xmax))
442 T* hold = (
T*)p->IsA()->New();
444 hold->SetDirectory(0);
448 if (p->fBinSumw2.fN) hold->Sumw2();
450 Int_t nbinsx = p->fXaxis.GetNbins();
451 Int_t nbinsy = p->fYaxis.GetNbins();
452 Int_t nbinsz = p->fZaxis.GetNbins();
458 Int_t ix, iy, iz, binx, biny, binz;
459 for (binz=1;binz<=nbinsz;binz++) {
460 bz = hold->GetZaxis()->GetBinCenter(binz);
461 iz = p->fZaxis.FindFixBin(bz);
462 for (biny=1;biny<=nbinsy;biny++) {
463 by = hold->GetYaxis()->GetBinCenter(biny);
464 iy = p->fYaxis.FindFixBin(by);
465 for (binx=1;binx<=nbinsx;binx++) {
466 bx = hold->GetXaxis()->GetBinCenter(binx);
467 ix = p->fXaxis.FindFixBin(bx);
469 Int_t sourceBin = hold->GetBin(binx,biny,binz);
470 Int_t destinationBin = p->GetBin(ix,iy,iz);
471 p->AddBinContent(destinationBin, hold->fArray[sourceBin]);
472 p->fBinEntries.fArray[destinationBin] += hold->fBinEntries.fArray[sourceBin];
473 p->fSumw2.fArray[destinationBin] += hold->fSumw2.fArray[sourceBin];
474 if (p->fBinSumw2.fN) p->fBinSumw2.fArray[destinationBin] += hold->fBinSumw2.fArray[sourceBin];
481 template <
typename T>
491 for (bin=0;bin<p->fN;bin++) {
492 p->fArray[bin] = c1*cu1[bin];
493 p->fSumw2.fArray[bin] = ac1*ac1*er1[bin];
494 p->fBinEntries.fArray[bin] = en1[bin];
498 template <
typename T>
511 if (p->fBinSumw2.fN > 0 ) p->fBinSumw2.Set(0);
515 if ( p->fBinSumw2.fN == p->fNcells) {
516 if (!p->fgDefaultSumw2)
517 Warning(
"Sumw2",
"Sum of squares of profile bin weights structure already created");
521 p->fBinSumw2.Set(p->fNcells);
524 for (
Int_t bin=0; bin<p->fNcells; bin++) {
525 p->fBinSumw2.fArray[bin] = p->fBinEntries.fArray[bin];
529 template <
typename T>
539 if (ax[0] ==
'y' || ax[0] ==
'Y') axis = p->GetYaxis();
540 if (ax[0] ==
'z' || ax[0] ==
'Z') axis = p->GetZaxis();
542 Error(
"TProfileHelper::LabelsDeflate",
"Invalid axis option %s",ax);
553 while ((obj =
next())) {
554 Int_t ibin = obj->GetUniqueID();
555 if (ibin > nbins) nbins = ibin;
557 if (nbins < 1) nbins = 1;
558 T *hold = (
T*)p->IsA()->New();;
567 axis->
Set(nbins,xmin,xmax);
568 p->SetBinsLength(-1);
569 p->fBinEntries.Set(p->fN);
570 p->fSumw2.Set(p->fN);
571 if (p->fBinSumw2.fN) p->fBinSumw2.Set(p->fN);
577 Int_t bin,binx,biny,binz;
578 for (bin =0; bin < hold->fN; ++bin)
580 hold->GetBinXYZ(bin, binx, biny, binz);
581 Int_t ibin = p->GetBin(binx, biny, binz);
582 p->fArray[ibin] += hold->fArray[bin];
583 p->fBinEntries.fArray[ibin] += hold->fBinEntries.fArray[bin];
584 p->fSumw2.fArray[ibin] += hold->fSumw2.fArray[bin];
585 if (p->fBinSumw2.fN) p->fBinSumw2.fArray[ibin] += hold->fBinSumw2.fArray[bin];
591 template <
typename T>
601 if (ax[0] ==
'y' || ax[0] ==
'Y') axis = p->GetYaxis();
602 T *hold = (
T*)p->IsA()->New();;
612 xmax = xmin + 2*(xmax-
xmin);
616 axis->
Set(nbins,xmin,xmax);
618 p->SetBinsLength(-1);
619 Int_t ncells = p->fN;
620 p->fBinEntries.Set(ncells);
621 p->fSumw2.Set(ncells);
622 if (p->fBinSumw2.fN) p->fBinSumw2.Set(ncells);
625 for (
Int_t ibin =0; ibin < p->fN; ibin++)
627 Int_t binx, biny, binz;
628 p->GetBinXYZ(ibin, binx, biny, binz);
629 Int_t bin = hold->GetBin(binx, biny, binz);
631 if (p->IsBinUnderflow(ibin) || p->IsBinOverflow(ibin)) {
632 p->UpdateBinContent(ibin, 0.0);
633 p->fBinEntries.fArray[ibin] = 0.0;
634 p->fSumw2.fArray[ibin] = 0.0;
635 if (p->fBinSumw2.fN) p->fBinSumw2.fArray[ibin] = 0.0;
637 p->fArray[ibin] = hold->fArray[bin];
638 p->fBinEntries.fArray[ibin] = hold->fBinEntries.fArray[bin];
639 p->fSumw2.fArray[ibin] = hold->fSumw2.fArray[bin];
640 if (p->fBinSumw2.fN) p->fBinSumw2.fArray[ibin] = hold->fBinSumw2.fArray[bin];
646 template <
typename T>
657 template <
typename T>
662 if (p->fBuffer) p->BufferEmpty();
664 if (bin < 0 || bin >= p->fNcells)
return 0;
666 Double_t sum = p->fBinEntries.fArray[bin];
667 Double_t err2 = p->fSumw2.fArray[bin];
668 Double_t neff = p->GetBinEffectiveEntries(bin);
669 if (sum == 0)
return 0;
689 if (err2 != 0 && neff < 5) test = eprim2*sum/err2;
691 if (p->fgApproximate && (test < 1.e-4 || eprim2 <= 0)) {
697 if (p->GetDimension() == 2) index = 7;
698 if (p->GetDimension() == 3) index = 11;
721 template <
typename T>
727 if (bin < 0 || bin >= p->fNcells)
return;
728 p->fBinEntries.fArray[bin] =
w;
729 if (p->fBinSumw2.fN) p->fBinSumw2.fArray[bin] =
w;
static void SetBinEntries(T *p, Int_t bin, Double_t w)
static void Scale(T *p, Double_t c1, Option_t *option)
virtual void SetLimits(Double_t xmin, Double_t xmax)
ClassImp(TSeqCollection) Int_t TSeqCollection TIter next(this)
Return index of object in collection.
static Double_t GetBinError(T *p, Int_t bin)
virtual void GetStats(Double_t *stats) const
fill the array stats from the contents of this histogram The array stats must be correctly dimensione...
static void SetErrorOption(T *p, Option_t *opt)
static void LabelsInflate(T *p, Option_t *)
virtual void AddFirst(TObject *obj)
Add object at the beginning of the list.
void ToLower()
Change string to lower-case.
static void BuildArray(T *p)
static Bool_t GetDefaultSumw2()
return kTRUE if TH1::Sumw2 must be called when creating new histograms.
virtual Bool_t IsEmpty() const
virtual void AddAll(const TCollection *col)
static double p2(double t, double a, double b, double c)
static void LabelsDeflate(T *p, Option_t *)
if(pyself &&pyself!=Py_None)
static T * ExtendAxis(T *p, Double_t x, TAxis *axis)
void Error(const char *location, const char *msgfmt,...)
virtual void SetRange(Int_t first=0, Int_t last=0)
Set the viewing range for the axis from bin first to last.
Class to manage histogram axis.
virtual Double_t GetBinUpEdge(Int_t bin) const
Return up edge of bin.
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Collection abstract base class.
Bool_t TestBit(UInt_t f) const
static double p1(double t, double a, double b)
void Warning(const char *location, const char *msgfmt,...)
static Bool_t Add(T *p, const TH1 *h1, const TH1 *h2, Double_t c1, Double_t c2=1)
virtual void SetDirectory(TDirectory *dir)
Change the tree's directory.
static Double_t GetBinEffectiveEntries(T *p, Int_t bin)
Mother of all ROOT objects.
THashList * GetLabels() const
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
virtual Int_t BufferEmpty(Int_t action=0)
Fill histogram with all entries in the buffer.
static void Sumw2(T *p, Bool_t flag)
Double_t Sqrt(Double_t x)
virtual void Set(Int_t nbins, Double_t xmin, Double_t xmax)
Initialize axis with fix bins.
static Long64_t Merge(T *p, TCollection *list)