Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
hist013_TH1_rebin.C
Go to the documentation of this file.
1/// \file
2/// \ingroup tutorial_hist
3/// \notebook -js
4/// Rebin a variable bin-width histogram.
5///
6/// This tutorial illustrates how to:
7/// - create a variable bin-width histogram with a binning such
8/// that the population per bin is about the same.
9/// - rebin a variable bin-width histogram into another one.
10///
11/// \macro_image
12/// \macro_code
13///
14/// \date July 2016
15/// \author Rene Brun
16
17#include "TH1.h"
18#include "TCanvas.h"
20{
21 // create a fix bin histogram
22 TH1F *h = new TH1F("h", "test rebin", 100, -3, 3);
23 Int_t nentries = 1000;
24 h->FillRandom("gaus", nentries);
25 Double_t xbins[1001];
26 Int_t k = 0;
27 TAxis *axis = h->GetXaxis();
28 for (Int_t i = 1; i <= 100; i++) {
29 Int_t y = (Int_t)h->GetBinContent(i);
30 if (y <= 0)
31 continue;
32 Double_t dx = axis->GetBinWidth(i) / y;
33 Double_t xmin = axis->GetBinLowEdge(i);
34 for (Int_t j = 0; j < y; j++) {
35 xbins[k] = xmin + j * dx;
36 k++;
37 }
38 }
39 xbins[k] = axis->GetXmax();
40 // create a variable bin-width histogram out of fix bin histogram
41 // new rebinned histogram should have about 10 entries per bin
42 TH1F *hnew = new TH1F("hnew", "rebinned", k, xbins);
43 hnew->FillRandom("gaus", 10 * nentries);
44
45 // rebin hnew keeping only 50% of the bins
46 Double_t xbins2[501];
47 Int_t kk = 0;
48 for (Int_t j = 0; j < k; j += 2) {
49 xbins2[kk] = xbins[j];
50 kk++;
51 }
52 xbins2[kk] = xbins[k];
53 TH1F *hnew2 = (TH1F *)hnew->Rebin(kk, "hnew2", xbins2);
54
55 // draw the 3 histograms
56 TCanvas *c1 = new TCanvas("c1", "c1", 800, 1000);
57 c1->Divide(1, 3);
58 c1->cd(1);
59 h->Draw();
60 c1->cd(2);
61 hnew->Draw();
62 c1->cd(3);
63 hnew2->Draw();
64}
#define h(i)
Definition RSha256.hxx:106
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
float xmin
int nentries
Class to manage histogram axis.
Definition TAxis.h:32
Double_t GetXmax() const
Definition TAxis.h:142
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition TAxis.cxx:513
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition TAxis.cxx:537
The Canvas class.
Definition TCanvas.h:23
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:645
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
const double xbins[xbins_n]