[hist] Improve precision of TAxis::FindFixBin / FindBin.

Due to the floating-point subtraction x - xMin, the bin index
occasionally flows over to the next bin. This is particularly
annoying when it goes into the overflow, although the coordinate is
strictly smaller than the max of the axis.
By adding a correction term that can detect this case, overflows such as
the one discussed in https://root-forum.cern.ch/t/bug-or-feature-in-ttree-draw/62862
can be avoided.

An error in the other direction is possible, too, and fixed in this
commit:
https://root-forum.cern.ch/t/floating-point-rounding-error-when-filling-the-histogram/35189

Microbenchmarking the changed lines showed them to be the same speed in
gcc, and 40% slower in clang. Both changes are by far outweighted by the
virtual call overhead, though.

This allowed for removing the cautionary note on rounding errors, added
in 1703c54, which is fixed now.

Fix #14091.

Author: hageboeck

hist/hist/src/TAxis.cxx

@@ -289,9 +289,6 @@ void TAxis::ExecuteEvent(Int_t event, Int_t px, Int_t py)
 /// whereas the overflow bin (`nbins+1`) is for any `x` greater or equal than `fXmax`,
 /// as well as for `NaN`. The first regular bin (`1`) is for any `x`
 /// greater or equal than `fXmin` and strictly smaller than `fXmin + binwidth`, and so on.
-/// @note The bin assignment equation uses doubles, thus rounding errors are
-/// expected to appear at the edges. For example: `TAxis(1, -1., 0.).FindBin(-1e-17)`
-/// returns the overflow bin (`2`) rather than the theoretically correct bin (`1`).
 
 Int_t TAxis::FindBin(Double_t x)
 {
@@ -314,7 +311,9 @@ Int_t TAxis::FindBin(Double_t x)
       return FindFixBin(x);
    } else {
       if (!fXbins.fN) {        //*-* fix bins
-         bin = 1 + int (fNbins*(x-fXmin)/(fXmax-fXmin) );
+         const double width = (fXmax-fXmin)/fNbins;
+         const int approxBin = (x-fXmin) / width; // Might be one unit too small/large due to limited precision of subtraction
+         bin = 1 + approxBin - (x < fXmin + width*approxBin) + (fXmin + width*(approxBin+1) <= x);
       } else {                  //*-* variable bin sizes
          //for (bin =1; x >= fXbins.fArray[bin]; bin++);
          bin = 1 + TMath::BinarySearch(fXbins.fN,fXbins.fArray,x);
@@ -429,7 +428,9 @@ Int_t TAxis::FindFixBin(Double_t x) const
       bin = fNbins+1;
    } else {
       if (!fXbins.fN) {        //*-* fix bins
-         bin = 1 + int (fNbins*(x-fXmin)/(fXmax-fXmin) );
+         const double width = (fXmax-fXmin)/fNbins;
+         const int approxBin = (x-fXmin) / width; // Might be one unit too small/large due to limited precision of subtraction
+         bin = 1 + approxBin - (x < fXmin + width*approxBin) + (fXmin + width*(approxBin+1) <= x);
       } else {                  //*-* variable bin sizes
          bin = 1 + TMath::BinarySearch(fXbins.fN,fXbins.fArray,x);
       }

hist/hist/test/test_TH1.cxx

@@ -6,7 +6,7 @@
 #include "TH1F.h"
 #include "THLimitsFinder.h"
 
-
+#include <random>
 #include <vector>
 
 // StatOverflows TH1