#include <iostream>
#include <fstream>
#include <cmath>
#include <vector.h>
#include <TROOT.h>
#include <TFile.h>
#include <TTree.h>
#include <TH1.h>
#include <TF1.h>
#include <TH2.h>
#include <TF2.h>
#include <TGraph.h>
#include <TMath.h>
#include <TChain.h>
#include <TMinuit.h>
#include <TVirtualFitter.h>
#include <TStopwatch.h>
#include <TLorentzVector.h>

#include"Truth.h"
#include"ZCandidates.h"
#include"HCandidates.h"
#include"Trigger.h"
#include"Muon.h"
#include"Electron.h"
#include"RawElectron.h"



//______________________________________________________________________________
int main()
{

  TFile hfile("output.root","RECREATE");

   TH1F *hMuonTrue_eta          = new TH1F("hMuonTrue_eta","MuonTrue_eta",25,0.,2.5);
   TH1F *hMuonReco_eta          = new TH1F("hMuonReco_eta","MuonReco_eta",25,0.,2.5);
   TH1F *hMuonCombined_eta      = new TH1F("hMuonCombined_eta","MuonCombined_eta",25,0.,2.5);

   TH1F *hMuonTrue_pt           = new TH1F("hMuonTrue_pt","MuonTrue_pt",55,5.,60.);
   TH1F *hMuonReco_pt           = new TH1F("hMuonReco_pt","MuonReco_pt",55,5.,60.);
   TH1F *hMuonCombined_pt       = new TH1F("hMuonCombined_pt","MuonCombined_pt",55,5.,60.);

   TH1F *hZeeee_mass            = new TH1F("hZeeee_mass","hZeeee_mass",120,5.,120.);
   TH1F *hZmmmm_mass            = new TH1F("hZmmmm_mass","hZmmmm_mass",120,5.,120.);
   TH1F *hZeemm_mass            = new TH1F("hZeemm_mass","hZeemm_mass",120,5.,120.);

   TH1F *hHeeee_mass            = new TH1F("hHeeee_mass","hHeeee_mass",100,80.,150.);
   TH1F *hHmmmm_mass            = new TH1F("hHmmmm_mass","hHmmmm_mass",100,80.,150.);
   TH1F *hHeemm_mass            = new TH1F("hHeemm_mass","hHeemm_mass",100,80.,150.);



  // instanciate chains
  TObjArray Chains;
  TChain * truth = new TChain("TruthTree"); Chains.Add(truth);
  TChain * rawelec = new TChain("RawElectronTree"); Chains.Add(rawelec);
  TChain * elec = new TChain("PreselectedElectronTree"); Chains.Add(elec);
  TChain * elecloose = new TChain("PreselectedLooseElectronTree"); Chains.Add(elecloose);
  TChain * elecmedium = new TChain("PreselectedMediumElectronTree"); Chains.Add(elecmedium);
  TChain * elecmediumcaloiso = new TChain("PreselectedMediumCALOISOElectronTree"); Chains.Add(elecmediumcaloiso);
  TChain * electight = new TChain("PreselectedTightElectronTree"); Chains.Add(electight);
  TChain * ztree = new TChain("ZTree"); Chains.Add(ztree);
  TChain * muon = new TChain("PreselectedMuonTree"); Chains.Add(muon);
//  TChain * heeee = new TChain("HeeeeCandidatesTree"); Chains.Add(heeee);
//  TChain * hmmmm = new TChain("HmmmmCandidatesTree"); Chains.Add(hmmmm);
//  TChain * heemm = new TChain("HeemmCandidatesTree"); Chains.Add(heemm);
  TChain * selelec = new TChain("SelectedMuonTree"); Chains.Add(selelec);
  TChain * selmuon = new TChain("SelectedMuonTree"); Chains.Add(selmuon);

  for (int i=0; i<Chains.GetEntries(); i++){
    TChain * chain = (TChain*)(Chains[i]);
    chain->Add("/hepstore/store2/dervan/TestNtupleOutput_130GeV.root"); 
//    chain->Add("/hepstore/store2/dervan/TestNtupleOutput_Zbb.root"); 
  }

   int TruthNum       =truth->GetEntries();
   int RawElectronNum    = rawelec->GetEntries();
   int ElectronNum    = elec->GetEntries();
   int ElectronLooseNum    = elecloose->GetEntries();
   int ElectronMediumNum    = elecmedium->GetEntries();
   int ElectronMediumCALOISONum    = elecmediumcaloiso->GetEntries();
   int ElectronTightNum    = electight->GetEntries();
   int MuonNum     = muon->GetEntries();
   int ZNum        = ztree->GetEntries();
//   int HeeeeNum    = heeee->GetEntries();
//   int HmmmmNum    = hmmmm->GetEntries();
//   int HeemmNum    = heemm->GetEntries();
   int SelElectronNum    = selelec->GetEntries();
   int SelMuonNum    = selmuon->GetEntries();
   cout << "Number of truth  events = " << TruthNum << endl;
   cout << "Number of Raw Electron events = " << RawElectronNum << endl;
   cout << "Number of Electron events = " << ElectronNum << endl;
   cout << "Number of Loose Electron events = " << ElectronLooseNum << endl;
   cout << "Number of Medium Electron events = " << ElectronMediumNum << endl;
   cout << "Number of Medium CALOISOElectron events = " << ElectronMediumCALOISONum << endl;
   cout << "Number of Tight Electron events = " << ElectronTightNum << endl;
   cout << "Number of Muon   events = " << MuonNum << endl;
   cout << "Number of Z      events = " << ZNum << endl;
//   cout << "Number of Heeee  events = " << HeeeeNum << endl;
//   cout << "Number of Hmmmm  events = " << HmmmmNum << endl;
//   cout << "Number of Heemm  events = " << HeemmNum << endl;
   cout << "Number of Selected Electron events = " << SelElectronNum << endl;
   cout << "Number of Selected Muon     events = " << SelMuonNum << endl;
   getchar();


   Truth *TruthTree;
   Muon *PreselectedMuonTree;
   RawElectron *RawElectronTree;
   Electron *PreselectedElectronTree;
   Electron *PreselectedLooseElectronTree;
   Electron *PreselectedMediumElectronTree;
   Electron *PreselectedMediumCALOISOElectronTree;
   Electron *PreselectedTightElectronTree;
   ZCandidates *ZTree;
//   HCandidates *HeeeeCandidatesTree;
//   HCandidates *HmmmmCandidatesTree;
//   HCandidates *HeemmCandidatesTree;
//   HeeeeCandidatesTree = new HCandidates(heeee);
//   HmmmmCandidatesTree = new HCandidates(hmmmm);
//   HeemmCandidatesTree = new HCandidates(heemm);

   

   TruthTree = new Truth(truth);
   PreselectedMuonTree = new Muon(muon);
   RawElectronTree = new RawElectron(rawelec);
   PreselectedElectronTree = new Electron(elec);
   PreselectedLooseElectronTree = new Electron(elecloose);
   PreselectedMediumElectronTree = new Electron(elecmedium);
   PreselectedMediumCALOISOElectronTree = new Electron(elecmediumcaloiso);
   PreselectedTightElectronTree = new Electron(electight);
   ZTree = new ZCandidates(ztree);
//   HeeeeCandidatesTree = new HCandidates(heeee);
//   HmmmmCandidatesTree = new HCandidates(hmmmm);
//   HeemmCandidatesTree = new HCandidates(heemm);



   int const NumEntries = TruthNum;

   int numberofhiggsevents = 0;
   int numberofZZevents = 0;
   int Z1ee = 0;
   int Z1mm = 0;
   int Z1tt = 0;
   int Z2ee = 0;
   int Z2mm = 0;
   int Z2tt = 0;
   int Zeeee = 0;
   int Zmmmm = 0;
   int Zeemm = 0;
   int Zeett = 0;
   int Zmmtt = 0;
   int Ztttt = 0;

   int Zeeee_pass_e22i = 0;
   int Zmmmm_pass_e22i = 0;
   int Zeemm_pass_e22i = 0;
   int Zeett_pass_e22i = 0;
   int Zmmtt_pass_e22i = 0;
   int Ztttt_pass_e22i = 0;

   int Zeeee_pass_2e15i = 0;
   int Zmmmm_pass_2e15i = 0;
   int Zeemm_pass_2e15i = 0;
   int Zeett_pass_2e15i = 0;
   int Zmmtt_pass_2e15i = 0;
   int Ztttt_pass_2e15i = 0;

   int Zeeee_pass_mu20 = 0;
   int Zmmmm_pass_mu20 = 0;
   int Zeemm_pass_mu20 = 0;
   int Zeett_pass_mu20 = 0;
   int Zmmtt_pass_mu20 = 0;
   int Ztttt_pass_mu20 = 0;

   int Zeeee_pass_mu20_e22i = 0;
   int Zmmmm_pass_mu20_e22i = 0;
   int Zeemm_pass_mu20_e22i = 0;
   int Zeett_pass_mu20_e22i = 0;
   int Zmmtt_pass_mu20_e22i = 0;
   int Ztttt_pass_mu20_e22i = 0;

   int truthmuons = 0;
   int truthetamuons = 0;
   int recomuons = 0;
   int recomuonscombined = 0;
   int recoetamuons = 0;
   int loosemuons = 0;
   int mediummuons = 0;
   int mediumcaloisomuons = 0;
   int tightmuons = 0;

   int numberofgoodzs = 0;
   int numberofgoodzeeee =0;
   int numberofgoodzmmmm =0;
   int numberofgoodzeemm =0;

   for (int jentry=0; jentry<NumEntries;jentry++) 
     {
       TruthTree->LoadTree(jentry);
       TruthTree->GetEntry(jentry);
       PreselectedMuonTree->LoadTree(jentry);
       PreselectedMuonTree->GetEntry(jentry);
       RawElectronTree->LoadTree(jentry);
       RawElectronTree->GetEntry(jentry);
       PreselectedElectronTree->LoadTree(jentry);
       PreselectedElectronTree->GetEntry(jentry);
       PreselectedLooseElectronTree->LoadTree(jentry);
       PreselectedLooseElectronTree->GetEntry(jentry);
       PreselectedMediumElectronTree->LoadTree(jentry);
       PreselectedMediumElectronTree->GetEntry(jentry);
       PreselectedMediumCALOISOElectronTree->LoadTree(jentry);
       PreselectedMediumCALOISOElectronTree->GetEntry(jentry);
       PreselectedTightElectronTree->LoadTree(jentry);
       PreselectedTightElectronTree->GetEntry(jentry);
       ZTree->LoadTree(jentry);
       ZTree->GetEntry(jentry);
//       HeeeeCandidatesTree->LoadTree(jentry);
//       HeeeeCandidatesTree->GetEntry(jentry);
//       HmmmmCandidatesTree->LoadTree(jentry);
//       HmmmmCandidatesTree->GetEntry(jentry);
//       HeemmCandidatesTree->LoadTree(jentry);
//       HeemmCandidatesTree->GetEntry(jentry);
     
       if(!(jentry%1000)) { 
	 cout << "event number "<<jentry 
	      << " (" << int(double(100.*jentry)/NumEntries) << "\% done)"
	      << endl;
       }

       
       
       if(TruthTree->truth_H_num != 0) {
	 numberofhiggsevents++;
       }
       
       if(TruthTree->truth_Z_num >= 2) {
	 numberofZZevents++;
       }

       
       int d1electrons = 0;
       int d1muons = 0;
       int d1taus = 0;

       int numetamuons = 0;
       
       for(int i=0; i< (*(TruthTree->truth_daughter1_pt)).size(); i++)
	 {
	   if(abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 13) 
	     {
	       d1muons++;
	       hMuonTrue_eta->Fill(fabs((*(TruthTree->truth_daughter1_eta))[i])); 
	       hMuonTrue_pt->Fill(fabs((*(TruthTree->truth_daughter1_pt))[i])); 
	       if(fabs((*(TruthTree->truth_daughter1_eta))[i]) < 2.5) numetamuons++;
	     }
	   if(abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 11) d1electrons++;
	   if(abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 15) d1taus++;
	 }
       
       if(d1electrons == 2) Z1ee++;
       if(d1muons     == 2) Z1mm++;
       if(d1taus      == 2) Z1tt++;
       
       
       int d2electrons = 0;
       int d2muons = 0;
       int d2taus = 0;
       
       for(int i=0; i< (*(TruthTree->truth_daughter2_pt)).size(); i++)
	 {
	   if(abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 13 )
	     {
	       d2muons++;
	       hMuonTrue_eta->Fill(fabs((*(TruthTree->truth_daughter2_eta))[i])); 
	       hMuonTrue_pt->Fill(fabs((*(TruthTree->truth_daughter2_pt))[i])); 	    
	       if(fabs((*(TruthTree->truth_daughter2_eta))[i]) < 2.5) numetamuons++;
	     }
	   if(abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 11) d2electrons++;
	   if(abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 15) d2taus++;
	 }
       
       if(d2electrons == 2) Z2ee++;
       if(d2muons     == 2) Z2mm++;
       if(d2taus      == 2) Z2tt++;
       
       if((d1electrons == 2) && (d2electrons == 2) ) Zeeee++;
       if((d1muons     == 2) && (d2muons     == 2) ) Zmmmm++;
       if((d1electrons == 2) && (d2muons     == 2) ) Zeemm++;
       if((d2electrons == 2) && (d1muons     == 2) ) Zeemm++;
       if((d1electrons == 2) && (d2taus      == 2) ) Zeett++;
       if((d2electrons == 2) && (d1taus      == 2) ) Zeett++;
       if((d1muons     == 2) && (d2taus      == 2) ) Zmmtt++;
       if((d2muons     == 2) && (d1taus      == 2) ) Zmmtt++;
       if((d1taus      == 2) && (d2taus      == 2) ) Ztttt++;

       printf("\n");
       printf("----------------------------------------------------\n");
       if((d1electrons == 2) && (d2electrons == 2) ) printf("h->eeee\n");
       if((d1muons     == 2) && (d2muons     == 2) ) printf("h->mmmm\n");
       if((d1electrons == 2) && (d2muons     == 2) ) printf("h->eemm\n");
       if((d2electrons == 2) && (d1muons     == 2) ) printf("h->mmee\n");
       if((d1electrons == 2) && (d2taus      == 2) ) printf("h->eett\n");
       if((d2electrons == 2) && (d1taus      == 2) ) printf("h->ttee\n");
       if((d1muons     == 2) && (d2taus      == 2) ) printf("h->mmtt\n");
       if((d2muons     == 2) && (d1taus      == 2) ) printf("h->ttmm\n");
       if((d1taus      == 2) && (d2taus      == 2) ) printf("h->tttt\n");
       printf("\n");

       vector<double>  Z_pt;
       vector<double>  Z_eta;
       vector<double>  Z_phi;
       vector<double>  Z_m;
       vector<int>     Z_charge;
       vector<int>     Z_pdgid;
       vector<double>  Z1_pt;
       vector<double>  Z1_eta;
       vector<double>  Z1_phi;
       vector<double>  Z1_m;
       vector<int>     Z1_charge;
       vector<int>     Z1_pdgid;
       vector<double>  Z2_pt;
       vector<double>  Z2_eta;
       vector<double>  Z2_phi;
       vector<double>  Z2_m;
       vector<int>     Z2_charge;
       vector<int>     Z2_pdgid;


       printf("Number of electrons (Medium+CALOISO)= %i\n",(*(PreselectedMediumCALOISOElectronTree->electron_pt)).size());
       printf("Number of Zee Candidates = %i\n",(*(ZTree->Zee_m)).size());
       printf("Number of Zee Candidates = %i\n",(*(ZTree->Zee_element_pt)).size());
       for (int i=0; i<(*(ZTree->Zee_element_pt)).size(); i++)
	 {	   
	   printf("Z  mass = %f\n",(*(ZTree->Zee_m))[i]);
	   printf("Number of Zee daughters = %i\n",(*(ZTree->Zee_element_pt))[i].size());
	   int chargesume = 0;
	   int pt7sume = 0;
	   int pt20sume = 0;
	   for (int j=0; j<(*(ZTree->Zee_element_pt))[i].size(); j++)
	     {
	       printf("Number of Zee daughter pdgid  = %i\n",(*(ZTree->Zee_element_pdgid))[i][j]);
	       printf("Number of Zee daughter charge = %i\n",(*(ZTree->Zee_element_charge))[i][j]);
	       printf("Number of Zee daughter pt    = %f\n",(*(ZTree->Zee_element_pt))[i][j]);
	       printf("Number of Zee daughter eta   = %f\n",(*(ZTree->Zee_element_eta))[i][j]);
	       chargesume = chargesume + (*(ZTree->Zee_element_charge))[i][j];
	       if((*(ZTree->Zee_element_pt))[i][j] > 7.) pt7sume++;
	       if((*(ZTree->Zee_element_pt))[i][j] > 20.) pt20sume++;
	     }
	   printf("chargesume = %i\n",chargesume);
	   printf("pt7sume     = %i\n",pt7sume);
	   if(chargesume == 0 && pt7sume ==2)
	     {
	       Z_pt.push_back( (*(ZTree->Zee_pt))[i] );
	       Z_eta.push_back( (*(ZTree->Zee_eta))[i] );
	       Z_phi.push_back( (*(ZTree->Zee_phi))[i] );
	       Z_m.push_back( (*(ZTree->Zee_m))[i] );
	       Z_charge.push_back( (*(ZTree->Zee_charge))[i] );
	       Z_pdgid.push_back( (*(ZTree->Zee_pdgid))[i] );
	       Z1_pt.push_back( (*(ZTree->Zee_element_pt))[i][0] );
	       Z1_eta.push_back( (*(ZTree->Zee_element_eta))[i][0] );
	       Z1_phi.push_back( (*(ZTree->Zee_element_phi))[i][0] );
	       Z1_m.push_back( (*(ZTree->Zee_element_m))[i][0] );
	       Z1_charge.push_back( (*(ZTree->Zee_element_charge))[i][0] );
	       Z1_pdgid.push_back( (*(ZTree->Zee_element_pdgid))[i][0] );
	       Z2_pt.push_back( (*(ZTree->Zee_element_pt))[i][1] );
	       Z2_eta.push_back( (*(ZTree->Zee_element_eta))[i][1] );
	       Z2_phi.push_back( (*(ZTree->Zee_element_phi))[i][1] );
	       Z2_m.push_back( (*(ZTree->Zee_element_m))[i][1] );
	       Z2_charge.push_back( (*(ZTree->Zee_element_charge))[i][1] );
	       Z2_pdgid.push_back( (*(ZTree->Zee_element_pdgid))[i][1] );
	     }

	 }
       int goodZee = Z_pt.size();
       printf("number of good Zee = %i\n",goodZee);

       printf("\n");
       printf("Number of muons (isCombined) = %i\n",(*(PreselectedMuonTree->muon_pt)).size());
       printf("Number of Zmumu Candidates = %i\n",(*(ZTree->Zmumu_m)).size());
       printf("Number of Zmumu Candidates = %i\n",(*(ZTree->Zmumu_element_pt)).size());
       for (int i=0; i<(*(ZTree->Zmumu_element_pt)).size(); i++)
	 {	   
	   printf("Z  mass = %f\n",(*(ZTree->Zmumu_m))[i]);
	   printf("Number of Zmumu daughters = %i\n",(*(ZTree->Zmumu_element_pt))[i].size());
	   int chargesummu = 0;
	   int pt7summu = 0;
	   int pt20summu = 0;
	   for (int j=0; j<(*(ZTree->Zmumu_element_pt))[i].size(); j++)
	     {
	       printf("Number of Zmumu daughter pdgid  = %i\n",(*(ZTree->Zmumu_element_pdgid))[i][j]);
	       printf("Number of Zmumu daughter charge = %i\n",(*(ZTree->Zmumu_element_charge))[i][j]);
	       printf("Number of Zmumu daughter pt    = %f\n",(*(ZTree->Zmumu_element_pt))[i][j]);
	       printf("Number of Zmumu daughter eta   = %f\n",(*(ZTree->Zmumu_element_eta))[i][j]);
	       chargesummu = chargesummu + (*(ZTree->Zmumu_element_charge))[i][j];
	       if((*(ZTree->Zmumu_element_pt))[i][j] > 7.) pt7summu++;
	       if((*(ZTree->Zmumu_element_pt))[i][j] > 20.) pt20summu++;
	     }
	   printf("chargesummu = %i\n",chargesummu);
	   printf("pt7summu     = %i\n",pt7summu);
	   if(chargesummu == 0 && pt7summu ==2)
	     {
	       Z_pt.push_back( (*(ZTree->Zmumu_pt))[i] );
	       Z_eta.push_back( (*(ZTree->Zmumu_eta))[i] );
	       Z_phi.push_back( (*(ZTree->Zmumu_phi))[i] );
	       Z_m.push_back( (*(ZTree->Zmumu_m))[i] );
	       Z_charge.push_back( (*(ZTree->Zmumu_charge))[i] );
	       Z_pdgid.push_back( (*(ZTree->Zmumu_pdgid))[i] );
	       Z1_pt.push_back( (*(ZTree->Zmumu_element_pt))[i][0] );
	       Z1_eta.push_back( (*(ZTree->Zmumu_element_eta))[i][0] );
	       Z1_phi.push_back( (*(ZTree->Zmumu_element_phi))[i][0] );
	       Z1_m.push_back( (*(ZTree->Zmumu_element_m))[i][0] );
	       Z1_charge.push_back( (*(ZTree->Zmumu_element_charge))[i][0] );
	       Z1_pdgid.push_back( (*(ZTree->Zmumu_element_pdgid))[i][0] );
	       Z2_pt.push_back( (*(ZTree->Zmumu_element_pt))[i][1] );
	       Z2_eta.push_back( (*(ZTree->Zmumu_element_eta))[i][1] );
	       Z2_phi.push_back( (*(ZTree->Zmumu_element_phi))[i][1] );
	       Z2_m.push_back( (*(ZTree->Zmumu_element_m))[i][1] );
	       Z2_charge.push_back( (*(ZTree->Zmumu_element_charge))[i][1] );
	       Z2_pdgid.push_back( (*(ZTree->Zmumu_element_pdgid))[i][1] );
	     }

	 }




       int goodZmumu = Z_pt.size() - goodZee;
       printf("number of good Zmumu = %i\n",goodZmumu);


       printf("----------------------------------------------------\n");

       printf("number of good Zll   = %i\n",Z_pt.size());
       int onshellz = -1;
       for(int i=0 ; i<Z_pt.size(); i++)
	 {
	   printf("Z mass = %f\n",Z_m[i]);
	   printf("\tZ1_pt = %f\n",Z1_pt[i]);
	   printf("\tZ2_pt = %f\n",Z2_pt[i]);
	   printf("%f\n",fabs(Z_m[i]-90.20));
	   if( fabs(Z_m[i]-90.20) < 15. )
	     {
	       onshellz = i;
	       printf("ON SHELL Z in position %i\n",onshellz);	      
	     }
	 }

       printf("----------------------------------------------------\n");

       float maxzmass = 0.;
       int otherz = -1;
       if(onshellz >= 0)
	 {
	   for(int i=0 ; i<Z_pt.size(); i++)
	     {
	       if(i != onshellz)
		 {
		   printf("Z mass = %f\n",Z_m[i]);
		   if(Z_m[i] > 20.)
		     {
		       printf("\tZ1_pt = %f\n",Z1_pt[i]);
		       printf("\tZ2_pt = %f\n",Z2_pt[i]);
		       if( Z_m[i] > maxzmass) 
			 {
			   maxzmass = Z_m[i];
			   otherz = i;
			   printf("maxzmass = %f\n",Z_m[i]);
			 }	       
		     }
		 }
	     }
	 }

       printf("----------------------------------------------------\n");

       TLorentzVector Heeee;
       TLorentzVector Hmmmm;
       TLorentzVector Heemm;

       if(onshellz >=0 && otherz >=0)
	 {
	   printf("SELECTED A GOOD PAIR OF Zs\n");
	   printf("Z 1 mass = %f\n",Z_m[onshellz]);
	   printf("\tZ1_pdgid = %i \tZ1_pt = %f\n",Z1_pdgid[onshellz],Z1_pt[onshellz]);
	   printf("\tZ2_pdgid = %i \tZ2_pt = %f\n",Z2_pdgid[onshellz],Z2_pt[onshellz]);
	   printf("Z 2 mass = %f\n",Z_m[otherz]);
	   printf("\tZ1_pdgid = %i \tZ1_pt = %f\n",Z1_pdgid[otherz],Z1_pt[otherz]);
	   printf("\tZ2_pdgid = %i \tZ2_pt = %f\n",Z2_pdgid[otherz],Z2_pt[otherz]);
	   numberofgoodzs++;
	   if(abs(Z1_pdgid[onshellz]) == 11 && abs(Z2_pdgid[otherz]) == 11) 
	     {
	       hZeeee_mass->Fill(Z_m[onshellz]);
	       hZeeee_mass->Fill(Z_m[otherz]);
	       numberofgoodzeeee++;
	       Heeee.SetPtEtaPhiM(Z1_pt[onshellz]+Z1_pt[otherz],Z1_eta[onshellz]+Z1_eta[otherz],Z1_phi[onshellz]+Z1_phi[otherz],Z1_m[onshellz]+Z1_m[otherz]),
	       hHeeee_mass->Fill(Heeee.M());
	     }
	   if(abs(Z1_pdgid[onshellz]) == 13 && abs(Z2_pdgid[otherz]) == 13) 
	     {
	       hZmmmm_mass->Fill(Z_m[onshellz]);
	       hZmmmm_mass->Fill(Z_m[otherz]);
	       numberofgoodzmmmm++;
	       Hmmmm.SetPtEtaPhiM(Z1_pt[onshellz]+Z1_pt[otherz],Z1_eta[onshellz]+Z1_eta[otherz],Z1_phi[onshellz]+Z1_phi[otherz],Z1_m[onshellz]+Z1_m[otherz]),
	       hHmmmm_mass->Fill(Hmmmm.M());
	     }
	   if(abs(Z1_pdgid[onshellz]) == 11 && abs(Z2_pdgid[otherz]) == 13) 
	     {
	       hZeemm_mass->Fill(Z_m[onshellz]);
	       hZeemm_mass->Fill(Z_m[otherz]);
	       numberofgoodzeemm++;
	       Heemm.SetPtEtaPhiM(Z1_pt[onshellz]+Z1_pt[otherz],Z1_eta[onshellz]+Z1_eta[otherz],Z1_phi[onshellz]+Z1_phi[otherz],Z1_m[onshellz]+Z1_m[otherz]),
	       hHeemm_mass->Fill(Heemm.M());
	     }
	   if(abs(Z1_pdgid[onshellz]) == 13 && abs(Z2_pdgid[otherz]) == 11) 
	     {
	       hZeemm_mass->Fill(Z_m[onshellz]);
	       hZeemm_mass->Fill(Z_m[otherz]);
	       numberofgoodzeemm++;
	       Heemm.SetPtEtaPhiM(Z1_pt[onshellz]+Z1_pt[otherz],Z1_eta[onshellz]+Z1_eta[otherz],Z1_phi[onshellz]+Z1_phi[otherz],Z1_m[onshellz]+Z1_m[otherz]),
	       hHeemm_mass->Fill(Heemm.M());
	     }
	 }

       cout << "Heemm" << endl;

       printf("----------------------------------------------------\n");
//       getchar();


       
















       truthmuons = truthmuons + d1muons+d2muons;
       truthetamuons = truthetamuons + numetamuons;

       for(int i=0; i< (*(PreselectedMuonTree->muon_pt)).size(); i++)
	 {
	   float deltaeta = 10000.;
	   float deltaphi = 10000.;
	   float deltaR = 10000.;
	   for(int j=0; j< (*(TruthTree->truth_daughter1_pt)).size(); j++)
	     {
	       deltaR = 10000.;
	       if(abs((*(TruthTree->truth_daughter1_pdgid))[j]) == 13 ) 
		 {
		   deltaeta = (*(TruthTree->truth_daughter1_eta))[j] - (*(PreselectedMuonTree->muon_eta))[i];
		   deltaphi = (*(TruthTree->truth_daughter1_phi))[j] - (*(PreselectedMuonTree->muon_phi))[i];
		   deltaR = sqrt(pow(deltaeta,2.) + pow(deltaphi,2.) );
		 }
	       if(deltaR < 0.02) 
		 {
		   hMuonReco_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
		   hMuonReco_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
		   recomuons++;	    
                   if( (*(PreselectedMuonTree->muon_isCombinedMuon))[i] == 1)
                     {
		        hMuonCombined_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
		        hMuonCombined_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
                        recomuonscombined++;
                     }
		 }
	     }
	   for(int j=0; j< (*(TruthTree->truth_daughter2_pt)).size(); j++)
	     {
	       deltaR = 10000.;
	       if(abs((*(TruthTree->truth_daughter2_pdgid))[j]) == 13 ) 
		 {
		   deltaeta = (*(TruthTree->truth_daughter2_eta))[j] - (*(PreselectedMuonTree->muon_eta))[i];
		   deltaphi = (*(TruthTree->truth_daughter2_phi))[j] - (*(PreselectedMuonTree->muon_phi))[i];
		   deltaR = sqrt(pow(deltaeta,2.) + pow(deltaphi,2.) );
		 }
	       if(deltaR < 0.02) 
		 {
		   hMuonReco_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
		   hMuonReco_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
		   recomuons++;	     
                   if( (*(PreselectedMuonTree->muon_isCombinedMuon))[i] == 1)
                     {
		        hMuonCombined_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
		        hMuonCombined_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
                        recomuonscombined++;
                     }
		 }
	     }
	 }


       
 //      int pass_mu20 = 0;
 //      int pass_e22i = 0;
 //      int pass_2e15i = 0;
	
// //	if( ((*(TrigResultsTree->Passed_mu20))[0] == 1) && ((*(TrigResultsTree->Passed_mu20))[1] == 1) && ((*(TrigResultsTree->Passed_mu20))[2] == 1) ) pass_mu20 = 1;
//       if( ((*(TrigResultsTree->Passed_mu20))[1] == 1) && ((*(TrigResultsTree->Passed_mu20))[2] == 1) ) pass_mu20 = 1;
//       if( ((*(TrigResultsTree->Passed_e22i_L33))[0] == 1) && ((*(TrigResultsTree->Passed_e22i_L33))[1] == 1) && ((*(TrigResultsTree->Passed_e22i_L33))[2] == 1) ) pass_e22i = 1;
//       if( ((*(TrigResultsTree->Passed_2e15i))[0] == 1) && ((*(TrigResultsTree->Passed_2e15i))[1] == 1) && ((*(TrigResultsTree->Passed_2e15i))[2] == 1) ) pass_2e15i = 1;
	

//       if((d1muons == 2) && (d2muons == 2) && (pass_e22i == 1) ) Zeeee_pass_e22i++;
//       if((d1muons     == 2) && (d2muons     == 2) && (pass_e22i == 1) ) Zmmmm_pass_e22i++;
//       if((d1muons == 2) && (d2muons     == 2) && (pass_e22i == 1) ) Zeemm_pass_e22i++;
//       if((d2muons == 2) && (d1muons     == 2) && (pass_e22i == 1) ) Zeemm_pass_e22i++;
//       if((d1muons == 2) && (d2taus      == 2) && (pass_e22i == 1) ) Zeett_pass_e22i++;
//       if((d2muons == 2) && (d1taus      == 2) && (pass_e22i == 1) ) Zeett_pass_e22i++;
//       if((d1muons     == 2) && (d2taus      == 2) && (pass_e22i == 1) ) Zmmtt_pass_e22i++;
//       if((d2muons     == 2) && (d1taus      == 2) && (pass_e22i == 1) ) Zmmtt_pass_e22i++;
//       if((d1taus      == 2) && (d2taus      == 2) && (pass_e22i == 1) ) Ztttt_pass_e22i++;
       
//       if((d1muons == 2) && (d2muons == 2) && (pass_2e15i == 1) ) Zeeee_pass_2e15i++;
//       if((d1muons     == 2) && (d2muons     == 2) && (pass_2e15i == 1) ) Zmmmm_pass_2e15i++;
//       if((d1muons == 2) && (d2muons     == 2) && (pass_2e15i == 1) ) Zeemm_pass_2e15i++;
//       if((d2muons == 2) && (d1muons     == 2) && (pass_2e15i == 1) ) Zeemm_pass_2e15i++;
//       if((d1muons == 2) && (d2taus      == 2) && (pass_2e15i == 1) ) Zeett_pass_2e15i++;
//       if((d2muons == 2) && (d1taus      == 2) && (pass_2e15i == 1) ) Zeett_pass_2e15i++;
//       if((d1muons     == 2) && (d2taus      == 2) && (pass_2e15i == 1) ) Zmmtt_pass_2e15i++;
//       if((d2muons     == 2) && (d1taus      == 2) && (pass_2e15i == 1) ) Zmmtt_pass_2e15i++;
//       if((d1taus      == 2) && (d2taus      == 2) && (pass_2e15i == 1) ) Ztttt_pass_2e15i++;
       
//       if((d1muons == 2) && (d2muons == 2) && (pass_mu20 == 1) ) Zeeee_pass_mu20++;
//       if((d1muons     == 2) && (d2muons     == 2) && (pass_mu20 == 1) ) Zmmmm_pass_mu20++;
//       if((d1muons == 2) && (d2muons     == 2) && (pass_mu20 == 1) ) Zeemm_pass_mu20++;
//       if((d2muons == 2) && (d1muons     == 2) && (pass_mu20 == 1) ) Zeemm_pass_mu20++;
//       if((d1muons == 2) && (d2taus      == 2) && (pass_mu20 == 1) ) Zeett_pass_mu20++;
//       if((d2muons == 2) && (d1taus      == 2) && (pass_mu20 == 1) ) Zeett_pass_mu20++;
//       if((d1muons     == 2) && (d2taus      == 2) && (pass_mu20 == 1) ) Zmmtt_pass_mu20++;
//       if((d2muons     == 2) && (d1taus      == 2) && (pass_mu20 == 1) ) Zmmtt_pass_mu20++;
//       if((d1taus      == 2) && (d2taus      == 2) && (pass_mu20 == 1) ) Ztttt_pass_mu20++;
        
//       if(((d1muons == 2) && (d2muons == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zeeee_pass_mu20_e22i++;
//       if(((d1muons     == 2) && (d2muons     == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zmmmm_pass_mu20_e22i++;
//       if(((d1muons == 2) && (d2muons     == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zeemm_pass_mu20_e22i++;
//       if(((d2muons == 2) && (d1muons     == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zeemm_pass_mu20_e22i++;
//       if(((d1muons == 2) && (d2taus      == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zeett_pass_mu20_e22i++;
//       if(((d2muons == 2) && (d1taus      == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zeett_pass_mu20_e22i++;
//       if(((d1muons     == 2) && (d2taus      == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zmmtt_pass_mu20_e22i++;
//       if(((d2muons     == 2) && (d1taus      == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Zmmtt_pass_mu20_e22i++;
//       if(((d1taus      == 2) && (d2taus      == 2)) && ((pass_mu20 == 1) || (pass_e22i == 1))) Ztttt_pass_mu20_e22i++;


//       for(int i=0; i< (*(TruthTree->truth_daughter1_pt)).size(); i++)
// 	 {
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 13) ) hMuonTrigger_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 11) ) hMuonTrigger_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
//
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 13) && ((*(TrigResultsTree->Passed_mu20))[0] == 1) ) hMuonTrigger_L1_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 13) && ((*(TrigResultsTree->Passed_mu20))[1] == 1) ) hMuonTrigger_L2_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 13) && ((*(TrigResultsTree->Passed_mu20))[2] == 1) ) hMuonTrigger_EF_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 11) && ((*(TrigResultsTree->Passed_e22i_L33))[0] == 1) ) hMuonTrigger_L1_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 11) && ((*(TrigResultsTree->Passed_e22i_L33))[1] == 1) ) hMuonTrigger_L2_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	   if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 11) && ((*(TrigResultsTree->Passed_e22i_L33))[2] == 1) ) hMuonTrigger_EF_pt->Fill((*(TruthTree->truth_daughter1_pt))[i]);
// 	 }

