// general
#include "TH1.h"
#include "TH2F.h"
#include "TFile.h"
#include "TTree.h"
#include "TCanvas.h"
#include "TSystem.h"
#include "TF1.h"
#include "TKey.h"
#include "TH1F.h"
#include "TStyle.h"
#include "TProfile.h"
#include "TLegend.h"
#include "TLine.h"
#include "TArrow.h"
#include "TLatex.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TFormula.h"
#include "TAxis.h"
#include "TRandom3.h"
#include "TMath.h"
#include "THStack.h"

#include <cmath>
#include <iostream>
#include <stdio.h>
#include <fstream>
#include <sstream>
#include <vector>
#include <algorithm>


void ReBin(TH1F * &h, TAxis * &ax,int entries_);
void ReBinTTbar(TH1F * &h);


void ReBin2()
{
  using namespace std;
  
  gStyle->SetOptStat(0);

  TString fileName = "filein.root";
  TFile *target = new TFile("fileout_Bin.root","recreate" );
  TFile *file = TFile::Open(fileName);

  if (!file) abort();

  TFile *first_source = TFile::Open(fileName);
  TString path ((char *) strstr (target->GetPath (), ":"));
  path.Remove (0, 2);

  //// Find in the file the TH1 = MVA_MLP_TT , make the custom binning and save a TAxis
  TAxis *axis ;
  TH1F * ht = (TH1F*) file->Get("MVA_MLP_TT");
  int entries = ht->GetEntries();
  
 /////////make the binning
  ReBinTTbar(ht);
  cout<<" OK, done with the TTJets......"<<endl; 

  ///create a new TAxis based on the previous binning
  TFile *taxisFile;
  taxisFile = new TFile ("TT_axis.root", "recreate");
  axis = ht->GetXaxis();
  taxisFile->cd ();
  axis->Write();
  taxisFile->Write ();
  ht->Write();
  cout<<"  before  "<<axis->GetNbins()<<"  "<<ht->GetNbinsX()<<endl;
  taxisFile->Close ();
  target->cd();
  ht->Write();
  
  //Now, you should have a file TT_axis.root with the TAxis and the TH1 of the Binned template coming from the the MVA_MLP_TT

  ///////Loop in the source for the rest of TH1 
  first_source->cd (path);
   
  TDirectory *current_sourcedir = gDirectory;
  //gain time, do not add the objects in the list in memory
  Bool_t status = TH1::AddDirectoryStatus ();
  TH1::AddDirectory (kFALSE);

  ///////// loop over all keys in this directory
  TChain *globChain = 0;
  TIter nextkey (current_sourcedir->GetListOfKeys ());
  TKey *key, *oldkey = 0;
  while ((key = (TKey *) nextkey ())) {

    //keep only the highest cycle number for each key
    if (oldkey && !strcmp (oldkey->GetName (), key->GetName ()))
      continue;

    first_source->cd (path);
    TObject *obj = key->ReadObj ();
   
    TLegend *legend_c1 = new TLegend (0.35, 0.90, 0.80, 0.70);
    

    if (obj->IsA ()->InheritsFrom ("TH1F")) {
      TH1F *h1 = (TH1F *) obj;
      string name=h1->GetName();

      /////////If TH1 is not MVA_MLP_TT rebin it by using the TAxis obtained from the previous step
      if (name!="MVA_MLP_TT"){
        ReBin(h1,axis,entries);
	target->cd();
        h1->Write();
      }

    }
  }
  target->Write();
  target->Close();


}
//////////////// This will make a binning  of the TT_MVA Classification- Default range is -05 to 2.
///We will first merge all bins below 0 and then scan the bins where x>0 - If one bin
// is found to have less than  10 entries all bins from this to the end will be merged to one bin.

  void ReBinTTbar(TH1F * &h){
  float binContent[100];
  binContent{100}=0;

  TH1F *histo;
  histo=h;

  int count_new_bins=0;
  int last_big_bin=0;
  int merge_up_to_this_bin=0;
  for (int nbbs=0;nbbs<int(h->GetNbinsX());nbbs++){
  }
  for  (unsigned int nb = 0;nb<int(h->GetNbinsX());nb++){
    
    ///////First find the first bin above x>0 which has more than 10 entries - merge_up_to_this_bin refers to this exact bin
  
     float  small_bins=0;	
    if ( h->GetBinContent(nb)>10 && h->GetBinContent(nb-1)<10 && h->GetBinLowEdge(nb-1)<0){
      merge_up_to_this_bin=nb-1;
    }

    //////// the last_big_bin counts for the last bin with more than 10 entries      	
    if (h->GetBinContent(nb)<10 && h->GetBinLowEdge(nb)>0){
      last_big_bin=nb-1;
      break;
    }
  }

  ///////merge the bin contents of the x<0 area
  for (int k=0;k<merge_up_to_this_bin+1;k++){
    small_bins+=h->GetBinContent(k);
  }

  ////////assign the new bin content
  for  (unsigned int n_big = merge_up_to_this_bin+1;n_big<last_big_bin+1;n_big++){

    count_new_bins++;
    binContent[count_new_bins] = float(h->GetBinContent(n_big));
  }

  binContent[1]+=small_bins;


  for  (unsigned int n_small = last_big_bin+1;n_small<int(h->GetNbinsX());n_small++){
    binContent[count_new_bins+1]+= float(h->GetBinContent(n_small));
  }


  ///now make a new TH1 which will hold  the new "binned" distribution
  TH1F *histo_binned;
  histo_binned = new TH1F("","",count_new_bins-1,0,2);

  /////get the same name,title from the input TH1
  histo_binned->SetName(h->GetName());
  histo_binned->SetTitle(h->GetTitle());
  
   ////now, make the new binning and return the new histogram

  histo_binned->SetBinContent(merge_up_to_this_bin+1,small_bins);
  for (unsigned int nbb=1;nbb<count_new_bins+1;nbb++){
    histo_binned->SetBinContent(nbb,binContent[nbb]);
    cout<<" Filling the new histo  "<<nbb<<" content "<<binContent[nbb]<<"  low edge   "<<h->GetBinLowEdge(nbb)<<endl;
  }
  histo_binned->Sumw2();
  h=histo_binned;
}

