Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
rebin.C File Reference

Detailed Description

View in nbviewer Open in SWAN
Rebin a variable bin-width histogram.

This tutorial illustrates how to:

  • create a variable bin-width histogram with a binning such that the population per bin is about the same.
  • rebin a variable bin-width histogram into another one.
#include "TH1.h"
#include "TCanvas.h"
void rebin() {
//create a fix bin histogram
TH1F *h = new TH1F("h","test rebin",100,-3,3);
Int_t nentries = 1000;
h->FillRandom("gaus",nentries);
Double_t xbins[1001];
Int_t k=0;
TAxis *axis = h->GetXaxis();
for (Int_t i=1;i<=100;i++) {
Int_t y = (Int_t)h->GetBinContent(i);
if (y <=0) continue;
Double_t dx = axis->GetBinWidth(i)/y;
for (Int_t j=0;j<y;j++) {
xbins[k] = xmin +j*dx;
k++;
}
}
xbins[k] = axis->GetXmax();
//create a variable bin-width histogram out of fix bin histogram
//new rebinned histogram should have about 10 entries per bin
TH1F *hnew = new TH1F("hnew","rebinned",k,xbins);
hnew->FillRandom("gaus",10*nentries);
//rebin hnew keeping only 50% of the bins
Double_t xbins2[501];
Int_t kk=0;
for (Int_t j=0;j<k;j+=2) {
xbins2[kk] = xbins[j];
kk++;
}
xbins2[kk] = xbins[k];
TH1F *hnew2 = (TH1F*)hnew->Rebin(kk,"hnew2",xbins2);
//draw the 3 histograms
TCanvas *c1 = new TCanvas("c1","c1",800,1000);
c1->Divide(1,3);
c1->cd(1);
h->Draw();
c1->cd(2);
hnew->Draw();
c1->cd(3);
hnew2->Draw();
}
#define h(i)
Definition RSha256.hxx:106
int Int_t
Definition RtypesCore.h:45
double Double_t
Definition RtypesCore.h:59
float xmin
int nentries
Class to manage histogram axis.
Definition TAxis.h:31
Double_t GetXmax() const
Definition TAxis.h:140
virtual Double_t GetBinLowEdge(Int_t bin) const
Return low edge of bin.
Definition TAxis.cxx:518
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width.
Definition TAxis.cxx:540
The Canvas class.
Definition TCanvas.h:23
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:621
virtual void FillRandom(const char *fname, Int_t ntimes=5000, TRandom *rng=nullptr)
Fill histogram following distribution in function fname.
Definition TH1.cxx:3519
void Draw(Option_t *option="") override
Draw this histogram with options.
Definition TH1.cxx:3066
virtual TH1 * Rebin(Int_t ngroup=2, const char *newname="", const Double_t *xbins=nullptr)
Rebin this histogram.
Definition TH1.cxx:6243
Double_t y[n]
Definition legend1.C:17
return c1
Definition legend1.C:41
const double xbins[xbins_n]
Author
Rene Brun

Definition in file rebin.C.