//       for(int i=0; i< (*(TruthTree->truth_daughter2_pt)).size(); i++)
// 	 {
// 	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 13) ) hMuonTrigger_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 11) ) hMuonTrigger_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 13) && ((*(TrigResultsTree->Passed_mu20))[0] == 1) ) hMuonTrigger_L1_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 13) && ((*(TrigResultsTree->Passed_mu20))[1] == 1) ) hMuonTrigger_L2_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 13) && ((*(TrigResultsTree->Passed_mu20))[2] == 1) ) hMuonTrigger_EF_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 11) && ((*(TrigResultsTree->Passed_e22i_L33))[0] == 1) ) hMuonTrigger_L1_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 11) && ((*(TrigResultsTree->Passed_e22i_L33))[1] == 1) ) hMuonTrigger_L2_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
//	   if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 11) && ((*(TrigResultsTree->Passed_e22i_L33))[2] == 1) ) hMuonTrigger_EF_pt->Fill((*(TruthTree->truth_daughter2_pt))[i]);
// 	 }

//       if((pass_mu20 == 1) || (pass_e22i == 1))
// 	 {
//	   cout << "PASSED TRIGGER" << endl;

//	   for(int i=0; i< (*(TruthTree->truth_daughter1_pt)).size(); i++)
//	     {
//	       if( (abs((*(TruthTree->truth_daughter1_pdgid))[i]) == 11) )
//	       {
//		 hMuonTrue_eta->Fill(fabs((*(TruthTree->truth_daughter1_eta))[i]));
//		 hMuonTrue_pt->Fill(fabs((*(TruthTree->truth_daughter1_pt))[i]));
//		 cout << "truth particle pt = " << fabs((*(TruthTree->truth_daughter1_pt))[i]) << " eta = " << fabs((*(TruthTree->truth_daughter1_eta))[i]) << endl;
//	       }
//	     }
//	   for(int i=0; i< (*(TruthTree->truth_daughter2_pt)).size(); i++)
///	     {
//	       if( (abs((*(TruthTree->truth_daughter2_pdgid))[i]) == 11) )
///	       {
//		 hMuonTrue_eta->Fill(fabs((*(TruthTree->truth_daughter2_eta))[i]));
//		 hMuonTrue_pt->Fill(fabs((*(TruthTree->truth_daughter2_pt))[i]));
//		 cout << "truth particle pt = " << fabs((*(TruthTree->truth_daughter2_pt))[i]) << " eta = " << fabs((*(TruthTree->truth_daughter2_eta))[i]) << endl;;
//	       }
//	     }

