Skip to content

Commit 655cdc7

Browse files
ferdymercuryferdymercury
committed
[hist] do not look at bin center if one edge is infinite
Fixes #20174 [hist] specifically dealing with edge cases
1 parent 34f9d42 commit 655cdc7

File tree

6 files changed

+117
-6
lines changed

6 files changed

+117
-6
lines changed

‎hist/hist/src/TAxis.cxx‎

Lines changed: 23 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -477,15 +477,38 @@ Int_t TAxis::GetLast() const
477477

478478
////////////////////////////////////////////////////////////////////////////////
479479
/// Return center of bin
480+
///
481+
/// If fXmin is -inf, it will return -inf even if fXmax is finite
482+
/// If fXmax is +inf, it will return +inf even if fXmin is finite
480483

481484
Double_t TAxis::GetBinCenter(Int_t bin) const
482485
{
483486
Double_t binwidth;
484487
if (!fXbins.fN || bin<1 || bin>fNbins) {
485488
binwidth = (fXmax - fXmin) / Double_t(fNbins);
489+
if (std::isinf(binwidth)) {
490+
if (std::isfinite(fXmin))
491+
return fXmax;
492+
else if (std::isfinite(fXmax))
493+
return fXmin;
494+
else if (fXmin == -std::numeric_limits<double>::infinity() && fXmax == std::numeric_limits<double>::infinity())
495+
return 0.;
496+
else
497+
return std::numeric_limits<double>::quiet_NaN();
498+
}
486499
return fXmin + (bin - 0.5) * binwidth;
487500
} else {
488501
binwidth = fXbins.fArray[bin] - fXbins.fArray[bin-1];
502+
if (std::isinf(binwidth)) {
503+
if (std::isfinite(fXbins.fArray[bin-1]))
504+
return fXbins.fArray[bin];
505+
else if (std::isfinite(fXbins.fArray[bin]))
506+
return fXbins.fArray[bin-1];
507+
else if (fXbins.fArray[bin-1] == -std::numeric_limits<double>::infinity() && fXbins.fArray[bin] == std::numeric_limits<double>::infinity())
508+
return 0.;
509+
else
510+
return std::numeric_limits<double>::quiet_NaN();
511+
}
489512
return fXbins.fArray[bin-1] + 0.5*binwidth;
490513
}
491514
}

‎hist/hist/src/TH1Merger.cxx‎

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -714,6 +714,17 @@ Bool_t TH1Merger::AutoP2Merge()
714714

715715
Double_t xu = hist->GetBinCenter(ibin);
716716
Int_t jbin = fH0->FindBin(xu);
717+
if (jbin == fH0->GetXaxis()->GetNbins() + 1 || jbin == 0) {
718+
auto eps = 1e-12;
719+
// if upper/lower edge is +/-infinite, the bin center is +/-infinite if the other edge is finite,
720+
// so FindBin is in overflow/underflow
721+
// Check close to the lower or upper edges instead of the bin center in these cases
722+
if (std::isinf(hist->GetXaxis()->GetBinUpEdge(ibin))) {
723+
jbin = fH0->GetXaxis()->FindBin(hist->GetXaxis()->GetBinLowEdge(ibin) + eps);
724+
} else if (std::isinf(hist->GetXaxis()->GetBinLowEdge(ibin))) {
725+
jbin = fH0->GetXaxis()->FindBin(hist->GetXaxis()->GetBinUpEdge(ibin) - eps);
726+
}
727+
}
717728

718729
fH0->AddBinContent(jbin, cu);
719730
if (fH0->fSumw2.fN)

‎hist/hist/src/TH1Merger.h‎

Lines changed: 15 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -16,6 +16,7 @@
1616
#include "TProfile2D.h"
1717
#include "TProfile3D.h"
1818
#include "TList.h"
19+
#include <cmath>
1920

