// ----------------------------------------------------------------------------
// $Id: 
// ----------------------------------------------------------------------------
//
//  Example file to handle CalHit.
//  
//  In this file, 1) read root event (output from GEANT4)
//                2) access CalHit
//                3) get E Theta Phi and R
//
//                                   Jun. 2001   M. Iwasaki
//
// :: How to run ::
// 
// 1) to call the root, type   
//    
//     root
//
// 2)
//     root [0]  .x rootg4_SeeCalHit.C
// 
// .. that's it.
//
// For example, if you like to run 3 events, type like
//     root [0]  .x rootg4_SeeCalHit(3)
//

#include <stdlib.h>

int rootg4_SeeCalHit(int nEvent=5) {
  gROOT->Reset();

  TString parfile_dir;
  parfile_dir += gSystem->Getenv("LCDROOT")   ;parfile_dir += "/";
  parfile_dir += "ParFiles/";

  // In order to refer to shareable libraries in this simple form,
  // include their paths in a .rootrc file.
  gSystem->Load("libLCDEvent");
  gSystem->Load("libLCDRootAppsUtil");
  //  gSystem->Load("libLCDFastMC");
  //  gSystem->Load("libLCDFullRecon");
  //  gSystem->Load("libLCDPhUtil");

  LCDEvent *event=new LCDEvent();

  // Input root file (GEANT4 output)
  Char_t* rootFile  = "test_geant4_SD.root";
  LCDreadRootFile source(rootFile,event);

  // Must use the same detector name as Full simulation
  Char_t* geomFile = 
    "SDMar01.par"; // SDMar01
  //"LDMar01.par"; // LDMar01
  TString s_geomFile(parfile_dir);
  s_geomFile += geomFile;
  FILE* m_geomfile =  fopen(s_geomFile.Data(),"r");
  LCDGetParameters gp(m_geomfile);
  fclose(m_geomfile);

  // Setup CalHitUtil
  LCDCalHitUtil chutil(&gp);

  // Define a histogram or two
  TH1F *EM_E    = new TH1F("EM_E"    ,"EM energy",50, 0.0, 2.5);
  TH1F *EM_phi  = new TH1F("EM_phi"  ,"EM phi"   ,50, 0.0, 6.3);
  TH1F *EM_theta= new TH1F("EM_theta","EM theta" ,50, -1., 1.0);
  TH1F *EM_r    = new TH1F("EM_r"    ,"EM r"     ,50, 0.0, 300.0);

  // Event loop
  Int_t iEvent;
  for (iEvent = 0; iEvent < nEvent; iEvent++) {

    if (!source.GetEvent()) break; // Get an event from the input

    LCDCalHit*  cE;
    Int_t       sysid;
    Double_t     Ecal, phi, theta, Rcyl;

    Int_t nCal = event->EMCAL()->GetEntries(); // Get EMCAL Hits.
    if (nCal > 0) {
      for (int i_em=0; i_em<nCal; i_em++) {
	cE = (LCDCalHit*)event->EMCAL()->At(i_em);

	Ecal  = chutil.GetEtot(cE);
	phi   = chutil.GetPhi(cE);
	theta = TMath::Cos(chutil.GetTheta(cE));	
	Rcyl  = chutil.GetRcyl(cE);

	sysid=cE->GetTower()->GetSystem();

	switch(sysid) {
	case 10: // EM Barrel
	  EM_E    ->Fill(Ecal);
	  EM_phi  ->Fill(phi);
	  EM_theta->Fill(theta);
	  EM_r    ->Fill(Rcyl);
	  break;
	case 11: // EM EndCap
	  EM_E    ->Fill(Ecal);
	  EM_phi  ->Fill(phi);
	  EM_theta->Fill(theta);
	  EM_r    ->Fill(Rcyl);
	  break;
	default:
	  break;
	}
      }
    }    
  }

  EM_E->SetFillColor(4);
  EM_phi->SetFillColor(4);
  EM_theta->SetFillColor(4);
  EM_r->SetFillColor(4);

  TCanvas *c1 = new TCanvas("c1", "c1", 600, 600);
  c1->Divide(2,2);
  c1->cd(1); EM_E->Draw("B");
  c1->cd(2); EM_phi->Draw("B");
  c1->cd(3); EM_theta->Draw("B");
  c1->cd(4); EM_r->Draw("B");

  return iEvent;
}