// 	   if((d1muons == 2) && (d2muons == 2) )
// 	     {
// 	       cout << "eeee event" << endl;
// 	     }
// 	   if((d1muons     == 2) && (d2muons     == 2) ) 
// 	     {
// 	       cout << "mmmm event" << endl;
// 	     }
// 	   if((d1muons == 2) && (d2muons     == 2) )
// 	     {
// 	       cout << "eemm event" << endl;
// 	     }
// 	   if((d2muons == 2) && (d1muons     == 2) ) 
// 	     {
// 	       cout << "eemm event" << endl;
// 	     }
//	   cout << "Num of preselected muons = " << (*(PreselectedMuonTree->muon_isloose)).size() << endl;
//	   for(int i=0; i< (*(PreselectedMuonTree->muon_pt)).size(); i++)
//	     {
//	       cout << i << " Loose/Medium/Tight  Muon = " << (*(PreselectedMuonTree->muon_isloose))[i] << " "<< (*(PreselectedMuonTree->muon_ismedium))[i] << " " << (*(PreselectedMuonTree->muon_istight))[i] << " pt " << fabs((*(PreselectedMuonTree->muon_pt))[i]) << " eta " << fabs((*(PreselectedMuonTree->muon_eta))[i]) <<endl;
//	       if((*(PreselectedMuonTree->muon_isloose))[i] == 1)
//		 {
//		   hMuonLoose_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
//		   hMuonLoose_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
//		 }
//	       if((*(PreselectedMuonTree->muon_ismedium))[i] == 1)
//		 {
//		   hMuonMedium_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
//		   hMuonMedium_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
//		 }
//	       if((*(PreselectedMuonTree->muon_istight))[i] == 1)
//		 {
//		   hMuonTight_eta->Fill(fabs((*(PreselectedMuonTree->muon_eta))[i])); 
//		   hMuonTight_pt->Fill(fabs((*(PreselectedMuonTree->muon_pt))[i])); 
//		 }

