#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"Trigger.h"
#include"Muon.h"



//______________________________________________________________________________
int main()
{

  TFile hfile("output_muons.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.);

  // instanciate chains
  TObjArray Chains;
  TChain * truth = new TChain("TruthTree"); Chains.Add(truth);
  TChain * muon = new TChain("PreselectedMuonTree"); Chains.Add(muon);

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

   int TruthNum       =truth->GetEntries();
   int MuonNum    = muon->GetEntries();
   cout << "Number of truth events = " << TruthNum << endl;
   cout << "Number of Muon events = " << MuonNum << endl;
   getchar();


   Truth *TruthTree;
   Muon *PreselectedMuonTree;

   

   TruthTree = new Truth(truth);
   PreselectedMuonTree = new Muon(muon);

	printf("HERE\n");

   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;

   for (int jentry=0; jentry<NumEntries;jentry++) 
     {
       TruthTree->LoadTree(jentry);
       TruthTree->GetEntry(jentry);
       PreselectedMuonTree->LoadTree(jentry);
       PreselectedMuonTree->GetEntry(jentry);
     
       if(!(jentry%1000)) { 
	 cout << "event number "<<jentry 
	      << " (" << int(double(100.*jentry)/NumEntries) << "\% done)"
	      << endl;
       }

//         printf("Number of muons = %i\n",(*(PreselectedMuonTree->muon_pt)).size());
//         printf("Number of loose muons = %i\n",(*(PreselectedLooseMuonTree->muon_pt)).size());
//         printf("Number of medium muons = %i\n",(*(PreselectedMediumMuonTree->muon_pt)).size());
//         printf("Number of medium CALOISO muons = %i\n",(*(PreselectedMediumCALOISOMuonTree->muon_pt)).size());
//         getchar();
       
       
       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++;

//        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");
//        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));
     } 
   
   
   hfile.Write();
   
   return 0;
}