2021
class TH1Merger {
2122

@@ -41,13 +42,25 @@ class TH1Merger {
4142
return outAxis.FindFixBin(inAxis.GetBinCenter(ibin));
4243
}
4344

44-
// find bin number estending the axis
45+
/// find bin number extending the axis
4546
static Int_t FindBinNumber(Int_t ibin, const TAxis & inAxis, TAxis & outAxis) {
4647
// should I ceck in case of underflow/overflow if underflow/overflow values of input axis
4748
// outside output axis ?
4849
if (ibin == 0 ) return 0; // return underflow
4950
if (ibin == inAxis.GetNbins()+1 ) return outAxis.GetNbins()+1; // return overflow
50-
return outAxis.FindBin(inAxis.GetBinCenter(ibin));
51+
auto binOut = outAxis.FindBin(inAxis.GetBinCenter(ibin));
52+
if (binOut == outAxis.GetNbins() + 1 || binOut == 0) {
53+
auto eps = 1e-12;
54+
// if upper/lower edge is +/-infinite, the bin center is +/-infinite if the other edge is finite,
55+
// so FindBin is in overflow/underflow
56+
// Check close to the lower or upper edges instead of the bin center in these cases
57+
if (std::isinf(inAxis.GetBinUpEdge(ibin))) {
58+
binOut = outAxis.FindBin(inAxis.GetBinLowEdge(ibin) + eps);
59+
} else if (std::isinf(inAxis.GetBinLowEdge(ibin))) {
60+
binOut = outAxis.FindBin(inAxis.GetBinUpEdge(ibin) - eps);
61+
}
62+
}
63+
return binOut;
5164
}
5265

5366
// Function to find if axis label list has duplicates

‎hist/hist/src/TH2.cxx‎

Lines changed: 12 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2293,7 +2293,18 @@ TH1D *TH2::DoProjection(bool onX, const char *name, Int_t firstbin, Int_t lastbi
22932293
}
22942294
}
22952295
// find corresponding bin number in h1 for outbin
2296-
Int_t binOut = h1->GetXaxis()->FindBin( outAxis->GetBinCenter(outbin) );
2296+
Int_t binOut = h1->GetXaxis()->FindBin(outAxis->GetBinCenter(outbin));
2297+
if (binOut == h1->GetXaxis()->GetNbins() + 1 || binOut == 0) {
2298+
auto eps = 1e-12;
2299+
// if upper/lower edge is +/-infinite, the bin center is +/-infinite if the other edge is finite,
2300+
// so FindBin is in overflow/underflow
2301+
// Check close to the lower or upper edges instead of the bin center in these cases
2302+
if (std::isinf(outAxis->GetBinUpEdge(outbin))) {
2303+
binOut = h1->GetXaxis()->FindBin(outAxis->GetBinLowEdge(outbin) + eps);
2304+
} else if (std::isinf(outAxis->GetBinLowEdge(outbin))) {
2305+
binOut = h1->GetXaxis()->FindBin(outAxis->GetBinUpEdge(outbin) - eps);
2306+
}
2307+
}
22972308
h1->SetBinContent(binOut ,cont);
22982309
if (computeErrors) h1->SetBinError(binOut,TMath::Sqrt(err2));
22992310
// sum all content

‎hist/hist/src/TH3.cxx‎

Lines changed: 35 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -2058,14 +2058,25 @@ TH1D *TH3::DoProject1D(const char* name, const char * title, const TAxis* projX,
20582058
}
20592059
}
20602060
}
2061-
Int_t ix = h1->FindBin( projX->GetBinCenter(ixbin) );
2061+
Int_t ix = h1->FindBin(projX->GetBinCenter(ixbin));
2062+
if (ix == h1->GetNbinsX() + 1 || ix == 0) {
2063+
auto eps = 1e-12;
2064+
// if upper/lower edge is +/-infinite, the bin center is +/-infinite if the other edge is finite,
2065+
// so FindBin is in overflow/underflow
2066+
// Check close to the lower or upper edges instead of the bin center in these cases
2067+
if (std::isinf(projX->GetBinUpEdge(ixbin))) {
2068+
ix = h1->FindBin(projX->GetBinLowEdge(ixbin) + eps);
2069+
} else if (std::isinf(projX->GetBinLowEdge(ixbin))) {
2070+
ix = h1->FindBin(projX->GetBinUpEdge(ixbin) - eps);
2071+
}
2072+
}
20622073
h1->SetBinContent(ix ,cont);
20632074
if (computeErrors) h1->SetBinError(ix, TMath::Sqrt(err2) );
20642075
// sum all content
20652076
totcont += cont;
20662077

20672078
}
2068-
if ( labels ) h1->GetXaxis()->SetCanExtend(extendable);
2079+
if (labels) h1->GetXaxis()->SetCanExtend(extendable);
20692080