//	     }
//	   cout << "Num of raw muons = " << (*(RawMuonTree->muon_isloose)).size() << endl;
//	   for(int i=0; i< (*(RawMuonTree->muon_pt)).size(); i++)
//	     {
//	       cout << i << " Loose/Medium/Tight  Muon = " << (*(RawMuonTree->muon_isloose))[i] << " "<< (*(RawMuonTree->muon_ismedium))[i] << " " << (*(RawMuonTree->muon_istight))[i] << " pt " << fabs((*(RawMuonTree->muon_pt))[i]) << " eta " << fabs((*(RawMuonTree->muon_eta))[i]) <<endl;
//	       if((*(RawMuonTree->muon_isloose))[i] == 1)
//		 {
//		   hMuonLoose_eta->Fill(fabs((*(RawMuonTree->muon_eta))[i])); 
//		 }

//	     }
//	   getchar();

//	 }
	
     }
   
   cout << "Number of Higgs events = " << numberofhiggsevents << endl;
   cout << "Number of ZZ    events = " << numberofZZevents << endl;
   cout << " " << endl;
   cout << "Number of Z1 ee events = " << Z1ee << endl;
   cout << "Number of Z1 mm events = " << Z1mm << endl;
   cout << "Number of Z1 tt events = " << Z1tt << endl;
   cout << "Number of Z1    events = " << Z1ee+Z1mm+Z1tt << endl;
   cout << " " << endl;
   cout << "Number of Z2 ee events = " << Z2ee << endl;
   cout << "Number of Z2 mm events = " << Z2mm << endl;
   cout << "Number of Z2 tt events = " << Z2tt << endl;
   cout << "Number of Z2    events = " << Z2ee+Z2mm+Z2tt << endl;
   cout << " " << endl;
   cout << "Number of Zeeee events = " << Zeeee << "\t (" << float(Zeeee)/NumEntries*100 << "%)" << endl;
   cout << "Number of Zmmmm events = " << Zmmmm << "\t (" << float(Zmmmm)/NumEntries*100 << "%)" << endl;
   cout << "Number of Zeemm events = " << Zeemm << "\t (" << float(Zeemm)/NumEntries*100 << "%)" << endl;
   cout << "Number of Zeett events = " << Zeett << "\t (" << float(Zeett)/NumEntries*100 << "%)" << endl;
   cout << "Number of Zmmtt events = " << Zmmtt << "\t (" << float(Zmmtt)/NumEntries*100 << "%)" << endl;
   cout << "Number of Ztttt events = " << Ztttt << "\t (" << float(Ztttt)/NumEntries*100 << "%)" << endl;
   cout << "Number of Zllll events = " << Zeeee+Zmmmm+Zeemm+Zeett+Zmmtt+Ztttt << endl;
   cout << " " << endl;
   cout << " " << endl;
   cout << " " << endl;
   cout << "Number of truth muons = " << truthmuons << endl;
   cout << "Number of truth muons (eta < 2.5) = " << truthetamuons << "\t (" << (float(truthetamuons)/float(truthmuons))*100. << "%)" << endl;
   cout << "Number of reco muons = " << recomuons << "\t (" << (float(recomuons)/float(truthmuons))*100. << "%)" << endl;
   cout << "Number of reco muons combined= " << recomuonscombined << "\t (" << (float(recomuonscombined)/float(truthmuons))*100. << "%)" << endl;