///////input is the initial histo, the TAxis and the netries from the template histo
void ReBin(TH1F * &h, TAxis * &axis,int entries){

  TH1F *histo;
  histo=h;

  
  cout<<"  In Rebin ==> NBins = "<<axis->GetNbins()<< " Min Axis = " << axis->GetXmin()
   << " Max Axis = " << axis->GetXmax() << " Hist Name = " << h->GetTitle()<< endl;
  /////// make a new histo, with TAxis 
//  TH1F *histo_binned = new TH1F(h->GetName(),h->GetTitle(),axis->GetNbins(), axis->GetXbins()->fArray );
   TH1F *histo_binned = new TH1F(h->GetName(),h->GetTitle(),axis->GetNbins(), axis->GetXmin(), axis->GetXmax() );

  vector <float > values;
  values.clear();

  Double_t *xbins = new Double_t[entries+1];
  Int_t k=0;

  ///get the values by combining the TAxis (for the N bins)  with  the enties of the input histo by using the BinLowEdge - dump the values in the values vector

  for (Int_t i=1;i<=axis->GetNbins();i++) {

    Int_t y = (Int_t)h->GetBinContent(i);
    if (y <=0) continue;
    Double_t dx = axis->GetBinWidth(i)/y;
    Double_t xmin = axis->GetBinLowEdge(i);
    for (Int_t j=0;j<y;j++) {
      xbins[k] = xmin +j*dx;
      values.push_back(xbins[k]);
      //printf("xbins[%d]= %g\n",k,xbins[k]);
      k++;
    }

  }
  ///fill the new binned histo
  for (int k=0;k<values.size();k++){
    histo_binned->Fill(values[k]);
  }

  ///return the binned histo
  h=histo_binned;  

}