2 /**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
17 //_________________________________________________________________________
18 // An analysis task to check the trigger data in ESD with MC data
19 // Creates an ntuple for 2x2 and NxN triggers
20 // Each ntuple connects the maximum trigger amplitudes
21 // and its positions with reconstructed clusters and MC
23 //*-- Yves Schutz (CERN) & Gustavo Conesa Balbastre (INFN-LNF)
24 //////////////////////////////////////////////////////////////////////////////
31 #include "AliAnalysisManager.h"
32 #include "AliMCEventHandler.h"
33 #include "AliMCEvent.h"
34 #include "AliAnaCaloTriggerMC.h"
35 #include "AliESDEvent.h"
36 #include "AliESDCaloCluster.h"
39 #include "TParticle.h"
41 //______________________________________________________________________________
42 AliAnaCaloTriggerMC::AliAnaCaloTriggerMC() :
51 // Default constructor.
54 //______________________________________________________________________________
55 AliAnaCaloTriggerMC::AliAnaCaloTriggerMC(const char *name) :
56 AliAnalysisTask(name, "AnaCaloTriggerMC"),
66 // Input slot #0 works with an Ntuple
67 DefineInput(0, TChain::Class());
68 // Output slot #0 writes into a TH1 container
69 DefineOutput(0, TObjArray::Class()) ;
71 //____________________________________________________________________________
72 AliAnaCaloTriggerMC::AliAnaCaloTriggerMC(const AliAnaCaloTriggerMC & ct) :
73 AliAnalysisTask(ct),fChain(ct.fChain), fESD(ct.fESD),
74 fOutputContainer(ct.fOutputContainer), fCalorimeter(ct. fCalorimeter),
75 fNtTrigger22(ct.fNtTrigger22), fNtTriggerNN(ct.fNtTriggerNN)
79 SetName (ct.GetName()) ;
80 SetTitle(ct.GetTitle()) ;
84 //_________________________________________________________________________
85 AliAnaCaloTriggerMC & AliAnaCaloTriggerMC::operator = (const AliAnaCaloTriggerMC & source)
87 // assignment operator
89 if(&source == this) return *this;
91 fChain = source.fChain ;
93 fOutputContainer = source.fOutputContainer ;
94 fCalorimeter = source. fCalorimeter ;
95 fNtTrigger22 = source.fNtTrigger22 ;
96 fNtTriggerNN = source.fNtTriggerNN ;
102 //______________________________________________________________________________
103 AliAnaCaloTriggerMC::~AliAnaCaloTriggerMC()
106 // fOutputContainer->Clear() ;
107 // delete fOutputContainer ;
112 //______________________________________________________________________________
113 void AliAnaCaloTriggerMC::ConnectInputData(const Option_t*)
115 // Initialisation of branch container and histograms
117 AliInfo(Form("*** Initialization of %s", GetName())) ;
120 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
122 AliError(Form("Input 0 for %s not found\n", GetName()));
126 fESD = new AliESDEvent();
127 fESD->ReadFromTree(fChain);
131 //________________________________________________________________________
133 void AliAnaCaloTriggerMC::CreateOutputObjects()
136 // Create the output container
140 fNtTrigger22 = new TNtuple(fCalorimeter+"trigger22", "Trigger data 2x2 patch", "a22:a220:enMax:phEnMax:eta22:phi22:etaMax:phiMax:phEtaMax:phPhiMax");
141 fNtTriggerNN = new TNtuple(fCalorimeter+"triggerNN", "Trigger data NxN patch", "aNN:aNN0:enMax:phEnMax:etaNN:phiNN:etaMax:phiMax:phEtaMax:phPhiMax");
143 // create output container
145 fOutputContainer = new TObjArray(2) ;
146 fOutputContainer->SetName(GetName()) ;
148 fOutputContainer->AddAt(fNtTrigger22, 0) ;
149 fOutputContainer->AddAt(fNtTriggerNN, 1) ;
153 //______________________________________________________________________________
154 void AliAnaCaloTriggerMC::Exec(Option_t *)
157 // Processing of one event
159 Long64_t entry = fChain->GetReadEntry() ;
162 AliError("fESD is not connected to the input!") ;
166 if ( !((entry-1)%100) )
167 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
170 AliStack* stack = 0x0;
171 AliMCEventHandler* mctruth = (AliMCEventHandler*)
172 ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
175 stack = mctruth->MCEvent()->Stack();
178 AliError("Stack not found") ;
182 // Get trigger information of fCalorimeter
183 TArrayF * triggerAmplitudes = 0x0 ;
184 TArrayF * triggerPosition = 0x0 ;
185 Int_t numberOfCaloClusters = 0 ;
187 if(fCalorimeter == "PHOS"){
188 triggerAmplitudes = fESD->GetPHOSTriggerAmplitudes();
189 triggerPosition = fESD->GetPHOSTriggerPosition();
191 else if(fCalorimeter == "EMCAL"){
192 triggerAmplitudes = fESD->GetEMCALTriggerAmplitudes();
193 triggerPosition = fESD->GetEMCALTriggerPosition();
196 if( triggerAmplitudes && triggerPosition ){
197 // trigger amplitudes
198 const Float_t ka22 = static_cast<Float_t>(triggerAmplitudes->At(0)) ;
199 const Float_t ka22O = static_cast<Float_t>(triggerAmplitudes->At(1)) ;
200 const Float_t kaNN = static_cast<Float_t>(triggerAmplitudes->At(2)) ;
201 const Float_t kaNNO = static_cast<Float_t>(triggerAmplitudes->At(3)) ;
204 const Float_t kx22 = static_cast<Float_t>(triggerPosition->At(0)) ;
205 const Float_t ky22 = static_cast<Float_t>(triggerPosition->At(1)) ;
206 const Float_t kz22 = static_cast<Float_t>(triggerPosition->At(2)) ;
207 const Float_t kxNN = static_cast<Float_t>(triggerPosition->At(3)) ;
208 const Float_t kyNN = static_cast<Float_t>(triggerPosition->At(4)) ;
209 const Float_t kzNN = static_cast<Float_t>(triggerPosition->At(5)) ;
212 Float_t phEnMax = 0. ;
213 Float_t etaMax = 0.5 ;
214 Float_t phiMax = 0. ;
215 Float_t phEtaMax = 0.5 ;
216 Float_t phPhiMax = 0. ;
218 TVector3 vpos22(kx22, ky22, kz22) ;
219 TVector3 vposNN(kxNN, kyNN, kzNN) ;
220 Float_t eta22 = vpos22.Eta() ;
221 Float_t phi22 = vpos22.Phi() * TMath::RadToDeg() + 360. ;
222 Float_t etaNN = vposNN.Eta() ;
223 Float_t phiNN = vposNN.Phi() * TMath::RadToDeg() + 360. ;
227 // loop over the Calorimeters Clusters
228 Float_t cluEnergy = 0;
230 for(icaloCluster = 0 ; icaloCluster < numberOfCaloClusters ; icaloCluster++) {
232 AliESDCaloCluster * cluster = fESD->GetCaloCluster(icaloCluster) ;
234 if (cluster && ( (fCalorimeter == "PHOS" && cluster->IsPHOS()) ||
235 (fCalorimeter == "EMCAL" && cluster->IsEMCAL()))) {
237 cluEnergy = cluster->E() ;
241 cluster->GetPosition( pos ) ;
243 if ( cluEnergy > enMax) {
245 vpos.SetXYZ(pos[0], pos[1], pos[2]) ;
246 etaMax = vpos.Eta() ;
247 phiMax = vpos.Phi() ;
248 labelmax = cluster->GetLabel();
251 Double_t * pid = cluster->GetPid() ;
253 if(pid[AliPID::kPhoton] > 0.9) {
254 if ( cluEnergy > phEnMax) {
255 phEnMax = cluEnergy ;
256 vpos.SetXYZ(pos[0], pos[1], pos[2]) ;
257 phEtaMax = vpos.Eta() ;
258 phPhiMax = vpos.Phi() ;
264 if(labelmax < stack->GetNtrack() && labelmax >= 0 ){
265 TParticle * particle = stack->Particle(labelmax);
266 Float_t ptgen = particle->Energy();
267 fNtTrigger22->Fill(ka22, ka22O, ptgen, enMax, phEnMax, eta22, phi22, etaMax, phiMax * TMath::RadToDeg() + 360., phEtaMax, phPhiMax * TMath::RadToDeg() + 360.);
268 fNtTriggerNN->Fill(kaNN, kaNNO, ptgen, enMax, phEnMax, etaNN, phiNN, etaMax, phiMax * TMath::RadToDeg() + 360., phEtaMax, phPhiMax * TMath::RadToDeg() + 360.);
270 else AliDebug(1, Form("Wrong label %d, ntrack %d, Emax %f ",labelmax, stack->GetNtrack(), phEnMax));
271 }//If trigger arrays filled
273 PostData(0, fOutputContainer);
278 //______________________________________________________________________________
279 void AliAnaCaloTriggerMC::Terminate(Option_t *)
281 // Processing when the event loop is ended