//   cout << "e22i" << endl;
//   cout << "Number of Zeeee events = " << Zeeee_pass_e22i << "\t (" << float(Zeeee_pass_e22i)/Zeeee*100 << "%)" << endl;
//   cout << "Number of Zmmmm events = " << Zmmmm_pass_e22i << "\t (" << float(Zmmmm_pass_e22i)/Zmmmm*100 << "%)" << endl;
//   cout << "Number of Zeemm events = " << Zeemm_pass_e22i << "\t (" << float(Zeemm_pass_e22i)/Zeemm*100 << "%)" << endl;
//   cout << "Number of Zeett events = " << Zeett_pass_e22i << "\t (" << float(Zeett_pass_e22i)/Zeett*100 << "%)" << endl;
//   cout << "Number of Zmmtt events = " << Zmmtt_pass_e22i << "\t (" << float(Zmmtt_pass_e22i)/Zmmtt*100 << "%)" << endl;
//   cout << "Number of Ztttt events = " << Ztttt_pass_e22i << "\t (" << float(Ztttt_pass_e22i)/Ztttt*100 << "%)" << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   cout << "2e15i" << endl;
//   cout << "Number of Zeeee events = " << Zeeee_pass_2e15i << "\t (" << float(Zeeee_pass_2e15i)/Zeeee*100 << "%)" << endl;
//   cout << "Number of Zmmmm events = " << Zmmmm_pass_2e15i << "\t (" << float(Zmmmm_pass_2e15i)/Zmmmm*100 << "%)" << endl;
//   cout << "Number of Zeemm events = " << Zeemm_pass_2e15i << "\t (" << float(Zeemm_pass_2e15i)/Zeemm*100 << "%)" << endl;
//   cout << "Number of Zeett events = " << Zeett_pass_2e15i << "\t (" << float(Zeett_pass_2e15i)/Zeett*100 << "%)" << endl;
//   cout << "Number of Zmmtt events = " << Zmmtt_pass_2e15i << "\t (" << float(Zmmtt_pass_2e15i)/Zmmtt*100 << "%)" << endl;
//   cout << "Number of Ztttt events = " << Ztttt_pass_2e15i << "\t (" << float(Ztttt_pass_2e15i)/Ztttt*100 << "%)" << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   cout << "mu20" << endl;
//   cout << "Number of Zeeee events = " << Zeeee_pass_mu20 << "\t (" << float(Zeeee_pass_mu20)/Zeeee*100 << "%)" << endl;
//  cout << "Number of Zmmmm events = " << Zmmmm_pass_mu20 << "\t (" << float(Zmmmm_pass_mu20)/Zmmmm*100 << "%)" << endl;
//  cout << "Number of Zeemm events = " << Zeemm_pass_mu20 << "\t (" << float(Zeemm_pass_mu20)/Zeemm*100 << "%)" << endl;
//  cout << "Number of Zeett events = " << Zeett_pass_mu20 << "\t (" << float(Zeett_pass_mu20)/Zeett*100 << "%)" << endl;
//  cout << "Number of Zmmtt events = " << Zmmtt_pass_mu20 << "\t (" << float(Zmmtt_pass_mu20)/Zmmtt*100 << "%)" << endl;
//   cout << "Number of Ztttt events = " << Ztttt_pass_mu20 << "\t (" << float(Ztttt_pass_mu20)/Ztttt*100 << "%)" << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   cout << "mu20 or e22i" << endl;
//   cout << "Number of Zeeee events = " << Zeeee_pass_mu20_e22i << "\t (" << (float(Zeeee_pass_mu20_e22i)/(Zeeee))*100 << "%)" << endl;
//   cout << "Number of Zmmmm events = " << Zmmmm_pass_mu20_e22i << "\t (" << (float(Zmmmm_pass_mu20_e22i)/(Zmmmm))*100 << "%)" << endl;
//   cout << "Number of Zeemm events = " << Zeemm_pass_mu20_e22i << "\t (" << (float(Zeemm_pass_mu20_e22i)/(Zeemm))*100 << "%)" << endl;
//   cout << "Number of Zeett events = " << Zeett_pass_mu20_e22i << "\t (" << (float(Zeett_pass_mu20_e22i)/Zeett)*100 << "%)" << endl;
//   cout << "Number of Zmmtt events = " << Zmmtt_pass_mu20_e22i << "\t (" << (float(Zmmtt_pass_mu20_e22i)/Zmmtt)*100 << "%)" << endl;
//   cout << "Number of Ztttt events = " << Ztttt_pass_mu20_e22i << "\t (" << (float(Ztttt_pass_mu20_e22i)/Ztttt)*100 << "%)" << endl;

