RandHisto.cxx
Go to the documentation of this file.
1 #include "RandHisto.h"
2 
3 #include "TMath.h"
4 
6 {}
7 
8 void RandHisto::GetRandom(TH1* histo, double& x, double& y, double& z)
9 {
10  int dim = histo->GetDimension();
11  int nbins = histo->GetNbinsX();
12  int nbinsy = 0;
13  int nbinsz = 0;
14  int nbinsx = histo->GetNbinsX();
15  if (dim > 1){
16  nbinsy = histo->GetNbinsY();
17  nbins *= nbinsy;
18  }
19  if (dim > 2){
20  nbinsz = histo->GetNbinsZ();
21  nbins *= nbinsz;
22  }
23 
24  int ibinX(0), ibinY(0), ibinZ(0);
25  //GetIntegral is normalized to 1 when it computes the bin integral array
26  int ibinG = TMath::BinarySearch(nbins, histo->GetIntegral(), fFlat.fire());
27  //histo->GetBinXYZ(ibinG, ibinX, ibinY, ibinZ);
28  //NOTE the GetIntegral computes its bin integral array without overflow/underflow bins (indexing starting at 0). So GetBinXYZ cannot be used to convert
29  //the global bin index. At best it will only be off by 1 for 1d histograms. At worst it will be very wrong for 2d and 3d histograms.
30  //To avoid this compute the indexing by hand.
31  if (dim == 1) ibinX = ibinG + 1;
32  if (dim == 2){
33  ibinY = (ibinG / nbinsx) + 1;
34  ibinX = (ibinG % nbinsx) + 1;
35  }
36  if (dim == 3){
37  ibinZ = ibinG / (nbinsx*nbinsy);
38  ibinY = ((ibinG - (nbinsx*nbinsy*ibinZ)) / nbinsx) + 1;
39  ibinX = ((ibinG - (nbinsx*nbinsy*ibinZ)) % nbinsx) + 1;
40  ibinZ++;
41  }
42 
43 
44  //x = histo->GetXaxis()->GetBinCenter(ibinX);
45  double xCenter = histo->GetXaxis()->GetBinCenter(ibinX);
46  double xWidth = histo->GetXaxis()->GetBinWidth(ibinX);
47  x = (xCenter-0.5*xWidth) + xWidth*fFlat.fire();
48 
49  if (dim > 1)
50  {
51  //y = histo->GetYaxis()->GetBinCenter(ibinY);
52  double yCenter = histo->GetYaxis()->GetBinCenter(ibinY);
53  double yWidth = histo->GetYaxis()->GetBinWidth(ibinY);
54  y = (yCenter-0.5*yWidth) + yWidth*fFlat.fire();
55  }
56 
57  if (dim > 2)
58  {
59  //z = histo->GetZaxis()->GetBinCenter(ibinZ);
60  double zCenter = histo->GetZaxis()->GetBinCenter(ibinZ);
61  double zWidth = histo->GetZaxis()->GetBinWidth(ibinZ);
62  z = (zCenter-0.5*zWidth) + zWidth*fFlat.fire();
63  }
64 
65  return;
66 }
67 
68 
69 
const int nbins
Definition: cellShifts.C:15
TH2D * histo
CLHEP::RandFlat fFlat
Definition: RandHisto.h:17
z
Definition: test.py:28
Int_t nbinsx
Definition: plot.C:23
void GetRandom(TH1 *histo, double &x, double &y, double &z)
Definition: RandHisto.cxx:8
RandHisto(CLHEP::HepRandomEngine &engine)
Definition: RandHisto.cxx:5