Rebin a TH1 plot with minEventsPerBin

Hi

Following my previous post
/viewtopic.php?f=10&t=10747

I am trying to re-process an existing TH1 and re-bin with a certain #entries/bin . What I do is this

 
 TH1F hnew;
  hnew = fin->Get("OldTH1");

int minEventsPerBin = 50;
  int nentries = hnew.GetEntries();
  for (int v=0; v<nentries;v++){

    if (v%minEventsPerBin==0){
 BinCenterObs.push_back(hnew.GetBinCenter(v));
 BinContentObs.push_back(hnew.GetBinContent(v)); 
     } 
   }

    int ssize =  BinCenterObs.size();
    double valsCenter[ssize], valsContent[ssize];

    for (int b=0;b< BinCenterObs.size();b++){
      valsCenter[b] = BinCenterObs[b];
  
  }
  
    for (int b=0;b< BinContentObs.size();b++){
  
      valsContent[b]=BinContentObs[b];
     
  }
 

 int nb = nentries/minEventsPerBin;
   TH1F * hnew2 = new TH1F ("hnew2", "AlljetsPt", nb, valsCenter);
  
   Double_t value;
 for(int b = 0; b <  BinCenterObs.size(); ++b) {
   value = valsCenter[b];//htemp->GetBinCenter( b );
   Double_t weight = valsContent[b];//htemp->GetBinContent( b );
         
    for(int w = 0; w < weight; ++w) {
      hnew2->Fill(w);
            }
      }
  
  hnew2->Write(); 

I know, I am doing something wrong at the end, as I cannot get the proper “final” histo…Can you possibly help out ??

Thanks in advance

Alternatevely I tried

[code]
TH1F hnew;
hnew = fin->Get(“OldTH1”);

int minEventsPerBin = 50;
int nentries = hnew.GetEntries();
for (int v=0; v<nentries;v++){

if (v%minEventsPerBin==0){

BinCenterObs.push_back(hnew.GetBinCenter(v));
BinContentObs.push_back(hnew.GetBinContent(v));
}
}

int ssize =  BinCenterObs.size();
double valsCenter[ssize], valsContent[ssize];

for (int b=0;b< BinCenterObs.size();b++){
  valsCenter[b] = BinCenterObs[b];

}

for (int b=0;b< BinContentObs.size();b++){

  valsContent[b]=BinContentObs[b];

}

for (int b=0;b< LowEdgeObs.size();b++){

  valsLowEdge[b]=LowEdgeObs[b];

}

int nb = nentries/50;

cout<<" nb “<<nb<<” valsLowEdge "<<LowEdgeObs.size()<<endl;
hnew.Rebin(nb,“hnew3”,valsLowEdge);
hnew.SetTitle(“Bins”);hnew.SetName(“Bins”);
hnew.Write();

   }

}[/code]

but, what I get, is not what I expect, ie a “flat” distribution and even the axis are wrong…

Any help is most welcome and appreciated




Hi,

At the very least this loop is weird/wrong:[code] int nentries = hnew.GetEntries();
for (int v=0; v<nentries;v++){

if (v%minEventsPerBin==0){

BinCenterObs.push_back(hnew.GetBinCenter(v));
BinContentObs.push_back(hnew.GetBinContent(v));
}
}[/code]as GetBinCenter does not expect an ‘entry index’ but a bin number …

Philippe.

Hi Philippe

So, could you possibly indicate what the loop should be like in order work as I need it ?? I know, that it starts to be annoying with all continuous posts from my side (also similar topic How to make custom binning ), but my “ultimate” goal as I already said, is to be able to bin a histo with custom#events/bin and maybe pass this to a TAxis for further processing / next plots of the same variable…and as the days pass by and I cannot get into this, I am starting to be really frustrated… So, if by any change you know to implement this (or even correct my existing scripts), in order to have this I would be more obliged

Thanks for your tolerance and time!

Hi,

I am probably not understand your goal, but it seems to me that you can not create this histogram only from the histogrammed version of the data (because there is missing information). You might vaguely approximate the result by calculating (pseudo code): double entries_per_bins = h->GetEntries() / h->GetNbins(); Double_t *xq = new Double_t[nbins+1]; Double_t *content = new Double_t[nbins+1]; xq[0] = h->GetXaxis()->GetXmin(); content[0] = h->GetBinContent(0); // underflow xq[nbins] = h->GetXaxis()->GetXmax(); content[nbins] = h->GetBinContent(nbins); // overflow double entries_so_far = 0; int new_i = 1; for(int i = 1; i < nbins; ++i) { entries_so_far += h->GetBinContent(1); if (entries_so_far >= entries_per_bins || (i==(nbins-1))); xq[ new_i ] = h->GetXaxis()->GetBinUpEdge( i ); content[new_i] = entries_so_far; entries_so_far = 0; } } etc ...

Philippe.

Hi Philippe

the truth is that I have all the original info also in TTree…So, could this solve the problem (how ??)

Thanks!

Hi,

as I said
a) it’s a very tricky problem,
b) I will post the code once I wrote it - probably tomorrow.

In the meantime I have learned how to solve it (n^2 in memory, n^3 in CPU time with n the number of bins of the original histogram - if anyone can find a better algorithm let me know!)

Cheers, Axel.

Ok Axel

I will eagerly wait for your script, and then maybe I could help on improving the code…

Thanks a lot!

Hi,

here you go! As I said: it’s complex :slight_smile: And thanks to Jakob Blomer for coming up with this algorithm! Run as
.x mergebins.C+
The function mergebins() creates an example histogram and then lets the algorithm do its job. The result is better for more bins (but the CPU time goes like nbins^3!) and the fewer target bins / the more target bin content you require. Please read the NOTE at the top of the file!

Enjoy!

Cheers, Axel.


mergebins.C (4.49 KB)

Hi Axel

Thanks a lot for this --it is really helpful…Nevertheless, can you directly point me what (or where :wink: I need to modify it in order to preserve the integral (ie making a binning according to a given #events/bin instead of #bins in the final histo)

Thanks again for all your help

Hi,

simple:

-      vBinContent.push_back(Content(idx));
+      vBinContent.push_back(Content(idx)/span);

Afterwards, bin width*bin content (i.e. the bin integral) is almost equal for all bins.

Cheers, Axel.

Hi Axel

Thanks for your prompt reply… Nevertheless I think as this loop

for (UInt_t i = 0, n = vbinx.size(); i < n; ++i) { hVarBins->SetBinContent(i +1, vbinc[i]); }

gives a plot where the 1rst bin is wrong as it look to include events, I think it should be like

for (UInt_t i = 1, n = vbinx.size(); i < n; ++i) { hVarBins->SetBinContent(i +1, vbinc[i]); }

in order to give the correct one… Correct me if I am wrong

Cheers




Hi,

sorry, I don’t understand why you think it’s wrong. I also don’t understand why you want to skip setting the content of the first bin.

Cheers, Axel.

Hi

Yes, on a second thought you are completely right… :wink: Also, I managed to write a script in order to have a certain binning according to #events per bin, ie the final plot should be a flat distrubution… Nevertheless, I am doing that while filling the histo (from the “original” data )and not from an already created TH1 object

All best

[quote=“Alkass”]Hi

Yes, on a second thought you are completely right… :wink: Also, I managed to write a script in order to have a certain binning according to #events per bin, ie the final plot should be a flat distrubution… Nevertheless, I am doing that while filling the histo (from the “original” data )and not from an already created TH1 object

All best[/quote]

Dear Alkass!

Could you please post your modified script that does the optimization during filling?

Thanks,

Sebastian