//   cout << " " << endl;
//   cout << " " << endl;
//   cout << " " << endl;
//   printf("Trigger Menu    |\t 4e   |\t 4mu  |\t2e2mu |\n"); 
//   printf("-----------------------------------------------------\n");
//   printf("1mu20           |\t%4.1f  |\t%4.1f  |\t%4.1f  |\n",float(Zeeee_pass_mu20)/Zeeee*100,float(Zmmmm_pass_mu20)/Zmmmm*100,float(Zeemm_pass_mu20)/Zeemm*100);
//   printf("1e22i           |\t%4.1f  |\t%4.1f  |\t%4.1f  |\n",float(Zeeee_pass_e22i)/Zeeee*100,float(Zmmmm_pass_e22i)/Zmmmm*100,float(Zeemm_pass_e22i)/Zeemm*100);
//   printf("2e15i           |\t%4.1f  |\t%4.1f  |\t%4.1f  |\n",float(Zeeee_pass_2e15i)/Zeeee*100,float(Zmmmm_pass_2e15i)/Zmmmm*100,float(Zeemm_pass_2e15i)/Zeemm*100);
//   printf("mu20 or e22i    |\t%4.1f  |\t%4.1f  |\t%4.1f  |\n",float(Zeeee_pass_mu20_e22i)/Zeeee*100,float(Zmmmm_pass_mu20_e22i)/Zmmmm*100,float(Zeemm_pass_mu20_e22i)/Zeemm*100);


