/*


*/

//------------------------------------------------------------------------------

struct TestPlots
{
  TH1 *fPTe0; 
  TH1 *fPTm0; 
  TH1 *fPTe1; 
  TH1 *fPTm1; 
  TH1 *fEtmiss; 
  TH1 *fMee; 
  TH1 *fMmumu; 
  TH1 *fMemu; 
};

//------------------------------------------------------------------------------

class ExRootResult;
class ExRootTreeReader;

//------------------------------------------------------------------------------

void BookHistograms(ExRootResult *result, TestPlots *plots)
{
  TLegend *legend;
  TPaveText *comment;


  plots->fPTe0= result-> AddHist1D("PTe0","PT^{e0}","P_{T}^{e0} GeV", "number o
f electron", 100, 0.0,300.); 

  plots->fPTm0= result-> AddHist1D("PTm0","PT^{mu0}","P_{T}^{mu0} GeV", "number of muon ", 100, 0.0,300.); 




  plots->fPTe1= result-> AddHist1D("PTe1","PT^{e1}","P_{T}^{e1} GeV", "number o
f electron", 100, 0.0,300.); 

  plots->fPTm1= result-> AddHist1D("PTm1","PT^{mu1}","P_{T}^{mu1} GeV", "number of muon ", 100, 0.0,300.); 


  plots->fEtmiss= result-> AddHist1D("Etmiss","Etmiss","Etmiss GeV", "number of muon ", 100, 0.0,500.); 
  plots->fMee= result-> AddHist1D("Mee","Mee","Mee GeV", "number ", 100, 0.0,500.); 
  plots->fMmumu= result-> AddHist1D("Mmumu","Mmumu","Mmumu GeV", "number ", 100, 0.0,500.); 
  plots->fMemu= result-> AddHist1D("Memu","Memu","Memu GeV", "number ", 100, 0.0,500.); 



}

//------------------------------------------------------------------------------

void AnalyseEvents(ExRootTreeReader *treeReader, TestPlots *plots)
{
  TClonesArray *branchElectron = treeReader->UseBranch("Electron");
  TClonesArray *branchMuon = treeReader->UseBranch("Muon");
  TClonesArray *branchMissingET = treeReader->UseBranch("MissingET");

  TFile *ftan = new TFile("small.root","RECREATE"); 
  TTree *tw = new TTree("wprim","parameter for wprim");  
  Int_t xnmu, xnele; 
  Double_t xpte0,xpte1; 
  Double_t xptm0,xptm1; 
  Double_t xetmiss; 
  Double_t xMmumu, xMee, xMemu; 
  tw->Branch("nele",&xnele); 
  tw->Branch("nmu",&xnmu); 
  tw->Branch("pte0",&xpte0); 
  tw->Branch("ptm0",&xptm0); 
  tw->Branch("pte1",&xpte1); 
  tw->Branch("ptm1",&xptm1); 
  tw->Branch("etmiss",&xetmiss); 
  tw->Branch("Mmumu",&xMmumu); 
  tw->Branch("Memu",&xMemu); 
  tw->Branch("Mmee",&xMee); 
  Long64_t allEntries = treeReader->GetEntries();
  
  cout << "** Chain contains " << allEntries << " events" << endl;
  
  TLorentzVector EMom0, EMom1;   
  TLorentzVector MuMom0, MuMom1;   

  Long64_t entry;   
  Int_t i, j;
  
  // Loop over all events
  for(entry = 0; entry < allEntries; ++entry)
    //for(entry = 0; entry < 100; ++entry)
    {
      // Load selected branches with data from specified event
      treeReader->ReadEntry(entry);
      //  cout<<"--  New event -- "<<endl;
      
      //counter of the electron/muon momentum bigger than 20 GeV
      xnmu=0; 
      xnele=0; 
      xpte0=0; 
      xptm0=0; 
      xpte1=0; 
      xptm1=0; 
      xetmiss=0;  
      if(entry%1000 == 0) cout << "Event number: "<< entry <<endl;
      
      Int_t nele=branchElectron-> GetEntriesFast(); 
      Int_t nmu=branchMuon-> GetEntriesFast(); 
      miss= (MissingET*) branchMissingET->At(0); 
      xetmiss=miss->MET; 
      plots-> fEtmiss-> Fill(xetmiss);
      
      //loop over electron 
      for(Int_t i=0; i<nele; ++i){
	ele=  (Electron*) branchElectron->At(i);
        cout<<"electron"<<endl; 
	Int_t echarge=ele->Charge;  
        cout<<ele->PT<<" "<<echarge <<endl; 
	if(i==0){
	  //This is how to use TLorentzVector  
	  EMom0=ele->P4();
	  //example to get charge information 
	  
          xpte0=EMom0.Pt(); 
	  cout<<xpte0<<" "<<ele->PT<<endl;  
	  plots-> fPTe0-> Fill(ele->PT);
	}
	if(i==1){
	  EMom1=ele->P4(); 
          xpte1=EMom1.Pt();
          cout<<xpte1<<" "<<ele->PT<<endl;  
	  plots-> fPTe1-> Fill(ele->PT);
	}
	if(ele->PT>20) xnele++;
      }
	//loop over muon 
      for(Int_t i=0; i<nmu; ++i){
	mu=  (Muon*) branchMuon->At(i);
	if(i==0)
	  {MuMom0=mu->P4();     
	    xptm0=MuMom0.Pt(); 
	    plots->fPTm0->Fill(mu->PT);
	  }
	if(i==1)
	  {MuMom1=mu->P4();  
	    xptm1=MuMom1.Pt();    
	    plots->fPTm1->Fill(mu->PT);
	  }
	if(mu->PT>20.) xnmu++;
      }
// inner products
      TLorentzVector dummy;
      xMmumu=0; 
      xMee=0; 
      xMemu=0; 
      if(nmu>0&& nele>0)
	{
	  dummy=EMom0+MuMom0; 
          xMemu=dummy.M(); 
	  plots->fMemu->Fill(xMemu);
	}
      
      if(nele>1)
	{
	  dummy=EMom0+EMom1; 
	  xMee=dummy.M();
	  plots->fMee->Fill(xMee);
	}
      
      if(nmu>1)
	{
	  dummy=MuMom0+MuMom1; 
	  xMmumu=dummy.M();
	  plots->fMmumu->Fill(xMmumu);
	}
      
      //place to fill the data to the output root file. call only once for a
      //evnet! 
      tw->Fill();
    } //loop over all events
  ftan->Write(); 

  
  
  //counts final number of events 
  
  
}
  

//------------------------------------------------------------------------------
//fill only relevant data to speedup analysis. 
void small2l(void)
{
  const char *outputFile="test.root"; 
  gSystem->Load("~/delphes_test/delphes2/libDelphes");

  TChain *chain = new TChain("Delphes");
  //chain->Add(inputFile);
  TString inputFileList="smallroot_input"; 
  std::ifstream infile( inputFileList); 
  if(!infile.is_open())
    { cout<<"file is not open"<<endl; 
      exit(); 
    }
  if(infile.is_open()){
    while(1){
      if(!infile.good()) break; 
      std::string buffer;
      infile>>buffer;  
      chain->Add(buffer.c_str()); 
      cout<<buffer<<endl;
    } 
  }

  ExRootTreeReader *treeReader = new ExRootTreeReader(chain);
  ExRootResult *result = new ExRootResult();

  TestPlots *plots = new TestPlots;

  BookHistograms(result, plots);

  AnalyseEvents(treeReader, plots);

  result->Write(outputFile);
  TCanvas *c1=new TCanvas("c1","test",841,594); 
  c1->Divide(2,2,0.0025,0.0025); 
  TVirtualPad *cutpad=c1-> cd(1);
  plots->fPTe0->Draw(); 
  c1->cd(2); 
  //plots->fPTm0->Draw();
  //c1->cd(3); 
  //  plots->fPTe1->Draw(); 
  //c1->cd(4);
  //plots->fPTm1->Draw(); 
  //c1->cd(2); 
  //plots-> fEtmiss->Draw(); 
  c1->cd(2); 
  plots-> fMee->Draw(); 
  c1->cd(3); 
  plots-> fMemu->Draw(); 
  c1->cd(4); 
  plots-> fMmumu->Draw(); 

  c1->Print("lepton.jpg"); 


  cout << "** Exiting..." << endl;

  //enable below if you want to clean up results 
  //delete plots;
  //delete result;
  delete treeReader;
  delete chain;
}

//------------------------------------------------------------------------------