ESD file is added to the list of proof output files. It is then automatically merged...
[u/mrichter/AliRoot.git] / PWG4 / AliAnaCaloTriggerMC.cxx
CommitLineData
884e9012 1
2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
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 **************************************************************************/
16
17//_________________________________________________________________________
f38b754b 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
884e9012 22//
f38b754b 23//*-- Yves Schutz (CERN) & Gustavo Conesa Balbastre (INFN-LNF)
884e9012 24//////////////////////////////////////////////////////////////////////////////
25
26#include <TChain.h>
27#include <TFile.h>
28#include <TNtuple.h>
29#include <TVector3.h>
30
31#include "AliAnalysisManager.h"
32#include "AliMCEventHandler.h"
4b707925 33#include "AliMCEvent.h"
884e9012 34#include "AliAnaCaloTriggerMC.h"
35#include "AliESDEvent.h"
36#include "AliESDCaloCluster.h"
37#include "AliLog.h"
38#include "AliStack.h"
39#include "TParticle.h"
40
41//______________________________________________________________________________
42AliAnaCaloTriggerMC::AliAnaCaloTriggerMC() :
43 fChain(0),
44 fESD(0),
45 fOutputContainer(0),
46 fCalorimeter("PHOS"),
47 fNtTrigger22(0),
48 fNtTriggerNN(0)
49
50{
51 // Default constructor.
52
53}
54//______________________________________________________________________________
55AliAnaCaloTriggerMC::AliAnaCaloTriggerMC(const char *name) :
56 AliAnalysisTask(name, "AnaCaloTriggerMC"),
57 fChain(0),
58 fESD(0),
59 fOutputContainer(0),
60 fCalorimeter("PHOS"),
61 fNtTrigger22(0),
62 fNtTriggerNN(0)
63
64{
65 // Constructor.
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()) ;
70}
71//____________________________________________________________________________
72AliAnaCaloTriggerMC::AliAnaCaloTriggerMC(const AliAnaCaloTriggerMC & ct) :
3bb2c538 73 AliAnalysisTask(ct),fChain(ct.fChain), fESD(ct.fESD),
884e9012 74 fOutputContainer(ct.fOutputContainer), fCalorimeter(ct. fCalorimeter),
75 fNtTrigger22(ct.fNtTrigger22), fNtTriggerNN(ct.fNtTriggerNN)
76{
77
78 // cpy ctor
79 SetName (ct.GetName()) ;
80 SetTitle(ct.GetTitle()) ;
81
82}
83
84//_________________________________________________________________________
85AliAnaCaloTriggerMC & AliAnaCaloTriggerMC::operator = (const AliAnaCaloTriggerMC & source)
86{
87 // assignment operator
88
89 if(&source == this) return *this;
90
91 fChain = source.fChain ;
92 fESD = source.fESD ;
93 fOutputContainer = source.fOutputContainer ;
94 fCalorimeter = source. fCalorimeter ;
95 fNtTrigger22 = source.fNtTrigger22 ;
96 fNtTriggerNN = source.fNtTriggerNN ;
97
98 return *this;
99
100}
101
102//______________________________________________________________________________
103AliAnaCaloTriggerMC::~AliAnaCaloTriggerMC()
104{
105 // dtor
3788def4 106 // fOutputContainer->Clear() ;
107 // delete fOutputContainer ;
108
884e9012 109}
110
111
112//______________________________________________________________________________
113void AliAnaCaloTriggerMC::ConnectInputData(const Option_t*)
114{
115 // Initialisation of branch container and histograms
116
117 AliInfo(Form("*** Initialization of %s", GetName())) ;
118
119 // Get input data
120 fChain = dynamic_cast<TChain *>(GetInputData(0)) ;
121 if (!fChain) {
122 AliError(Form("Input 0 for %s not found\n", GetName()));
123 return ;
124 }
125
126 fESD = new AliESDEvent();
127 fESD->ReadFromTree(fChain);
128
129}
130
131//________________________________________________________________________
132
133void AliAnaCaloTriggerMC::CreateOutputObjects()
134{
135
136 // Create the output container
137 OpenFile(0);
138
139 // create histograms
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");
142
143 // create output container
144
145 fOutputContainer = new TObjArray(2) ;
146 fOutputContainer->SetName(GetName()) ;
147
148 fOutputContainer->AddAt(fNtTrigger22, 0) ;
149 fOutputContainer->AddAt(fNtTriggerNN, 1) ;
150
151}
152
153//______________________________________________________________________________
154void AliAnaCaloTriggerMC::Exec(Option_t *)
155{
156
157 // Processing of one event
158
159 Long64_t entry = fChain->GetReadEntry() ;
160
161 if (!fESD) {
162 AliError("fESD is not connected to the input!") ;
163 return ;
164 }
165
166 if ( !((entry-1)%100) )
167 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ;
168
169 //Get MC data
170 AliStack* stack = 0x0;
171 AliMCEventHandler* mctruth = (AliMCEventHandler*)
172 ((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
173
174 if(mctruth)
4b707925 175 stack = mctruth->MCEvent()->Stack();
884e9012 176
177 if (!stack) {
178 AliError("Stack not found") ;
179 return ;
180 }
181
182 // Get trigger information of fCalorimeter
183 TArrayF * triggerAmplitudes = 0x0 ;
184 TArrayF * triggerPosition = 0x0 ;
185 Int_t numberOfCaloClusters = 0 ;
186
187 if(fCalorimeter == "PHOS"){
188 triggerAmplitudes = fESD->GetPHOSTriggerAmplitudes();
189 triggerPosition = fESD->GetPHOSTriggerPosition();
190 }
191 else if(fCalorimeter == "EMCAL"){
192 triggerAmplitudes = fESD->GetEMCALTriggerAmplitudes();
193 triggerPosition = fESD->GetEMCALTriggerPosition();
194 }
195
196 if( triggerAmplitudes && triggerPosition ){
197 // trigger amplitudes
f38b754b 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)) ;
884e9012 202
203 // trigger position
f38b754b 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)) ;
884e9012 210
211 Float_t enMax = 0. ;
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. ;
217
f38b754b 218 TVector3 vpos22(kx22, ky22, kz22) ;
219 TVector3 vposNN(kxNN, kyNN, kzNN) ;
884e9012 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. ;
224
225 Int_t icaloCluster ;
226
227 // loop over the Calorimeters Clusters
228 Float_t cluEnergy = 0;
229 Int_t labelmax = -5;
230 for(icaloCluster = 0 ; icaloCluster < numberOfCaloClusters ; icaloCluster++) {
231
232 AliESDCaloCluster * cluster = fESD->GetCaloCluster(icaloCluster) ;
233
234 if (cluster && ( (fCalorimeter == "PHOS" && cluster->IsPHOS()) ||
235 (fCalorimeter == "EMCAL" && cluster->IsEMCAL()))) {
236
237 cluEnergy = cluster->E() ;
238 Float_t pos[3] ;
239 TVector3 vpos ;
240
241 cluster->GetPosition( pos ) ;
242
243 if ( cluEnergy > enMax) {
244 enMax = cluEnergy ;
245 vpos.SetXYZ(pos[0], pos[1], pos[2]) ;
246 etaMax = vpos.Eta() ;
247 phiMax = vpos.Phi() ;
248 labelmax = cluster->GetLabel();
249 }
250
4b707925 251 Double_t * pid = cluster->GetPid() ;
884e9012 252
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() ;
259 }
260 }
261 }//if cluster
262 }//CaloCluster loop
263
264 if(labelmax < stack->GetNtrack() && labelmax >= 0 ){
265 TParticle * particle = stack->Particle(labelmax);
266 Float_t ptgen = particle->Energy();
f38b754b 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.);
884e9012 269 }
270 else AliDebug(1, Form("Wrong label %d, ntrack %d, Emax %f ",labelmax, stack->GetNtrack(), phEnMax));
271 }//If trigger arrays filled
272
273 PostData(0, fOutputContainer);
274
275}
276
277
278//______________________________________________________________________________
279void AliAnaCaloTriggerMC::Terminate(Option_t *)
280{
281 // Processing when the event loop is ended
282
283}