//   hMuonTrigger_L1_pt->Divide(hMuonTrigger_pt);
//   hMuonTrigger_L1_pt->Divide(hMuonTrigger_pt);

//   for(int i=0; i<96; i++)
//     {
//       cout << hMuonTrigger_L1_pt->GetBinContent(i) << endl;
//       cout << hMuonTrigger_L2_pt->GetBinContent(i) << endl;
//       cout << hMuonTrigger_EF_pt->GetBinContent(i) << endl;
//       cout << hMuonTrigger_pt->GetBinContent(i) << endl;
//     }

//   TH1F *hMuonTrigger_L1_eff = (TH1F*)hMuonTrigger_L1_pt->Clone("hMuonTrigger_L1_eff");
//   hMuonTrigger_L1_eff->Divide(hMuonTrigger_pt);
//   TH1F *hMuonTrigger_L2_eff = (TH1F*)hMuonTrigger_L2_pt->Clone("hMuonTrigger_L2_eff");
//   hMuonTrigger_L2_eff->Divide(hMuonTrigger_pt);
//   TH1F *hMuonTrigger_EF_eff = (TH1F*)hMuonTrigger_EF_pt->Clone("hMuonTrigger_EF_eff");
//   hMuonTrigger_EF_eff->Divide(hMuonTrigger_pt);

