//C++ code

#include <stdio.h>
#include <stdlib.h>
#include <iostream>
#include <iomanip>
#include <assert.h>
#include <time.h>
#include <fstream>
#include <string>

//Root Code

#include <TH1.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TPad.h>
#include <TMath.h>
#include <TLine.h>
#include <TPaveText.h>
#include <TString.h>
#include <TTree.h>
#include <TCut.h>
#include <TGraphErrors.h>
#include <TLeaf.h>

//roofit


#ifndef __CINT__
#include "RooGlobalFunc.h"
#endif
#include <RooRealVar.h>
#include <RooDataSet.h>
#include <RooGaussian.h>
#include <RooChebychev.h>
#include <RooArgusBG.h>
#include <RooAddPdf.h>
#include <TCanvas.h>
#include <TAxis.h>
#include <RooPlot.h>
#include <RooHist.h>
#include <RooDataHist.h>
#include <RooExponential.h>
using namespace RooFit ;

//etc

using namespace std;

//supporting methods


TTree* makeTTree() ;

//main test code

void testdatahist()
{
  //Test Data hist and such.
  TTree * mytree = makeTTree();
  TCanvas *can = new TCanvas("c1");  
  TH1F* test_hist[10];

  Double_t chi2cutval = 1.0;
  Double_t massmin = 0;
  Double_t massmax = 40;
  
  TCut chi2_test_cut[10];
  TCut mass_cut = Form("mass>%lf && mass<%lf",massmin, massmax);

  //Check cuts for root histogram only;
  for(int i=0; i<10; ++i){
    chi2cutval+=0.5;
    
    chi2_test_cut[i] = Form("max(track0_chi2,track1_chi2)<%lf",chi2cutval);
    mytree->Draw("mass>>h1",chi2_test_cut[i]);
    test_hist[i] = (TH1F*)gPad->GetPrimitive("h1");
    cout<<"test_hist[i]->Integral() = "<<test_hist[i]->Integral()<<endl;

  }

  //Now use roofit

  //observables from ttree
  RooRealVar mymass("mass","mass",massmin,massmax);
  RooRealVar track0_chi2("track0_chi2","track0_chi2",0.0,20);
  RooRealVar track1_chi2("track1_chi2","track1_chi2",0.0,20);

  //get the dataset
  RooDataSet *data_unbinned = new RooDataSet("mass_example","mass example", mytree,RooArgSet(mymass,track0_chi2,track1_chi2));
  RooDataHist* data = data_unbinned->binnedClone("data");
  RooDataHist* reduced_data;
  RooDataSet *data_unbinned_reduced;
  
  for(int i=0; i<10;++i){
    data_unbinned_reduced = (RooDataSet*)data_unbinned->reduce(Cut(chi2_test_cut[i]));
    reduced_data = (RooDataHist*)data->reduce(Cut(chi2_test_cut[i]));

    cout<<"reduced_data->sumEntries() = "<<reduced_data->sumEntries()<<endl;
    cout<<"data_ubinned_reduced->sumEntries() = "<<data_unbinned_reduced->sumEntries()<<endl;
  }
  


}


TTree* makeTTree() 
{
  //modified from rf102_dataimport.C
  //Make some distributions for mass and chi2
  //Fill the tree

  TTree* tree = new TTree("tree","tree") ;
  Double_t* mass_x = new Double_t;
  Double_t* track0_chi2_x = new Double_t;
  Double_t* track1_chi2_x = new Double_t;
  
  //Double_t* px = new Double_t ;
  //Double_t* py = new Double_t ;

  //tree->Branch("x",px,"x/D") ;
  //tree->Branch("y",py,"y/D") ;

  tree->Branch("track0_chi2",track0_chi2_x,"track0_chi2/D");
  tree->Branch("track1_chi2",track1_chi2_x,"track1_chi2/D");
  tree->Branch("mass",mass_x,"mass/D");
  for (int i=0 ; i<10000 ; i++) {
    *track0_chi2_x = gRandom->Landau(1,0.5);
    *track1_chi2_x = gRandom->Landau(1,0.5);
    *mass_x = gRandom->Gaus(20,3);
    tree->Fill() ;
  }
  return tree ;
}