20702081
// since we use a combination of fill and SetBinError we need to reset and recalculate the statistics
20712082
// for weighted histograms otherwise sumw2 will be wrong.
@@ -2267,11 +2278,32 @@ TH2D *TH3::DoProject2D(const char* name, const char * title, const TAxis* projX,
22672278

22682279
for (ixbin=0;ixbin<=1+projX->GetNbins();ixbin++) {
22692280
if ( projX->TestBit(TAxis::kAxisRange) && ( ixbin < ixmin || ixbin > ixmax )) continue;
2270-
Int_t ix = h2->GetYaxis()->FindBin( projX->GetBinCenter(ixbin) );
2281+
Int_t ix = h2->GetYaxis()->FindBin(projX->GetBinCenter(ixbin));
2282+
auto eps = 1e-12;
2283+
if (ix == h2->GetYaxis()->GetNbins() + 1 || ix == 0) {
2284+
// if upper/lower edge is +/-infinite, the bin center is +/-infinite if the other edge is finite,
2285+
// so FindBin is in overflow/underflow
2286+
// Check close to the lower or upper edges instead of the bin center in these cases
2287+
if (std::isinf(projX->GetBinUpEdge(ixbin))) {
2288+
ix = h2->GetYaxis()->FindBin(projX->GetBinLowEdge(ixbin) + eps);
2289+
} else if (std::isinf(projX->GetBinLowEdge(ixbin))) {
2290+
ix = h2->GetYaxis()->FindBin(projX->GetBinUpEdge(ixbin) - eps);
2291+
}
2292+
}
22712293

22722294
for (iybin=0;iybin<=1+projY->GetNbins();iybin++) {
22732295
if ( projY->TestBit(TAxis::kAxisRange) && ( iybin < iymin || iybin > iymax )) continue;
22742296
Int_t iy = h2->GetXaxis()->FindBin( projY->GetBinCenter(iybin) );
2297+
if (iy == h2->GetXaxis()->GetNbins() + 1 || iy == 0) {
2298+
// if upper/lower edge is +/-infinite, the bin center is +/-infinite if the other edge is finite,
2299+
// so FindBin is in overflow/underflow
2300+
// Check close to the lower or upper edges instead of the bin center in these cases
2301+
if (std::isinf(projY->GetBinUpEdge(iybin))) {
2302+
iy = h2->GetXaxis()->FindBin(projY->GetBinLowEdge(iybin) + eps);
2303+
} else if (std::isinf(projY->GetBinLowEdge(iybin))) {
2304+
iy = h2->GetXaxis()->FindBin(projY->GetBinUpEdge(iybin) - eps);
2305+
}
2306+
}
22752307

22762308
Double_t cont = 0;
22772309
Double_t err2 = 0;

‎hist/hist/test/test_TH1.cxx‎

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@
1010
#include <cstddef>
1111
#include <random>
1212
#include <vector>
13+
#include <limits>
1314

1415
#include "ROOT/TestSupport.hxx"
1516

@@ -320,3 +321,23 @@ TEST(TAxis, EqualBinEdges)
320321
{
321322
ROOT_EXPECT_ERROR(TAxis _({1, 1}), "TAxis::Set", "bins must be in increasing order");
322323
}
324+
325+
TEST(TH2, ProjectionYInfiniteUpperEdge)
326+
{
327+
Double_t xedges[] = {0., 1.};
328+
Double_t yedges[] = {1., std::numeric_limits<double>::infinity()};
329+
TH2D h("h_inf", "h_inf;X;Y", 1, xedges, 1, yedges);
330+
h.SetBinContent(1, 1, 11.);
331+
auto projY = h.ProjectionY();
332+
EXPECT_EQ(projY->GetBinContent(1), h.Integral(1, 1, 1, 1));
333+
}
334+
335+
TEST(TH2, ProjectionYInfiniteLowerEdge)
336+
{
337+
Double_t xedges[] = {0, 1.};
338+
Double_t yedges[] = {-std::numeric_limits<double>::infinity(), 2};
339+
TH2D h("h_inf", "h_inf;X;Y", 1, xedges, 1, yedges);
340+
h.SetBinContent(1, 1, 11.);
341+
auto projY = h.ProjectionY();
342+
EXPECT_EQ(projY->GetBinContent(1), h.Integral(1, 1, 1, 1));
343+
}

0 commit comments

Comments
 (0)