//   TH1F *hMuonLoose_eta_eff = (TH1F*)hMuonLoose_eta->Clone("hMuonLoose_eta_eff");
//   hMuonLoose_eta_eff->Divide(hMuonTrue_eta);
//   TH1F *hMuonMedium_eta_eff = (TH1F*)hMuonMedium_eta->Clone("hMuonMedium_eta_eff");
//   hMuonMedium_eta_eff->Divide(hMuonTrue_eta);
//   TH1F *hMuonTight_eta_eff = (TH1F*)hMuonTight_eta->Clone("hMuonTight_eta_eff");
//   hMuonTight_eta_eff->Divide(hMuonTrue_eta);

//   TH1F *hMuonLoose_pt_eff = (TH1F*)hMuonLoose_pt->Clone("hMuonLoose_pt_eff");
//   hMuonLoose_pt_eff->Divide(hMuonTrue_pt);
//   TH1F *hMuonMedium_pt_eff = (TH1F*)hMuonMedium_pt->Clone("hMuonMedium_pt_eff");
//   hMuonMedium_pt_eff->Divide(hMuonTrue_pt);
//   TH1F *hMuonTight_pt_eff = (TH1F*)hMuonTight_pt->Clone("hMuonTight_pt_eff");
//   hMuonTight_pt_eff->Divide(hMuonTrue_pt);


   TH1F *hMuonReco_pt_test       = new TH1F("hMuonReco_pt_test","Muon_pt_test",55,5.,60.);
   TH1F *hMuonCombined_pt_test   = new TH1F("hMuonCombined_pt_test","MuonCombined_pt_test",55,5.,60.);
   TH1F *hMuonReco_eta_test      = new TH1F("hMuonReco_eta_test","Muon_eta_test",25,0.,2.5);
   TH1F *hMuonCombined_eta_test  = new TH1F("hMuonCombined_eta_test","MuonCombined_eta_test",25,0.,2.5);
   for(int i=1; i<56; i++)
     {
       hMuonReco_pt_test->SetBinContent(i,hMuonReco_pt->GetBinContent(i)/hMuonTrue_pt->GetBinContent(i));
       hMuonReco_pt_test->SetBinError(i,((hMuonReco_pt->GetBinContent(i)/hMuonTrue_pt->GetBinContent(i))*((sqrt(hMuonReco_pt->GetBinContent(i))/hMuonReco_pt->GetBinContent(i))+(sqrt(hMuonTrue_pt->GetBinContent(i))/hMuonTrue_pt->GetBinContent(i)))));
       
       hMuonCombined_pt_test->SetBinContent(i,hMuonCombined_pt->GetBinContent(i)/hMuonTrue_pt->GetBinContent(i));
       hMuonCombined_pt_test->SetBinError(i,((hMuonCombined_pt->GetBinContent(i)/hMuonTrue_pt->GetBinContent(i))*((sqrt(hMuonCombined_pt->GetBinContent(i))/hMuonCombined_pt->GetBinContent(i))+(sqrt(hMuonTrue_pt->GetBinContent(i))/hMuonTrue_pt->GetBinContent(i)))));
       
       printf("PT %i %f %f %f  \n",i,hMuonCombined_pt->GetBinContent(i),hMuonTrue_pt->GetBinContent(i),hMuonCombined_pt->GetBinContent(i)/hMuonTrue_pt->GetBinContent(i));
     }

   printf("\n\n\n\n\n");

   for(int i=1; i<26; i++)	
     {
       hMuonReco_eta_test->SetBinContent(i,hMuonReco_eta->GetBinContent(i)/hMuonTrue_eta->GetBinContent(i));
       hMuonReco_eta_test->SetBinError(i,((hMuonReco_eta->GetBinContent(i)/hMuonTrue_eta->GetBinContent(i))*((sqrt(hMuonReco_eta->GetBinContent(i))/hMuonReco_eta->GetBinContent(i))+(sqrt(hMuonTrue_eta->GetBinContent(i))/hMuonTrue_eta->GetBinContent(i)))));
       
       hMuonCombined_eta_test->SetBinContent(i,hMuonCombined_eta->GetBinContent(i)/hMuonTrue_eta->GetBinContent(i));
       hMuonCombined_eta_test->SetBinError(i,((hMuonCombined_eta->GetBinContent(i)/hMuonTrue_eta->GetBinContent(i))*((sqrt(hMuonCombined_eta->GetBinContent(i))/hMuonCombined_eta->GetBinContent(i))+(sqrt(hMuonTrue_eta->GetBinContent(i))/hMuonTrue_eta->GetBinContent(i)))));
     
       printf("ETA %i %f %f %f  \n",i,hMuonCombined_eta->GetBinContent(i),hMuonTrue_eta->GetBinContent(i),hMuonCombined_eta->GetBinContent(i)/hMuonTrue_eta->GetBinContent(i));
     } 

   cout << "NUMBER OF GOOD Zeeee = " << numberofgoodzeeee << "\t (" <<  float(numberofgoodzeeee*100)/float(Zeeee) << "%)" << endl; 
   cout << "NUMBER OF GOOD Zmmmm = " << numberofgoodzmmmm << "\t (" <<  float(numberofgoodzmmmm*100)/float(Zmmmm) << "%)" << endl; 
   cout << "NUMBER OF GOOD Zeemm = " << numberofgoodzeemm << "\t (" <<  float(numberofgoodzeemm*100)/float(Zeemm) << "%)" << endl; 

   
   
   hfile.Write();
   
   return 0;
}