[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskTrgContam.cxx

// $Id$
//
// Task to compute the trigger contamination by exotics.
//
// Authors: M.Cosentino

#include "TChain.h"
#include "TTree.h"
#include "TH1F.h"
#include "TH2F.h"
#include "TCanvas.h"

#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"
#include "AliESDEvent.h"
#include "AliESDHeader.h"
#include "AliESDUtils.h"
#include "AliESDInputHandler.h"
#include "AliESDpid.h"
#include "AliKFParticle.h"
#include "AliMCEventHandler.h"
#include "AliMCEvent.h"
#include "AliStack.h"
#include "AliESDtrackCuts.h"
#include "AliESDv0.h"
#include "AliV0vertexer.h"
#include "AliESDCaloCluster.h"
#include "AliESDCaloCells.h"
#include "AliEMCALGeometry.h"
#include "AliEMCALRecoUtils.h"
#include "TLorentzVector.h"
#include "AliVCluster.h"
#include "AliAnalysisTaskTrgContam.h"
#include "TFile.h"

ClassImp(AliAnalysisTaskTrgContam)

AliAnalysisTaskTrgContam::AliAnalysisTaskTrgContam() : 
  AliAnalysisTaskSE(), 
  fCaloClusters(0),
  fEMCalCells(0),
  fGeom(0x0),
  fGeoName("EMCAL_COMPLETEV1"),
  fPeriod("LHC11c"),
  fIsTrain(0),
  fTrigThresh(4.8),
  fExoticCut(0.97),
  fESD(0),
  fOutputList(0),
  fEvtSel(0),
  fClusEt(0),
  fClusEtTM(0),
  fClusEtLead(0),
  fClusEtSubLead(0),
  fClusEtLeadTM(0),
  fClusEtSubLeadTM(0),
  fClusEtExotic(0), 
  fClusEtExoticTM(0),
  fClusEtSingleExotic(0),
  fCellEnergy(0),
  fM02Et(0),
  fM02EtTM(0),
  fM02EtExot(0),
  fM02EtExotTM(0)
{
  // Default constructor.
}

//________________________________________________________________________
AliAnalysisTaskTrgContam::AliAnalysisTaskTrgContam(const char *name) : 
  AliAnalysisTaskSE(name), 
  fCaloClusters(0),
  fEMCalCells(0),
  fGeom(0x0),
  fGeoName("EMCAL_COMPLETEV1"),
  fPeriod("LHC11c"),
  fIsTrain(0),
  fTrigThresh(4.8),
  fExoticCut(0.97),
  fESD(0),
  fOutputList(0),
  fEvtSel(0),
  fClusEt(0),
  fClusEtTM(0),
  fClusEtLead(0),
  fClusEtSubLead(0),
  fClusEtLeadTM(0),
  fClusEtSubLeadTM(0),
  fClusEtExotic(0),
  fClusEtExoticTM(0), 
  fClusEtSingleExotic(0),
  fCellEnergy(0),
  fM02Et(0),
  fM02EtTM(0),
  fM02EtExot(0),
  fM02EtExotTM(0)
{
  // Constructor

  // Define input and output slots here
  // Input slot #0 works with a TChain
  DefineInput(0, TChain::Class());
  // Output slot #0 id reserved by the base class for AOD
  // Output slot #1 writes into a TH1 container
  DefineOutput(1, TList::Class());
}

//________________________________________________________________________
void AliAnalysisTaskTrgContam::UserCreateOutputObjects()
{
  // Create histograms, called once.
    
  fCaloClusters = new TRefArray();
  
  fOutputList = new TList();
  fOutputList->SetOwner();// Container cleans up all histos (avoids leaks in merging) 
  
  fGeom = AliEMCALGeometry::GetInstance(fGeoName.Data());
  
  fEvtSel = new TH1F("hEvtSel","Event selection counter (0=all trg, 1=pvz cut) ;evt cut ;dN/dcut}",2,0,2);
  fOutputList->Add(fEvtSel);
  
  fClusEt = new TH1F("hClusEt","Clusters E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEt);
  
  fClusEtTM = new TH1F("hClusEtTM","Clusters (track-matched) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtTM);
  
  fClusEtLead = new TH1F("hClusEtLead","Clusters (leading-trig) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtLead);
  
  fClusEtSubLead = new TH1F("hClusEtSubLead","Clusters (subleading-trig) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtSubLead);
  
  fClusEtLeadTM = new TH1F("hClusEtLeadTM","Clusters (leading-trig, TM) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtLeadTM);
  
  fClusEtSubLeadTM = new TH1F("hClusEtSubLeadTM","Clusters (subleading-trig, TM) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtSubLeadTM);
  
  fClusEtExotic = new TH1F("hClusEtExotic","Exotic trigger clusters  E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtExotic);
  
  fClusEtExoticTM = new TH1F("hClusEtExoticTM","Exotic trigger clusters (TM)  E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtExoticTM);
  
  fClusEtSingleExotic = new TH1F("hClusEtSingleExotic","Exotic trigger clusters (only this above thrs.)  E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
  fOutputList->Add(fClusEtSingleExotic);

  fCellEnergy = new TH1F("hCellE","cell energy spectrum;E_{cell} (GeV);entries",200,0,20);
  fOutputList->Add(fCellEnergy);
  
  fM02Et = new TH2F("hM02Et","#lambda_{0}^{2} vs. E_{T} for trigger clusters ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
  fOutputList->Add(fM02Et);
  
  fM02EtTM = new TH2F("hM02EtTM","#lambda_{0}^{2} vs. E_{T} for trigger clusters(TM) ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
  fOutputList->Add(fM02EtTM);
  
  fM02EtExot = new TH2F("hM02EtExot","#lambda_{0}^{2} vs. E_{T} for trigger clusters(Exotic) ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
  fOutputList->Add(fM02EtExot);

  fM02EtExotTM = new TH2F("hM02EtExotTM","#lambda_{0}^{2} vs. E_{T} for trigger clusters(TM+Exotic) ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
  fOutputList->Add(fM02EtExotTM);
  
  PostData(1, fOutputList);
}

//________________________________________________________________________
void AliAnalysisTaskTrgContam::UserExec(Option_t *) 
{
  // User exec. Called once per event.

  Bool_t isSelected = 0;
  if(fPeriod.Contains("11a"))
    isSelected =  (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kEMC1);
  else
    isSelected =  ((((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kCentral) ||
		   (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kSemiCentral));
  if(!isSelected )
    return; 

  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
  if (!fESD) {
    printf("ERROR: fESD not available\n");
    return;
  }

  fEvtSel->Fill(0);
  if (0) {
    AliESDVertex *pv = (AliESDVertex*)fESD->GetPrimaryVertex();
    if(!pv) 
      return;
    if(TMath::Abs(pv->GetZ())>15)
      return;
  }
  fEvtSel->Fill(1);
  if (!fIsTrain) {
    for(Int_t mod=0; mod < (fGeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
      if(fGeoName=="EMCAL_FIRSTYEARV1" && mod>3)
        break;
      fGeom->SetMisalMatrix(fESD->GetEMCALMatrix(mod), mod);
    }
  }
  fESD->GetEMCALClusters(fCaloClusters);
  fEMCalCells = fESD->GetEMCALCells();
  for(int i=0;i<fEMCalCells->GetNumberOfCells();i++){
    Double_t e = fEMCalCells->GetCellAmplitude(TMath::Abs(fEMCalCells->GetAmplitude(i)));
    fCellEnergy->Fill(e);
  }
  FillClusHists(); 
  
  fCaloClusters->Clear();
  PostData(1, fOutputList);
}
 
//________________________________________________________________________
void AliAnalysisTaskTrgContam::FillClusHists()
{
  // Fill cluster histograms.

  if(!fCaloClusters)
    return;
  const Int_t nclus = fCaloClusters->GetEntries();
  if(nclus==0)
    return;
  Double_t EtArray[nclus];
  Bool_t isTM[nclus];
  Bool_t isEx[nclus];
  Int_t index[nclus];
  Int_t nthresholds = 0;
  for(Int_t ic=0;ic<nclus;ic++){
    EtArray[ic]=0;
    isTM[ic] = 0;
    isEx[ic] = 0;
    index[ic]=0;
    AliESDCaloCluster *c = static_cast<AliESDCaloCluster*>(fCaloClusters->At(ic));
    if(!c)
      continue;
    if(!c->IsEMCAL())
      continue;
    if(c->E()<fTrigThresh)
      continue;
    nthresholds++;
    Short_t id;
    Double_t Emax = GetMaxCellEnergy( c, id);
    Double_t Ecross = GetCrossEnergy( c, id);
    if((1-Ecross/Emax)>fExoticCut)
      isEx[ic] = 1;
    Float_t clsPos[3] = {0,0,0};
    c->GetPosition(clsPos);
    TVector3 clsVec(clsPos);
    Double_t Et = c->E()*TMath::Sin(clsVec.Theta());
    EtArray[ic] = Et;
    fClusEt->Fill(Et);
    fM02Et->Fill(Et, c->GetM02());
    if(isEx[ic]){
      fClusEtExotic->Fill(Et);
      fM02EtExot->Fill(Et,c->GetM02()); 
     }
    Double_t dR = TMath::Sqrt(pow(c->GetTrackDx(),2)+pow(c->GetTrackDz(),2));
    if(dR<0.025){
      isTM[ic]=1;
      fClusEtTM->Fill(Et);
      fM02EtTM->Fill(Et, c->GetM02());
      if(isEx[ic]){
	fClusEtExoticTM->Fill(Et);
	fM02EtExotTM->Fill(Et,c->GetM02());
      }
    }
  }
  TMath::Sort(nclus,EtArray,index, kTRUE);
  if(EtArray[index[0]]>0){
    fClusEtLead->Fill(EtArray[index[0]]);
    if(nthresholds==1 && isEx[index[0]])
       fClusEtSingleExotic->Fill(EtArray[index[0]]);
  }
  if(nclus>1)if(EtArray[index[1]]>0)
    fClusEtSubLead->Fill(EtArray[index[1]]);
  if(isTM[index[0]] && EtArray[index[0]]>0)
    fClusEtLeadTM->Fill(EtArray[index[0]]);
  if(nclus>1)if(isTM[index[1]] && EtArray[index[1]]>0)
    fClusEtSubLeadTM->Fill(EtArray[index[1]]);
} 

//________________________________________________________________________
Double_t AliAnalysisTaskTrgContam::GetCrossEnergy(const AliVCluster *cluster, Short_t &idmax)
{
  // Calculate the energy of cross cells around the leading cell.

  AliVCaloCells *cells = 0;
  cells = fESD->GetEMCALCells();
  if (!cells)
    return 0;

  if (!fGeom)
    return 0;

  Int_t iSupMod = -1;
  Int_t iTower  = -1;
  Int_t iIphi   = -1;
  Int_t iIeta   = -1;
  Int_t iphi    = -1;
  Int_t ieta    = -1;
  Int_t iphis   = -1;
  Int_t ietas   = -1;

  Double_t crossEnergy = 0;

  fGeom->GetCellIndex(idmax,iSupMod,iTower,iIphi,iIeta);
  fGeom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphis,ietas);

  Int_t ncells = cluster->GetNCells();
  for (Int_t i=0; i<ncells; i++) {
    Int_t cellAbsId = cluster->GetCellAbsId(i);
    fGeom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
    fGeom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
    Int_t aphidiff = TMath::Abs(iphi-iphis);
    if (aphidiff>1)
      continue;
    Int_t aetadiff = TMath::Abs(ieta-ietas);
    if (aetadiff>1)
      continue;
    if ( (aphidiff==1 && aetadiff==0) ||
	(aphidiff==0 && aetadiff==1) ) {
      crossEnergy += cells->GetCellAmplitude(cellAbsId);
    }
  }

  return crossEnergy;
}

//________________________________________________________________________
Double_t AliAnalysisTaskTrgContam ::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const
{
  // Get maximum energy of attached cell.

  id = -1;

  AliVCaloCells *cells = 0;
  cells = fESD->GetEMCALCells();
  if (!cells)
    return 0;

  Double_t maxe = 0;
  Int_t ncells = cluster->GetNCells();
  for (Int_t i=0; i<ncells; i++) {
    Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
    if (e>maxe) {
      maxe = e;
      id   = cluster->GetCellAbsId(i);
    }
  }
  return maxe;
}

//________________________________________________________________________
void AliAnalysisTaskTrgContam::Terminate(Option_t *) 
{
  // Called once at the end of the query
}
Commit	Line	Data
6176cf5e	1	// $Id$
70d53162	2	//
	3	// Task to compute the trigger contamination by exotics.
	4	//
	5	// Authors: M.Cosentino
6176cf5e	6
eb7f0edc	7	#include "TChain.h"
	8	#include "TTree.h"
	9	#include "TH1F.h"
	10	#include "TH2F.h"
	11	#include "TCanvas.h"
c363681d	12
eb7f0edc	13	#include "AliAnalysisTask.h"
eb7f0edc	14	#include "AliAnalysisManager.h"
eb7f0edc	15	#include "AliESDEvent.h"
	16	#include "AliESDHeader.h"
	17	#include "AliESDUtils.h"
	18	#include "AliESDInputHandler.h"
	19	#include "AliESDpid.h"
	20	#include "AliKFParticle.h"
eb7f0edc	21	#include "AliMCEventHandler.h"
	22	#include "AliMCEvent.h"
	23	#include "AliStack.h"
eb7f0edc	24	#include "AliESDtrackCuts.h"
	25	#include "AliESDv0.h"
	26	#include "AliV0vertexer.h"
	27	#include "AliESDCaloCluster.h"
	28	#include "AliESDCaloCells.h"
	29	#include "AliEMCALGeometry.h"
	30	#include "AliEMCALRecoUtils.h"
	31	#include "TLorentzVector.h"
	32	#include "AliVCluster.h"
eb7f0edc	33	#include "AliAnalysisTaskTrgContam.h"
	34	#include "TFile.h"
	35
eb7f0edc	36	ClassImp(AliAnalysisTaskTrgContam)
eb7f0edc	37
6176cf5e	38	AliAnalysisTaskTrgContam::AliAnalysisTaskTrgContam() :
	39	AliAnalysisTaskSE(),
	40	fCaloClusters(0),
	41	fEMCalCells(0),
f912f9a9	42	fGeom(0x0),
6176cf5e	43	fGeoName("EMCAL_COMPLETEV1"),
	44	fPeriod("LHC11c"),
	45	fIsTrain(0),
	46	fTrigThresh(4.8),
	47	fExoticCut(0.97),
6176cf5e	48	fESD(0),
	49	fOutputList(0),
	50	fEvtSel(0),
	51	fClusEt(0),
	52	fClusEtTM(0),
	53	fClusEtLead(0),
	54	fClusEtSubLead(0),
	55	fClusEtLeadTM(0),
	56	fClusEtSubLeadTM(0),
	57	fClusEtExotic(0),
	58	fClusEtExoticTM(0),
	59	fClusEtSingleExotic(0),
	60	fCellEnergy(0),
	61	fM02Et(0),
	62	fM02EtTM(0),
	63	fM02EtExot(0),
	64	fM02EtExotTM(0)
	65	{
	66	// Default constructor.
	67	}
	68
eb7f0edc	69	//________________________________________________________________________
6176cf5e	70	AliAnalysisTaskTrgContam::AliAnalysisTaskTrgContam(const char *name) :
6176cf5e	71	AliAnalysisTaskSE(name),
eb7f0edc	72	fCaloClusters(0),
db0b19b8	73	fEMCalCells(0),
eb7f0edc	74	fGeom(0x0),
	75	fGeoName("EMCAL_COMPLETEV1"),
	76	fPeriod("LHC11c"),
	77	fIsTrain(0),
	78	fTrigThresh(4.8),
	79	fExoticCut(0.97),
eb7f0edc	80	fESD(0),
eb7f0edc	81	fOutputList(0),
eb7f0edc	82	fEvtSel(0),
eb7f0edc	83	fClusEt(0),
	84	fClusEtTM(0),
	85	fClusEtLead(0),
	86	fClusEtSubLead(0),
	87	fClusEtLeadTM(0),
	88	fClusEtSubLeadTM(0),
	89	fClusEtExotic(0),
	90	fClusEtExoticTM(0),
	91	fClusEtSingleExotic(0),
c363681d	92	fCellEnergy(0),
eb7f0edc	93	fM02Et(0),
	94	fM02EtTM(0),
	95	fM02EtExot(0),
	96	fM02EtExotTM(0)
eb7f0edc	97	{
	98	// Constructor
	99
	100	// Define input and output slots here
	101	// Input slot #0 works with a TChain
	102	DefineInput(0, TChain::Class());
	103	// Output slot #0 id reserved by the base class for AOD
	104	// Output slot #1 writes into a TH1 container
	105	DefineOutput(1, TList::Class());
	106	}
6176cf5e	107
eb7f0edc	108	//________________________________________________________________________
	109	void AliAnalysisTaskTrgContam::UserCreateOutputObjects()
	110	{
6176cf5e	111	// Create histograms, called once.
eb7f0edc	112
	113	fCaloClusters = new TRefArray();
	114
	115	fOutputList = new TList();
	116	fOutputList->SetOwner();// Container cleans up all histos (avoids leaks in merging)
	117
f912f9a9	118	fGeom = AliEMCALGeometry::GetInstance(fGeoName.Data());
eb7f0edc	119
	120	fEvtSel = new TH1F("hEvtSel","Event selection counter (0=all trg, 1=pvz cut) ;evt cut ;dN/dcut}",2,0,2);
	121	fOutputList->Add(fEvtSel);
eb7f0edc	122
	123	fClusEt = new TH1F("hClusEt","Clusters E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	124	fOutputList->Add(fClusEt);
	125
	126	fClusEtTM = new TH1F("hClusEtTM","Clusters (track-matched) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	127	fOutputList->Add(fClusEtTM);
	128
	129	fClusEtLead = new TH1F("hClusEtLead","Clusters (leading-trig) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	130	fOutputList->Add(fClusEtLead);
	131
	132	fClusEtSubLead = new TH1F("hClusEtSubLead","Clusters (subleading-trig) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	133	fOutputList->Add(fClusEtSubLead);
	134
	135	fClusEtLeadTM = new TH1F("hClusEtLeadTM","Clusters (leading-trig, TM) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	136	fOutputList->Add(fClusEtLeadTM);
	137
	138	fClusEtSubLeadTM = new TH1F("hClusEtSubLeadTM","Clusters (subleading-trig, TM) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	139	fOutputList->Add(fClusEtSubLeadTM);
	140
	141	fClusEtExotic = new TH1F("hClusEtExotic","Exotic trigger clusters E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	142	fOutputList->Add(fClusEtExotic);
	143
	144	fClusEtExoticTM = new TH1F("hClusEtExoticTM","Exotic trigger clusters (TM) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	145	fOutputList->Add(fClusEtExoticTM);
	146
	147	fClusEtSingleExotic = new TH1F("hClusEtSingleExotic","Exotic trigger clusters (only this above thrs.) E_{T} ;E_{T} ;dN/dE_{T}",400,0,200);
	148	fOutputList->Add(fClusEtSingleExotic);
c363681d	149
	150	fCellEnergy = new TH1F("hCellE","cell energy spectrum;E_{cell} (GeV);entries",200,0,20);
	151	fOutputList->Add(fCellEnergy);
eb7f0edc	152
	153	fM02Et = new TH2F("hM02Et","#lambda_{0}^{2} vs. E_{T} for trigger clusters ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
	154	fOutputList->Add(fM02Et);
	155
	156	fM02EtTM = new TH2F("hM02EtTM","#lambda_{0}^{2} vs. E_{T} for trigger clusters(TM) ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
	157	fOutputList->Add(fM02EtTM);
	158
	159	fM02EtExot = new TH2F("hM02EtExot","#lambda_{0}^{2} vs. E_{T} for trigger clusters(Exotic) ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
	160	fOutputList->Add(fM02EtExot);
	161
	162	fM02EtExotTM = new TH2F("hM02EtExotTM","#lambda_{0}^{2} vs. E_{T} for trigger clusters(TM+Exotic) ;E_{T} ;#lambda_{0}^{2}",400,0,200, 400,0,4);
	163	fOutputList->Add(fM02EtExotTM);
	164
	165	PostData(1, fOutputList);
	166	}
	167
	168	//________________________________________________________________________
	169	void AliAnalysisTaskTrgContam::UserExec(Option_t *)
	170	{
6176cf5e	171	// User exec. Called once per event.
6176cf5e	172
eb7f0edc	173	Bool_t isSelected = 0;
	174	if(fPeriod.Contains("11a"))
	175	isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kEMC1);
	176	else
c363681d	177	isSelected = ((((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kCentral) \|\|
c363681d	178	(((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kSemiCentral));
eb7f0edc	179	if(!isSelected )
c363681d	180	return;
c363681d	181
eb7f0edc	182	fESD = dynamic_cast<AliESDEvent*>(InputEvent());
6176cf5e	183	if (!fESD) {
eb7f0edc	184	printf("ERROR: fESD not available\n");
	185	return;
	186	}
c363681d	187
6176cf5e	188	fEvtSel->Fill(0);
	189	if (0) {
	190	AliESDVertex pv = (AliESDVertex)fESD->GetPrimaryVertex();
	191	if(!pv)
	192	return;
	193	if(TMath::Abs(pv->GetZ())>15)
	194	return;
	195	}
	196	fEvtSel->Fill(1);
	197	if (!fIsTrain) {
	198	for(Int_t mod=0; mod < (fGeom->GetEMCGeometry())->GetNumberOfSuperModules(); mod++){
	199	if(fGeoName=="EMCAL_FIRSTYEARV1" && mod>3)
	200	break;
	201	fGeom->SetMisalMatrix(fESD->GetEMCALMatrix(mod), mod);
	202	}
	203	}
	204	fESD->GetEMCALClusters(fCaloClusters);
	205	fEMCalCells = fESD->GetEMCALCells();
	206	for(int i=0;i<fEMCalCells->GetNumberOfCells();i++){
	207	Double_t e = fEMCalCells->GetCellAmplitude(TMath::Abs(fEMCalCells->GetAmplitude(i)));
	208	fCellEnergy->Fill(e);
eb7f0edc	209	}
6176cf5e	210	FillClusHists();
	211
	212	fCaloClusters->Clear();
	213	PostData(1, fOutputList);
c363681d	214	}
6176cf5e	215
eb7f0edc	216	//________________________________________________________________________
	217	void AliAnalysisTaskTrgContam::FillClusHists()
	218	{
6176cf5e	219	// Fill cluster histograms.
6176cf5e	220
eb7f0edc	221	if(!fCaloClusters)
	222	return;
	223	const Int_t nclus = fCaloClusters->GetEntries();
	224	if(nclus==0)
	225	return;
c363681d	226	Double_t EtArray[nclus];
	227	Bool_t isTM[nclus];
	228	Bool_t isEx[nclus];
	229	Int_t index[nclus];
eb7f0edc	230	Int_t nthresholds = 0;
	231	for(Int_t ic=0;ic<nclus;ic++){
	232	EtArray[ic]=0;
	233	isTM[ic] = 0;
	234	isEx[ic] = 0;
	235	index[ic]=0;
	236	AliESDCaloCluster c = static_cast<AliESDCaloCluster>(fCaloClusters->At(ic));
	237	if(!c)
	238	continue;
	239	if(!c->IsEMCAL())
	240	continue;
	241	if(c->E()<fTrigThresh)
	242	continue;
	243	nthresholds++;
	244	Short_t id;
	245	Double_t Emax = GetMaxCellEnergy( c, id);
	246	Double_t Ecross = GetCrossEnergy( c, id);
	247	if((1-Ecross/Emax)>fExoticCut)
	248	isEx[ic] = 1;
	249	Float_t clsPos[3] = {0,0,0};
	250	c->GetPosition(clsPos);
	251	TVector3 clsVec(clsPos);
	252	Double_t Et = c->E()*TMath::Sin(clsVec.Theta());
	253	EtArray[ic] = Et;
	254	fClusEt->Fill(Et);
	255	fM02Et->Fill(Et, c->GetM02());
	256	if(isEx[ic]){
	257	fClusEtExotic->Fill(Et);
	258	fM02EtExot->Fill(Et,c->GetM02());
	259	}
	260	Double_t dR = TMath::Sqrt(pow(c->GetTrackDx(),2)+pow(c->GetTrackDz(),2));
	261	if(dR<0.025){
	262	isTM[ic]=1;
	263	fClusEtTM->Fill(Et);
	264	fM02EtTM->Fill(Et, c->GetM02());
	265	if(isEx[ic]){
	266	fClusEtExoticTM->Fill(Et);
	267	fM02EtExotTM->Fill(Et,c->GetM02());
	268	}
	269	}
	270	}
	271	TMath::Sort(nclus,EtArray,index, kTRUE);
	272	if(EtArray[index[0]]>0){
	273	fClusEtLead->Fill(EtArray[index[0]]);
	274	if(nthresholds==1 && isEx[index[0]])
	275	fClusEtSingleExotic->Fill(EtArray[index[0]]);
	276	}
	277	if(nclus>1)if(EtArray[index[1]]>0)
	278	fClusEtSubLead->Fill(EtArray[index[1]]);
	279	if(isTM[index[0]] && EtArray[index[0]]>0)
	280	fClusEtLeadTM->Fill(EtArray[index[0]]);
	281	if(nclus>1)if(isTM[index[1]] && EtArray[index[1]]>0)
	282	fClusEtSubLeadTM->Fill(EtArray[index[1]]);
	283	}
6176cf5e	284
eb7f0edc	285	//________________________________________________________________________
	286	Double_t AliAnalysisTaskTrgContam::GetCrossEnergy(const AliVCluster *cluster, Short_t &idmax)
	287	{
	288	// Calculate the energy of cross cells around the leading cell.
	289
	290	AliVCaloCells *cells = 0;
	291	cells = fESD->GetEMCALCells();
	292	if (!cells)
	293	return 0;
	294
eb7f0edc	295	if (!fGeom)
	296	return 0;
	297
	298	Int_t iSupMod = -1;
	299	Int_t iTower = -1;
	300	Int_t iIphi = -1;
	301	Int_t iIeta = -1;
	302	Int_t iphi = -1;
	303	Int_t ieta = -1;
	304	Int_t iphis = -1;
	305	Int_t ietas = -1;
	306
	307	Double_t crossEnergy = 0;
	308
	309	fGeom->GetCellIndex(idmax,iSupMod,iTower,iIphi,iIeta);
	310	fGeom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphis,ietas);
	311
	312	Int_t ncells = cluster->GetNCells();
	313	for (Int_t i=0; i<ncells; i++) {
	314	Int_t cellAbsId = cluster->GetCellAbsId(i);
	315	fGeom->GetCellIndex(cellAbsId,iSupMod,iTower,iIphi,iIeta);
	316	fGeom->GetCellPhiEtaIndexInSModule(iSupMod,iTower,iIphi, iIeta,iphi,ieta);
	317	Int_t aphidiff = TMath::Abs(iphi-iphis);
	318	if (aphidiff>1)
	319	continue;
	320	Int_t aetadiff = TMath::Abs(ieta-ietas);
	321	if (aetadiff>1)
	322	continue;
	323	if ( (aphidiff==1 && aetadiff==0) \|\|
	324	(aphidiff==0 && aetadiff==1) ) {
	325	crossEnergy += cells->GetCellAmplitude(cellAbsId);
	326	}
	327	}
	328
	329	return crossEnergy;
	330	}
	331
eb7f0edc	332	//________________________________________________________________________
	333	Double_t AliAnalysisTaskTrgContam ::GetMaxCellEnergy(const AliVCluster *cluster, Short_t &id) const
	334	{
	335	// Get maximum energy of attached cell.
	336
	337	id = -1;
	338
	339	AliVCaloCells *cells = 0;
	340	cells = fESD->GetEMCALCells();
	341	if (!cells)
	342	return 0;
	343
	344	Double_t maxe = 0;
	345	Int_t ncells = cluster->GetNCells();
	346	for (Int_t i=0; i<ncells; i++) {
	347	Double_t e = cells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
	348	if (e>maxe) {
	349	maxe = e;
	350	id = cluster->GetCellAbsId(i);
	351	}
	352	}
	353	return maxe;
	354	}
6176cf5e	355
eb7f0edc	356	//________________________________________________________________________
	357	void AliAnalysisTaskTrgContam::Terminate(Option_t *)
	358	{
eb7f0edc	359	// Called once at the end of the query
eb7f0edc	360	}