]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/totEt/AliAnalysisEtMonteCarlo.cxx
correct inputHistogramFileName (Prabhat Ranjan Pujahari <p.pujahari@cern.ch>)
[u/mrichter/AliRoot.git] / PWGLF / totEt / AliAnalysisEtMonteCarlo.cxx
CommitLineData
cf6522d1 1//_________________________________________________________________________
2// Utility Class for transverse energy studies
3// Base class for MC analysis
4// - MC output
5// implementation file
6//
7//*-- Authors: Oystein Djuvsland (Bergen), David Silvermyr (ORNL)
e2ee5727 8//_________________________________________________________________________
cf6522d1 9
2fbf38ac 10#include "AliAnalysisEtMonteCarlo.h"
11#include "AliAnalysisEtCuts.h"
f61cec2f 12#include "AliAnalysisEtSelector.h"
13#include "AliAnalysisEtSelector.h"
13b0d3c1 14#include "AliESDtrack.h"
2fbf38ac 15#include "AliStack.h"
0651f6b4 16#include "AliVEvent.h"
2fbf38ac 17#include "AliMCEvent.h"
0651f6b4 18#include "AliESDEvent.h"
13b0d3c1 19#include "TH2F.h"
e2ee5727 20#include "TH3F.h"
cf6522d1 21#include "TParticle.h"
ce546038 22#include "AliGenHijingEventHeader.h"
23#include "AliGenPythiaEventHeader.h"
0651f6b4 24#include "TList.h"
25#include "AliESDCaloCluster.h"
0f6416f3 26#include "AliLog.h"
e2ee5727 27#include <iostream>
28#include <AliCentrality.h>
ef647350 29#include "AliPHOSGeoUtils.h"
30#include "AliPHOSGeometry.h"
31#include "TFile.h"
c79c36cc 32#include "AliESDtrackCuts.h"
e2ee5727 33using namespace std;
16abb579 34
35ClassImp(AliAnalysisEtMonteCarlo);
36
37
ce546038 38// ctor
0651f6b4 39AliAnalysisEtMonteCarlo::AliAnalysisEtMonteCarlo():AliAnalysisEt()
c79c36cc 40 ,fIsMC(kTRUE)
f61cec2f 41 ,fImpactParameter(0)
42 ,fNcoll(0)
43 ,fNpart(0)
44 ,fPrimaryTree(0)
45 ,fTotEtWithSecondaryRemoved(0)
46 ,fTotEtSecondaryFromEmEtPrimary(0)
47 ,fTotEtSecondary(0)
48 ,fPrimaryCode(0)
49 ,fPrimaryCharge(0)
50 ,fPrimaryE(0)
51 ,fPrimaryEt(0)
52 ,fPrimaryPx(0)
53 ,fPrimaryPy(0)
54 ,fPrimaryPz(0)
55 ,fPrimaryVx(0)
56 ,fPrimaryVy(0)
57 ,fPrimaryVz(0)
58 ,fPrimaryAccepted(0)
b2c10007 59 ,fPrimaryMatched(0)
f61cec2f 60 ,fDepositedCode(0)
b2c10007 61 ,fDepositedE(0)
f61cec2f 62 ,fDepositedEt(0)
63 ,fDepositedCharge(0)
64 ,fDepositedVx(0)
65 ,fDepositedVy(0)
66 ,fDepositedVz(0)
b2c10007 67 ,fSecondary(kFALSE)
68 ,fReconstructedE(0)
69 ,fReconstructedEt(0)
70 ,fTotPx(0)
71 ,fTotPy(0)
72 ,fTotPz(0)
73 ,fClusterMult(0)
f61cec2f 74 ,fHistDecayVertexNonRemovedCharged(0)
75 ,fHistDecayVertexRemovedCharged(0)
76 ,fHistDecayVertexNonRemovedNeutral(0)
77 ,fHistDecayVertexRemovedNeutral(0)
78
79 ,fHistRemovedOrNot(0)
80 ,fHistEtNonRemovedProtons(0)
81 ,fHistEtNonRemovedAntiProtons(0)
82 ,fHistEtNonRemovedPiPlus(0)
83 ,fHistEtNonRemovedPiMinus(0)
84 ,fHistEtNonRemovedKaonPlus(0)
85 ,fHistEtNonRemovedKaonMinus(0)
86 ,fHistEtNonRemovedK0s(0)
87 ,fHistEtNonRemovedK0L(0)
88 ,fHistEtNonRemovedLambdas(0)
89 ,fHistEtNonRemovedElectrons(0)
90 ,fHistEtNonRemovedPositrons(0)
91 ,fHistEtNonRemovedMuPlus(0)
92 ,fHistEtNonRemovedMuMinus(0)
93 ,fHistEtNonRemovedNeutrons(0)
94 ,fHistEtNonRemovedAntiNeutrons(0)
95 ,fHistEtNonRemovedGammas(0)
96 ,fHistEtNonRemovedGammasFromPi0(0)
97 ,fHistEtRemovedGammas(0)
98 ,fHistEtRemovedNeutrons(0)
99 ,fHistEtRemovedAntiNeutrons(0)
100 ,fHistEtRemovedCharged(0)
101 ,fHistEtRemovedNeutrals(0)
102 ,fHistEtNonRemovedCharged(0)
103 ,fHistEtNonRemovedNeutrals(0)
104 ,fHistMultNonRemovedProtons(0)
105 ,fHistMultNonRemovedAntiProtons(0)
106 ,fHistMultNonRemovedPiPlus(0)
107 ,fHistMultNonRemovedPiMinus(0)
108 ,fHistMultNonRemovedKaonPlus(0)
109 ,fHistMultNonRemovedKaonMinus(0)
110 ,fHistMultNonRemovedK0s(0)
111 ,fHistMultNonRemovedK0L(0)
112 ,fHistMultNonRemovedLambdas(0)
113 ,fHistMultNonRemovedElectrons(0)
114 ,fHistMultNonRemovedPositrons(0)
115 ,fHistMultNonRemovedMuPlus(0)
116 ,fHistMultNonRemovedMuMinus(0)
117 ,fHistMultNonRemovedNeutrons(0)
118 ,fHistMultNonRemovedAntiNeutrons(0)
119 ,fHistMultNonRemovedGammas(0)
120 ,fHistMultRemovedGammas(0)
121 ,fHistMultRemovedNeutrons(0)
122 ,fHistMultRemovedAntiNeutrons(0)
123 ,fHistMultRemovedCharged(0)
124 ,fHistMultRemovedNeutrals(0)
125 ,fHistMultNonRemovedCharged(0)
126 ,fHistMultNonRemovedNeutrals(0)
127 ,fHistTrackMultvsNonRemovedCharged(0)
128 ,fHistTrackMultvsNonRemovedNeutral(0)
129 ,fHistTrackMultvsRemovedGamma(0)
130 ,fHistClusterMultvsNonRemovedCharged(0)
131 ,fHistClusterMultvsNonRemovedNeutral(0)
132 ,fHistClusterMultvsRemovedGamma(0)
133 ,fHistMultvsNonRemovedChargedE(0)
134 ,fHistMultvsNonRemovedNeutralE(0)
03bfd268 135 ,fHistMultvsRemovedGammaE(0)
0861cc1f 136 ,fCalcForKaonCorrection(kFALSE)
137 ,fHistK0EDepositsVsPtInAcceptance(0)
138 ,fHistK0EGammaVsPtInAcceptance(0)
139 ,fHistK0EDepositsVsPtOutOfAcceptance(0)
140 ,fHistK0EGammaVsPtOutOfAcceptance(0)
141 ,fHistSimKaonsInAcceptance(0)
142 ,fHistSimK0SInAcceptance(0)
143 ,fHistSimKPlusInAcceptance(0)
144 ,fHistSimKMinusInAcceptance(0)
145 ,fHistSimK0LInAcceptance(0)
146 ,fHistSimKaonsOutOfAcceptance(0)
147 ,fHistSimKaonsInAcceptanceWithDepositsPrimaries(0)
148 ,fHistSimKaonsOutOfAcceptanceWithDepositsSecondaries(0)
149 ,fHistSimKaonsOutOfAcceptanceWithDepositsPrimaries(0)
f61cec2f 150 ,fEtNonRemovedProtons(0)
151 ,fEtNonRemovedAntiProtons(0)
152 ,fEtNonRemovedPiPlus(0)
153 ,fEtNonRemovedPiMinus(0)
154 ,fEtNonRemovedKaonPlus(0)
155 ,fEtNonRemovedKaonMinus(0)
156 ,fEtNonRemovedK0S(0)
157 ,fEtNonRemovedK0L(0)
158 ,fEtNonRemovedLambdas(0)
159 ,fEtNonRemovedElectrons(0)
160 ,fEtNonRemovedPositrons(0)
161 ,fEtNonRemovedMuMinus(0)
162 ,fEtNonRemovedMuPlus(0)
163 ,fEtNonRemovedGammas(0)
164 ,fEtNonRemovedGammasFromPi0(0)
165 ,fEtNonRemovedNeutrons(0)
166 ,fEtNonRemovedAntiNeutrons(0)
167 ,fEtRemovedProtons(0)
168 ,fEtRemovedAntiProtons(0)
169 ,fEtRemovedPiPlus(0)
170 ,fEtRemovedPiMinus(0)
171 ,fEtRemovedKaonPlus(0)
172 ,fEtRemovedKaonMinus(0)
173 ,fEtRemovedK0s(0)
174 ,fEtRemovedK0L(0)
175 ,fEtRemovedLambdas(0)
176 ,fEtRemovedElectrons(0)
177 ,fEtRemovedPositrons(0)
178 ,fEtRemovedMuMinus(0)
179 ,fEtRemovedMuPlus(0)
180 ,fEtRemovedGammasFromPi0(0)
181 ,fEtRemovedGammas(0)
182 ,fEtRemovedNeutrons(0)
183 ,fEtRemovedAntiNeutrons(0)
184 ,fMultNonRemovedProtons(0)
185 ,fMultNonRemovedAntiProtons(0)
186 ,fMultNonRemovedPiPlus(0)
187 ,fMultNonRemovedPiMinus(0)
188 ,fMultNonRemovedKaonPlus(0)
189 ,fMultNonRemovedKaonMinus(0)
190 ,fMultNonRemovedK0s(0)
191 ,fMultNonRemovedK0L(0)
192 ,fMultNonRemovedLambdas(0)
193 ,fMultNonRemovedElectrons(0)
194 ,fMultNonRemovedPositrons(0)
195 ,fMultNonRemovedMuMinus(0)
196 ,fMultNonRemovedMuPlus(0)
197 ,fMultNonRemovedGammas(0)
198 ,fMultNonRemovedNeutrons(0)
199 ,fMultNonRemovedAntiNeutrons(0)
200 ,fMultRemovedProtons(0)
201 ,fMultRemovedAntiProtons(0)
202 ,fMultRemovedPiPlus(0)
203 ,fMultRemovedPiMinus(0)
204 ,fMultRemovedKaonPlus(0)
205 ,fMultRemovedKaonMinus(0)
206 ,fMultRemovedK0s(0)
207 ,fMultRemovedK0L(0)
208 ,fMultRemovedLambdas(0)
209 ,fMultRemovedElectrons(0)
210 ,fMultRemovedPositrons(0)
211 ,fMultRemovedMuMinus(0)
212 ,fMultRemovedMuPlus(0)
213 ,fMultRemovedGammas(0)
214 ,fMultRemovedNeutrons(0)
215 ,fMultRemovedAntiNeutrons(0)
216 ,fTrackMultInAcc(0)
217 ,fHistDxDzNonRemovedCharged(0)
218 ,fHistDxDzRemovedCharged(0)
219 ,fHistDxDzNonRemovedNeutral(0)
220 ,fHistDxDzRemovedNeutral(0)
221 ,fHistPiPlusMult(0)
222 ,fHistPiMinusMult(0)
223 ,fHistPiZeroMult(0)
224 ,fHistPiPlusMultAcc(0)
225 ,fHistPiMinusMultAcc(0)
226 ,fHistPiZeroMultAcc(0)
b2c10007 227// ,fPiPlusMult(0)
228// ,fPiMinusMult(0)
f61cec2f 229 ,fPiZeroMult(0)
230 ,fPiPlusMultAcc(0)
231 ,fPiMinusMultAcc(0)
232 ,fPiZeroMultAcc(0)
233 ,fNeutralRemoved(0)
234 ,fChargedRemoved(0)
235 ,fChargedNotRemoved(0)
236 ,fNeutralNotRemoved(0)
b2c10007 237 ,fGammaRemoved(0)
238 ,fSecondaryNotRemoved(0)
f61cec2f 239 ,fEnergyNeutralRemoved(0)
240 ,fEnergyChargedRemoved(0)
241 ,fEnergyChargedNotRemoved(0)
242 ,fEnergyNeutralNotRemoved(0)
b2c10007 243 ,fEnergyGammaRemoved(0)
f61cec2f 244 ,fNClusters(0)
245 ,fTotNeutralEtAfterMinEnergyCut(0)
0861cc1f 246 ,fCalcTrackMatchVsMult(kFALSE)
247 ,fHistGammasFound(0)
9ef6f13f 248 ,fHistGammasGenerated(0)
5881f036 249 ,fHistGammasFoundCent(0)
250 ,fHistGammasGeneratedCent(0)
9ef6f13f 251 ,fHistChargedTracksCut(0)
252 ,fHistChargedTracksAccepted(0)
253 ,fHistGammasCut(0)
254 ,fHistGammasAccepted(0)
0861cc1f 255 ,fHistChargedTracksCutMult(0)
256 ,fHistChargedTracksAcceptedMult(0)
ac610b08 257 ,fHistChargedTracksAcceptedLowPtCent(0)
5881f036 258 ,fHistChargedTracksAcceptedLowPtCent500MeV(0)
ac610b08 259 ,fHistChargedTracksAcceptedLowPtCentNoAntiProtons(0)
260 ,fHistChargedTracksAcceptedLowPtCentAntiProtons(0)
0861cc1f 261 ,fHistGammasCutMult(0)
262 ,fHistGammasAcceptedMult(0)
9ef6f13f 263 ,fHistBadTrackMatches(0)
264 ,fHistMatchedTracksEvspTBkgd(0)
265 ,fHistMatchedTracksEvspTSignal(0)
266 ,fHistMatchedTracksEvspTBkgdPeripheral(0)
267 ,fHistMatchedTracksEvspTSignalPeripheral(0)
ac610b08 268 ,fHistMatchedTracksEvspTBkgdvsCent(0)
269 ,fHistMatchedTracksEvspTSignalvsCent(0)
270 ,fHistMatchedTracksEvspTBkgdvsCentEffCorr(0)
271 ,fHistMatchedTracksEvspTSignalvsCentEffCorr(0)
0861cc1f 272
9ef6f13f 273 ,fHistChargedTracksCutPeripheral(0)
274 ,fHistChargedTracksAcceptedPeripheral(0)
275 ,fHistGammasCutPeripheral(0)
276 ,fHistGammasAcceptedPeripheral(0)
277 ,fHistBadTrackMatchesdPhidEta(0)
278 ,fHistGoodTrackMatchesdPhidEta(0)
c79c36cc 279 ,fHistHadronDepositsAll(0)
280 ,fHistHadronDepositsReco(0)
ac610b08 281 ,fHistHadronDepositsAllCent(0)
5881f036 282 ,fHistHadronDepositsAllCent500MeV(0)
ac610b08 283 ,fHistHadronDepositsRecoCent(0)
284 ,fHistHadronsAllCent(0)
d3ce32b8 285 ,fHistMultChVsSignalVsMult(0)
ce546038 286{
287}
288
289// dtor
e2ee5727 290AliAnalysisEtMonteCarlo::~AliAnalysisEtMonteCarlo()
f61cec2f 291{ //Destructor
311c6540 292
293 if(fPrimaryTree){
294 fPrimaryTree->Clear();
295 delete fPrimaryTree;
296 }
f61cec2f 297 delete fHistDecayVertexNonRemovedCharged; // Decay vertex for non-removed charged particles
298 delete fHistDecayVertexRemovedCharged; // Decay vertex for non-removed charged particles
299 delete fHistDecayVertexNonRemovedNeutral; // Decay vertex for non-removed charged particles
300 delete fHistDecayVertexRemovedNeutral; // Decay vertex for non-removed charged particles
301
302 delete fHistRemovedOrNot; // If charged/neutral particles were removed or not
303
304 delete fHistEtNonRemovedProtons; // enter comment here
305 delete fHistEtNonRemovedAntiProtons; // enter comment here
306 delete fHistEtNonRemovedPiPlus; // enter comment here
307 delete fHistEtNonRemovedPiMinus; // enter comment here
308 delete fHistEtNonRemovedKaonPlus; // enter comment here
309 delete fHistEtNonRemovedKaonMinus; // enter comment here
310 delete fHistEtNonRemovedK0s; // enter comment here
311 delete fHistEtNonRemovedK0L; // enter comment here
312 delete fHistEtNonRemovedLambdas; // enter comment here
313 delete fHistEtNonRemovedElectrons; // enter comment here
314 delete fHistEtNonRemovedPositrons; // enter comment here
315 delete fHistEtNonRemovedMuPlus; // enter comment here
316 delete fHistEtNonRemovedMuMinus; // enter comment here
317 delete fHistEtNonRemovedNeutrons; // enter comment here
318 delete fHistEtNonRemovedAntiNeutrons; // enter comment here
319 delete fHistEtNonRemovedGammas; // enter comment here
320 delete fHistEtNonRemovedGammasFromPi0; // enter comment here
321
322 delete fHistEtRemovedGammas; // enter comment here
323 delete fHistEtRemovedNeutrons; // enter comment here
324 delete fHistEtRemovedAntiNeutrons; // enter comment here
325
326
327 delete fHistMultNonRemovedProtons; // enter comment here
328 delete fHistMultNonRemovedAntiProtons; // enter comment here
329 delete fHistMultNonRemovedPiPlus; // enter comment here
330 delete fHistMultNonRemovedPiMinus; // enter comment here
331 delete fHistMultNonRemovedKaonPlus; // enter comment here
332 delete fHistMultNonRemovedKaonMinus; // enter comment here
333 delete fHistMultNonRemovedK0s; // enter comment here
334 delete fHistMultNonRemovedK0L; // enter comment here
335 delete fHistMultNonRemovedLambdas; // enter comment here
336 delete fHistMultNonRemovedElectrons; // enter comment here
337 delete fHistMultNonRemovedPositrons; // enter comment here
338 delete fHistMultNonRemovedMuPlus; // enter comment here
339 delete fHistMultNonRemovedMuMinus; // enter comment here
340 delete fHistMultNonRemovedNeutrons; // enter comment here
341 delete fHistMultNonRemovedAntiNeutrons; // enter comment here
342 delete fHistMultNonRemovedGammas; // enter comment here
343
344 delete fHistMultRemovedGammas; // enter comment here
345 delete fHistMultRemovedNeutrons; // enter comment here
346 delete fHistMultRemovedAntiNeutrons; // enter comment here
347
348 delete fHistTrackMultvsNonRemovedCharged; // enter comment here
349 delete fHistTrackMultvsNonRemovedNeutral; // enter comment here
350 delete fHistTrackMultvsRemovedGamma; // enter comment here
351
352 delete fHistClusterMultvsNonRemovedCharged; // enter comment here
353 delete fHistClusterMultvsNonRemovedNeutral; // enter comment here
354 delete fHistClusterMultvsRemovedGamma; // enter comment here
355
356 delete fHistMultvsNonRemovedChargedE; // enter comment here
357 delete fHistMultvsNonRemovedNeutralE; // enter comment here
358 delete fHistMultvsRemovedGammaE; // enter comment here
0861cc1f 359 delete fHistK0EDepositsVsPtInAcceptance;//enter comment here
360 delete fHistK0EGammaVsPtInAcceptance;//enter comment here
361 delete fHistK0EDepositsVsPtOutOfAcceptance;
362 delete fHistK0EGammaVsPtOutOfAcceptance;
363
364 delete fHistSimKaonsInAcceptance;
365 delete fHistSimK0SInAcceptance;
366 delete fHistSimKPlusInAcceptance;
367 delete fHistSimKMinusInAcceptance;
368 delete fHistSimK0LInAcceptance;
369 delete fHistSimKaonsOutOfAcceptance;
370 delete fHistSimKaonsInAcceptanceWithDepositsPrimaries;
371 delete fHistSimKaonsOutOfAcceptanceWithDepositsSecondaries;
372 delete fHistSimKaonsOutOfAcceptanceWithDepositsPrimaries;
2aab9269 373
f61cec2f 374 delete fHistDxDzNonRemovedCharged; // enter comment here
375 delete fHistDxDzRemovedCharged; // enter comment here
376 delete fHistDxDzNonRemovedNeutral; // enter comment here
377 delete fHistDxDzRemovedNeutral; // enter comment here
378
379 delete fHistPiPlusMult; // enter comment here
380 delete fHistPiMinusMult; // enter comment here
381 delete fHistPiZeroMult; // enter comment here
382
383 delete fHistPiPlusMultAcc; // enter comment here
384 delete fHistPiMinusMultAcc; // enter comment here
385 delete fHistPiZeroMultAcc; // enter comment here
4503e29d 386 delete fHistGammasFound; // enter comment here
387 delete fHistGammasGenerated; // enter comment here
5881f036 388 delete fHistGammasFoundCent; // enter comment here
389 delete fHistGammasGeneratedCent; // enter comment here
9ef6f13f 390 delete fHistChargedTracksCut;
391 delete fHistChargedTracksAccepted;
392 delete fHistGammasCut;
393 delete fHistGammasAccepted;
0861cc1f 394 delete fHistChargedTracksCutMult;
395 delete fHistChargedTracksAcceptedMult;
ac610b08 396 delete fHistChargedTracksAcceptedLowPtCent;
5881f036 397 delete fHistChargedTracksAcceptedLowPtCent500MeV;
ac610b08 398 delete fHistChargedTracksAcceptedLowPtCentNoAntiProtons;
399 delete fHistChargedTracksAcceptedLowPtCentAntiProtons;
0861cc1f 400 delete fHistGammasCutMult;
401 delete fHistGammasAcceptedMult;
9ef6f13f 402 delete fHistBadTrackMatches;
403 delete fHistMatchedTracksEvspTBkgd;
404 delete fHistMatchedTracksEvspTSignal;
405 delete fHistMatchedTracksEvspTBkgdPeripheral;
406 delete fHistMatchedTracksEvspTSignalPeripheral;
ac610b08 407 delete fHistMatchedTracksEvspTBkgdvsCent;
408 delete fHistMatchedTracksEvspTSignalvsCent;
409 delete fHistMatchedTracksEvspTBkgdvsCentEffCorr;
410 delete fHistMatchedTracksEvspTSignalvsCentEffCorr;
9ef6f13f 411 delete fHistChargedTracksCutPeripheral;
412 delete fHistChargedTracksAcceptedPeripheral;
413 delete fHistGammasCutPeripheral;
414 delete fHistGammasAcceptedPeripheral;
415 delete fHistBadTrackMatchesdPhidEta;
416 delete fHistGoodTrackMatchesdPhidEta;
c79c36cc 417 delete fHistHadronDepositsAll;
418 delete fHistHadronDepositsReco;
ac610b08 419 delete fHistHadronDepositsAllCent;
5881f036 420 delete fHistHadronDepositsAllCent500MeV;
ac610b08 421 delete fHistHadronDepositsRecoCent;
422 delete fHistHadronsAllCent;
d3ce32b8 423 delete fHistMultChVsSignalVsMult;
ce546038 424}
2fbf38ac 425
426Int_t AliAnalysisEtMonteCarlo::AnalyseEvent(AliVEvent* ev)
f61cec2f 427{ // analyse MC event
c79c36cc 428 if(!fIsMC) return 0;
e2ee5727 429 ResetEventValues();
430
b2c10007 431
e2ee5727 432 fPiPlusMult = 0;
433 fPiMinusMult = 0;
434 fPiZeroMult = 0;
435 if (fCentrality)
436 {
ac610b08 437 fCentClass = fCentrality->GetCentralityClass5(fCentralityMethod);
e2ee5727 438
439 }
e2ee5727 440
441 // Get us an mc event
442 if (!ev) {
443 AliFatal("ERROR: Event does not exist");
444 return 0;
445 }
446 AliMCEvent *event = dynamic_cast<AliMCEvent*>(ev);
447 if (!event) {
448 AliFatal("ERROR: MC Event does not exist");
449 return 0;
450 }
451
452 Double_t protonMass =fgProtonMass;
453
454 // Hijing header
455 AliGenEventHeader* genHeader = event->GenEventHeader();
456 if (!genHeader) {
457 Printf("ERROR: Event generation header does not exist");
458 return 0;
459 }
460 AliGenHijingEventHeader* hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
461 if (hijingGenHeader) {
462 fImpactParameter = hijingGenHeader->ImpactParameter();
463 fNcoll = hijingGenHeader->HardScatters(); // or should this be some combination of NN() NNw() NwN() NwNw() ?
464 fNpart = hijingGenHeader->ProjectileParticipants() + hijingGenHeader->TargetParticipants();
f61cec2f 465 }
466
e2ee5727 467 // Let's play with the stack!
468 AliStack *stack = event->Stack();
469
470 Int_t nPrim = stack->GetNtrack();
471
472 Int_t partCount = 0;
473 for (Int_t iPart = 0; iPart < nPrim; iPart++)
2fbf38ac 474 {
e2ee5727 475
476 TParticle *part = stack->Particle(iPart);
b2c10007 477
478
e2ee5727 479
480 if (!part)
481 {
482 Printf("ERROR: Could not get particle %d", iPart);
483 continue;
484 }
e2ee5727 485 TParticlePDG *pdg = part->GetPDG(0);
e2ee5727 486
f61cec2f 487 if (!pdg)
2fbf38ac 488 {
f61cec2f 489 Printf("ERROR: Could not get particle PDG %d", iPart);
e2ee5727 490 continue;
2fbf38ac 491 }
492
e2ee5727 493 Double_t particleMassPart = 0; //The mass part in the Et calculation for this particle
f61cec2f 494 Int_t code = pdg->PdgCode();
b2c10007 495
496 if(stack->IsPhysicalPrimary(iPart))
497 {
498 fTotPx += part->Px();
499 fTotPy += part->Py();
500 fTotPz += part->Pz();
501 }
502
e2ee5727 503 // Check for reasonable (for now neutral and singly charged) charge on the particle
86e7d5db 504 if (fSelector->IsNeutralMcParticle(iPart,*stack,*pdg))
e2ee5727 505 {
506
507 fMultiplicity++;
b2c10007 508// PrintFamilyTree(iPart, stack);
f61cec2f 509//
e2ee5727 510 if (TMath::Abs(part->Eta()) < fCuts->GetCommonEtaCut())
2fbf38ac 511 {
f61cec2f 512 //Printf("Particle with eta: %f, pid: %d", part->Eta(), code);
e2ee5727 513 // calculate E_T
514 if (
ef647350 515 TMath::Abs(code) == fgProtonCode ||
516 TMath::Abs(code) == fgNeutronCode ||
517 TMath::Abs(code) == fgLambdaCode ||
518 TMath::Abs(code) == fgXiCode ||
519 TMath::Abs(code) == fgXi0Code ||
520 TMath::Abs(code) == fgOmegaCode
e2ee5727 521 )
522 {
ef647350 523 if (code > 0) {
e2ee5727 524 particleMassPart = - protonMass;
525 }
ef647350 526 if (code < 0) {
e2ee5727 527 particleMassPart = protonMass;
528 }
529 }
530 Double_t et = part->Energy() * TMath::Sin(part->Theta()) + particleMassPart;
531
e2ee5727 532
6a152780 533// // Fill up total E_T counters for each particle species
534// if (code == fgProtonCode || code == fgAntiProtonCode)
535// {
536// }
537// if (code == fgPiPlusCode || code == fgPiMinusCode)
538// {
539// if (code == fgPiPlusCode)
540// {
541// }
542// else
543// {
544// }
545// }
546// if (code == fgGammaCode)
547// {
548// }
549// if (code == fgKPlusCode)
550// {
551// }
552// if(code == fgKMinusCode)
553// {
554// }
555// if (code == fgMuPlusCode || code == fgMuMinusCode)
556// {
557// }
558// if (code == fgEPlusCode || code == fgEMinusCode)
559// {
560// }
561// // some neutrals also
562// if (code == fgNeutronCode)
563// {
564// }
565// if (code == fgAntiNeutronCode)
566// {
567// }
568// if (code == fgGammaCode)
569// {
570// }
e2ee5727 571
572 // Neutral particles
ef647350 573 //if (TMath::Abs(pdg->Charge() - fCuts->GetMonteCarloNeutralParticle()) <1e-3 )
f61cec2f 574
ef647350 575 if(code == fgGammaCode || code == fgPi0Code || code == fgEtaCode)
e2ee5727 576 {
b2c10007 577 // PrintFamilyTree(iPart,stack);
f61cec2f 578 //Printf("Gamma, phi: %f, eta: %f, phi cut min: %f, phi cut max: %f, eta cut: %f", part->Phi(), part->Eta(), fPhiMinCutAcc, fPhiMaxCutAcc, fEtaCutAcc);
579 //if (et > fCuts->GetCommonClusterEnergyCut()) fTotNeutralEt += et;
e2ee5727 580
581 // inside EMCal acceptance
f61cec2f 582
583 //if (TMath::Abs(part->Eta()) < fEtaCutAcc && part->Phi() < fPhiMaxCutAcc && part->Phi() > fPhiMinCutAcc)
584
b2c10007 585 if(fSelector->CutGeometricalAcceptance(*part) )
e2ee5727 586 {
587 fNeutralMultiplicity++;
f61cec2f 588 fTotNeutralEt += et;
86e7d5db 589 if(fSelector->PassMinEnergyCut(*part))
f61cec2f 590 {
591 fTotNeutralEtAfterMinEnergyCut += et;
592 }
e2ee5727 593 if (part->Energy() > 0.05) partCount++;
e2ee5727 594 }
595 }
596 //Charged particles
597 else if (TMath::Abs( pdg->Charge() - fCuts->GetMonteCarloNeutralParticle())>1e-3 )
598 {
e2ee5727 599
600 // inside EMCal acceptance
f61cec2f 601 if (fSelector->CutGeometricalAcceptance(*part))
e2ee5727 602 {
f61cec2f 603
604 fChargedMultiplicity++;
605
606 fTotChargedEt += et;
e2ee5727 607
ef647350 608 if (code == fgProtonCode || code == fgAntiProtonCode)
e2ee5727 609 {
e2ee5727 610 }
ef647350 611 if (code == fgPiPlusCode || code == fgPiMinusCode)
e2ee5727 612 {
e2ee5727 613 }
ef647350 614 if (code == fgKPlusCode || code == fgKMinusCode)
e2ee5727 615 {
e2ee5727 616 }
ef647350 617 if (code == fgMuPlusCode || code == fgMuMinusCode)
e2ee5727 618 {
e2ee5727 619 }
f61cec2f 620
ef647350 621 if (code == fgEPlusCode || code == fgEMinusCode)
e2ee5727 622 {
f61cec2f 623 fTotNeutralEt += et; // calling electrons neutral
624 fTotChargedEt -= et;
e1fa1966 625 }
e2ee5727 626 } // inside EMCal acceptance
627
e2ee5727 628 if (TrackHitsCalorimeter(part)) // magnetic field info not filled?
629 {
630 if (pdg->Charge() > 0) fHistPhivsPtPos->Fill(part->Phi(),part->Pt());
631 else fHistPhivsPtNeg->Fill(part->Phi(), part->Pt());
632 }
633 }
634 }
635
636 }
2fbf38ac 637 }
b2c10007 638 // std::cout << "Total: p_x = " << fTotPx << ", p_y = " << fTotPy << ", p_z = " << fTotPz << std::endl;
e2ee5727 639 fTotEt = fTotChargedEt + fTotNeutralEt;
f61cec2f 640 //fTotEtAcc = fTotChargedEtAcc + fTotNeutralEtAcc;//
e2ee5727 641 //std::cout << "Event done! # of particles: " << partCount << std::endl;
642 return 0;
2fbf38ac 643}
0651f6b4 644//Int_t AliAnalysisEtMonteCarlo::AnalyseEvent(AliMCEvent* mcEvent,AliESDEvent* realEvent)
645Int_t AliAnalysisEtMonteCarlo::AnalyseEvent(AliVEvent* ev,AliVEvent* ev2)
f61cec2f 646{ // analyse MC and real event info
e2ee5727 647 //if(!mcEvent || !realEvent){
c79c36cc 648 if(!fIsMC) return 0;
e2ee5727 649 if (!ev || !ev2) {
650 AliFatal("ERROR: Event does not exist");
651 return 0;
652 }
f61cec2f 653 AliAnalysisEt::AnalyseEvent(ev);
e2ee5727 654 AliMCEvent *mcEvent = dynamic_cast<AliMCEvent*>(ev);
655 AliESDEvent *realEvent = dynamic_cast<AliESDEvent*>(ev2);
656 if (!mcEvent || !realEvent) {
657 AliFatal("ERROR: mcEvent or realEvent does not exist");
b2c10007 658
e2ee5727 659 }
660
b2c10007 661 std::vector<Int_t> foundGammas;
662
f61cec2f 663 fSelector->SetEvent(realEvent);
664
665 AnalyseEvent(ev);
666
e2ee5727 667 AliStack *stack = mcEvent->Stack();
668
9ef6f13f 669 // get all detector clusters
f61cec2f 670 // TRefArray* caloClusters = new TRefArray();
671
672// if (fDetector == fCuts->GetDetectorEmcal()) realEvent->GetEMCALClusters( caloClusters );
673 //else if (fDetector == fCuts->GetDetectorPhos()) realEvent->GetPHOSClusters( caloClusters );
674 //else {
675 //AliFatal("Detector ID has not been specified");
676 //return -1;
677// }
678
4503e29d 679//Note that this only returns clusters for the selected detector. fSelector actually calls the right GetClusters... for the detector
cfe23ff0 680//It does not apply any cuts on these clusters
f61cec2f 681 TRefArray *caloClusters = fSelector->GetClusters();
e2ee5727 682
e2ee5727 683 Int_t nCluster = caloClusters->GetEntries();
b2c10007 684 fClusterMult = nCluster;
ef647350 685 fNClusters = 0;
d3ce32b8 686 Int_t fClusterMultChargedTracks = 0;
687 Int_t fClusterMultGammas = 0;
e2ee5727 688 // loop the clusters
689 for (int iCluster = 0; iCluster < nCluster; iCluster++ )
690 {
f61cec2f 691 Int_t cf = 0;
e2ee5727 692 AliESDCaloCluster* caloCluster = ( AliESDCaloCluster* )caloClusters->At( iCluster );
2aab9269 693 //Float_t caloE = caloCluster->E()
5881f036 694 if (!fSelector->CutGeometricalAcceptance(*caloCluster)) continue;
f61cec2f 695 fNClusters++;
b2c10007 696 const UInt_t iPart = (UInt_t)TMath::Abs(caloCluster->GetLabel());
ef647350 697 TParticle *part = stack->Particle(iPart);
e2ee5727 698
699 if (!part)
0651f6b4 700 {
e2ee5727 701 Printf("No MC particle %d", iCluster);
702 continue;
703 }
9ef6f13f 704
f61cec2f 705 int primIdx = iPart;
e2ee5727 706 if (!stack->IsPhysicalPrimary(iPart)) // check whether particle is primary. we keep secondary electron and gamma for testing.
707 {
f61cec2f 708 primIdx = GetPrimMother(iPart, stack);
b2c10007 709 if(primIdx != stack->GetPrimary(iPart))
710 {
711 //std::cout << primIdx << " = " << stack->GetPrimary(iPart) << std::endl;
712 //PrintFamilyTree(iPart, stack);
713 }
2aab9269 714 //if it is from a K0S
f61cec2f 715 if(primIdx < 0)
716 {
ac610b08 717 //std::cout << "What!? No primary?" << std::endl;
718 //PrintFamilyTree(iPart, stack);
9ef6f13f 719 //continue;
720 //This is a work around to fix a bug. For the EMCal when you use the tender supply, the parent particle ID gets messed up.
721 primIdx = iPart;
f61cec2f 722 }
723
e2ee5727 724 } // end of primary particle check
a0dcb3c1 725 //const int primCode = stack->Particle(primIdx)->GetPdgCode();
ef647350 726 TParticlePDG *pdg = part->GetPDG();
f61cec2f 727 if (!pdg)
ef647350 728 {
729 Printf("ERROR: Could not get particle PDG %d", iPart);
730 continue;
f61cec2f 731 }
e2ee5727 732
a0dcb3c1 733 //Int_t code = pdg->PdgCode();
53ad5ec1 734// if(primCode == fgGammaCode)
735// {
736
d3ce32b8 737
738
739 Bool_t nottrackmatched = kTRUE;//default to no track matched
6a152780 740 Float_t matchedTrackp = 0.0;
741 Float_t matchedTrackpt = 0.0;
d3ce32b8 742 nottrackmatched = fSelector->PassTrackMatchingCut(*caloCluster);
743 //by default ALL matched tracks are accepted, whether or not the match is good. So we check to see if the track is good.
584f2478 744 if(!nottrackmatched){//if the track is trackmatched
d3ce32b8 745 Int_t trackMatchedIndex = caloCluster->GetTrackMatchedIndex();
746 if(trackMatchedIndex < 0) nottrackmatched=kTRUE;
747 AliESDtrack *track = realEvent->GetTrack(trackMatchedIndex);
6a152780 748 matchedTrackp = track->P();
749 matchedTrackpt = track->Pt();
d3ce32b8 750 //if this is a good track, accept track will return true. The track matched is good, so not track matched is false
584f2478 751 nottrackmatched = !(fEsdtrackCutsTPC->AcceptTrack(track));//if the track is bad, this is track matched
6a152780 752 //if(!nottrackmatched) cout<<"Matched track p: "<<matchedTrackp<<" sim "<<part->P()<<endl;
d3ce32b8 753 }
754
755
a0dcb3c1 756 for(UInt_t i = 0; i < caloCluster->GetNLabels(); i++)
b2c10007 757 {
53ad5ec1 758 Int_t pIdx = caloCluster->GetLabelAt(i);
5881f036 759
0861cc1f 760 //TParticle *p = stack->Particle(pIdx);
761
53ad5ec1 762 if(!stack->IsPhysicalPrimary(pIdx))
763 {
0861cc1f 764// PrintFamilyTree(pIdx, stack);
53ad5ec1 765 pIdx = GetPrimMother(pIdx, stack);
766 }
0861cc1f 767 if(fSelector->PassDistanceToBadChannelCut(*caloCluster))//&&fSelector->CutGeometricalAcceptance(*(stack->Particle(primIdx))))
768 {
d3ce32b8 769
0861cc1f 770// std::cout << "Gamma primary: " << primIdx << std::endl;
771// foundGammas.push_back(primIdx);
d3ce32b8 772 if(nottrackmatched){
773 foundGammas.push_back(pIdx);
774 }
0861cc1f 775 }
b2c10007 776 }
53ad5ec1 777 fCutFlow->Fill(cf++);
86e7d5db 778 if(!fSelector->PassDistanceToBadChannelCut(*caloCluster)) continue;
5881f036 779 Double_t clEt = CorrectForReconstructionEfficiency(*caloCluster,fCentClass);
ef647350 780// if(code == fgK0SCode) std::cout << "K0 energy: " << caloCluster->E() << std::endl;
86e7d5db 781 if(!fSelector->PassMinEnergyCut(*caloCluster)) continue;
53ad5ec1 782
783
f61cec2f 784 fCutFlow->Fill(cf++);
0861cc1f 785 Float_t pos[3];
786 //PrintFamilyTree(
787 caloCluster->GetPosition(pos);
788 TVector3 cp(pos);
789
f61cec2f 790 TParticle *primPart = stack->Particle(primIdx);
791 fPrimaryCode = primPart->GetPdgCode();
03bfd268 792 fPrimaryCharge = (Int_t) primPart->GetPDG()->Charge();
f61cec2f 793
794 fPrimaryE = primPart->Energy();
795 fPrimaryEt = primPart->Energy()*TMath::Sin(primPart->Theta());
f61cec2f 796 fPrimaryPx = primPart->Px();
797 fPrimaryPy = primPart->Py();
798 fPrimaryPz = primPart->Pz();
cfe23ff0 799 //cout<<"I have a cluster and it's good energy "<<caloCluster->E()<<" simulated "<<fPrimaryE<<endl;
f61cec2f 800 fPrimaryVx = primPart->Vx();
801 fPrimaryVy = primPart->Vy();
802 fPrimaryVz = primPart->Vz();
803
804 fPrimaryAccepted = false;
b2c10007 805 fPrimaryMatched = false;
f61cec2f 806
807 fDepositedCode = part->GetPdgCode();
b2c10007 808 fDepositedE = part->Energy();
809 fDepositedEt = part->Energy()*TMath::Sin(part->Theta());
03bfd268 810 fDepositedCharge = (Int_t) part->GetPDG()->Charge();
f61cec2f 811
812 fDepositedVx = part->Vx();
813 fDepositedVy = part->Vy();
814 fDepositedVz = part->Vz();
b2c10007 815
816 //fSecondary = fSelector->FromSecondaryInteraction(*primPart, *stack);
817 fSecondary =fSelector->FromSecondaryInteraction(*part, *stack);
9ef6f13f 818// if(fSecondary)
819// {
820// //std::cout << "Have secondary!" << std::endl;
821// //PrintFamilyTree(iPart, stack);
822// }
b2c10007 823 fReconstructedE = caloCluster->E();
824 fReconstructedEt = caloCluster->E()*TMath::Sin(cp.Theta());
b2c10007 825
826 pdg = primPart->GetPDG(0);
a0dcb3c1 827 //Int_t code = primPart->GetPdgCode();
e2ee5727 828
9ef6f13f 829 Bool_t written = kFALSE;
830
d3ce32b8 831// Bool_t nottrackmatched = kTRUE;//default to no track matched
832// nottrackmatched = fSelector->PassTrackMatchingCut(*caloCluster);
833// //by default ALL matched tracks are accepted, whether or not the match is good. So we check to see if the track is good.
834// if(!nottrackmatched){
835// Int_t trackMatchedIndex = caloCluster->GetTrackMatchedIndex();
836// if(trackMatchedIndex < 0) nottrackmatched=kTRUE;
837// AliESDtrack *track = realEvent->GetTrack(trackMatchedIndex);
838// //if this is a good track, accept track will return true. The track matched is good, so not track matched is false
839// nottrackmatched = !(fEsdtrackCutsTPC->AcceptTrack(track));
840// }
9ef6f13f 841
842 if(fSecondary){//all particles from secondary interactions
843 written = kTRUE;
844 if(nottrackmatched){//secondaries not removed
845 if (fDepositedCharge != 0){//charged track not removed
846 fChargedNotRemoved++;
847 fEnergyChargedNotRemoved += clEt;
848 fHistRemovedOrNot->Fill(2.0, fCentClass);
849 fHistDxDzNonRemovedCharged->Fill(caloCluster->GetTrackDx(), caloCluster->GetTrackDz());
850 }
851 else{
852 fSecondaryNotRemoved++;
853 }
854 }
855 else{//secondaries removed
856 if (fDepositedCharge != 0){
857 fChargedRemoved++;
858 fEnergyChargedRemoved += clEt;
859 fHistRemovedOrNot->Fill(0.0, fCentClass);
860 fHistDxDzRemovedCharged->Fill(caloCluster->GetTrackDx(), caloCluster->GetTrackDz());
861 fHistChargedTracksCut->Fill(fDepositedEt);
0861cc1f 862 if(fCalcTrackMatchVsMult){
863 fHistChargedTracksCutMult->Fill(fDepositedEt,fClusterMult);
864 if(fClusterMult<25){fHistChargedTracksCutPeripheral->Fill(fDepositedEt);}
865 }
6a152780 866 fHistMatchedTracksEvspTBkgd->Fill(matchedTrackp,fReconstructedE);
0861cc1f 867 if(fCalcTrackMatchVsMult){
6a152780 868 if(fClusterMult<25){fHistMatchedTracksEvspTBkgdPeripheral->Fill(matchedTrackp,fReconstructedEt);}
ac610b08 869 fHistMatchedTracksEvspTBkgdvsCent->Fill(matchedTrackp,fReconstructedEt, fCentClass);
870 fHistMatchedTracksEvspTBkgdvsCentEffCorr->Fill(matchedTrackp,clEt, fCentClass);//Fill with the efficiency corrected energy
0861cc1f 871 }
9ef6f13f 872 Int_t trackindex = (caloCluster->GetLabelsArray())->At(1);
873 if(caloCluster->GetLabel()!=trackindex){
874 fHistBadTrackMatches->Fill(part->Pt(),fReconstructedE);
875 fHistBadTrackMatchesdPhidEta->Fill(caloCluster->GetTrackDx(),caloCluster->GetTrackDz());
876 //cout<<"Track matched, label cluster "<<caloCluster->GetLabel()<<" track "<<trackindex<<endl;
877 }
878 else{
879 fHistGoodTrackMatchesdPhidEta->Fill(caloCluster->GetTrackDx(),caloCluster->GetTrackDz());
880 }
881 }
882 else{//neutral energy removed
883 fNeutralRemoved++;
884 fEnergyNeutralRemoved += clEt;
885 fHistRemovedOrNot->Fill(1.0, fCentClass);
886 fHistDxDzRemovedNeutral->Fill(caloCluster->GetTrackDx(), caloCluster->GetTrackDz());
887 }
888 }
889 }
890 else{
891
892 if (fDepositedCharge != 0 && fDepositedCode!=fgEMinusCode && fDepositedCode!=fgEPlusCode){//if the particle hitting the calorimeter is pi/k/p/mu
893 written = kTRUE;
d3ce32b8 894 fClusterMultChargedTracks++;
9ef6f13f 895 if(nottrackmatched){//not removed but should be
c79c36cc 896 fHistHadronDepositsAll->Fill(part->Pt());
ac610b08 897 fHistHadronDepositsAllCent->Fill(part->Pt(), fCentClass);
5881f036 898 if(fReconstructedEt>0.5) fHistHadronDepositsAllCent500MeV->Fill(part->Pt(), fCentClass);
9ef6f13f 899 fChargedNotRemoved++;
900 fEnergyChargedNotRemoved += clEt;
901 fHistRemovedOrNot->Fill(2.0, fCentClass);
902 fHistDxDzNonRemovedCharged->Fill(caloCluster->GetTrackDx(), caloCluster->GetTrackDz());
903 fHistChargedTracksAccepted->Fill(fDepositedEt);
0861cc1f 904 if(fCalcTrackMatchVsMult){
6a152780 905 if(matchedTrackpt<0.5){//if we could never have matched this because of its pt, how much energy did it deposit?
ac610b08 906 fHistChargedTracksAcceptedLowPtCent->Fill(fDepositedEt, fCentClass);
5881f036 907 if(fDepositedEt>=0.5) fHistChargedTracksAcceptedLowPtCent500MeV->Fill(fDepositedEt, fCentClass);
ac610b08 908 if(pdg->PdgCode()!=fgAntiProtonCode){
909 fHistChargedTracksAcceptedLowPtCentNoAntiProtons->Fill(fDepositedEt, fCentClass);
910 }
911 else{
912 fHistChargedTracksAcceptedLowPtCentAntiProtons->Fill(fDepositedEt, fCentClass);
913 }
6a152780 914 }
0861cc1f 915 fHistChargedTracksAcceptedMult->Fill(fDepositedEt,fClusterMult);
916 if(fClusterMult<25){fHistChargedTracksAcceptedPeripheral->Fill(fDepositedEt);}
917 }
9ef6f13f 918 }
919 else{//removed and should have been
920 Int_t trackindex = (caloCluster->GetLabelsArray())->At(0);
c79c36cc 921 fHistHadronDepositsReco->Fill(part->Pt());
ac610b08 922 fHistHadronDepositsRecoCent->Fill(part->Pt(), fCentClass);
c79c36cc 923 fHistHadronDepositsAll->Fill(part->Pt());
ac610b08 924 fHistHadronDepositsAllCent->Fill(part->Pt(), fCentClass);
5881f036 925 if(fReconstructedEt>0.5) fHistHadronDepositsAllCent500MeV->Fill(part->Pt(), fCentClass);
9ef6f13f 926 if(caloCluster->GetLabel()!=trackindex){
927 fHistBadTrackMatches->Fill(part->Pt(),fReconstructedE);
928 fHistBadTrackMatchesdPhidEta->Fill(caloCluster->GetTrackDx(),caloCluster->GetTrackDz());
929 //cout<<"Track matched, label cluster "<<caloCluster->GetLabel()<<" track "<<trackindex<<endl;
930 }
931 else{
932 fHistGoodTrackMatchesdPhidEta->Fill(caloCluster->GetTrackDx(),caloCluster->GetTrackDz());
933 }
934 fChargedRemoved++;
935 fEnergyChargedRemoved += clEt;
936 fHistRemovedOrNot->Fill(0.0, fCentClass);
937 fHistDxDzRemovedCharged->Fill(caloCluster->GetTrackDx(), caloCluster->GetTrackDz());
938 fHistChargedTracksCut->Fill(fDepositedEt);
0861cc1f 939 if(fCalcTrackMatchVsMult){
940 fHistChargedTracksCutMult->Fill(fDepositedEt,fClusterMult);
941 if(fClusterMult<25){fHistChargedTracksCutPeripheral->Fill(fDepositedEt);}
942 }
6a152780 943 fHistMatchedTracksEvspTBkgd->Fill(matchedTrackp,fReconstructedE);
0861cc1f 944 if(fCalcTrackMatchVsMult){
6a152780 945 if(fClusterMult<25){fHistMatchedTracksEvspTBkgdPeripheral->Fill(matchedTrackp,fReconstructedEt);}
ac610b08 946 fHistMatchedTracksEvspTBkgdvsCent->Fill(matchedTrackp,fReconstructedEt, fCentClass);
947 fHistMatchedTracksEvspTBkgdvsCentEffCorr->Fill(matchedTrackp,clEt, fCentClass);//fill with the efficiency corrected energy
0861cc1f 948 }
9ef6f13f 949 }
950 }
951 //K0L and any neutral particles from the decay of K+/- or K0S
952 if(!written && (fPrimaryCode==fgKPlusCode || fPrimaryCode==fgKMinusCode || fPrimaryCode==fgK0SCode ||fPrimaryCode==fgK0LCode)){
953 written = kTRUE;//At this point we are not tracking them but we don't count them as neutrals accidentally removed
cfe23ff0 954 }
9ef6f13f 955
956 if(!written && (fDepositedCode==fgGammaCode || fDepositedCode==fgEMinusCode || fDepositedCode ==fgEPlusCode)){//if the particle hitting the calorimeter is gamma, electron and not from a kaon
d3ce32b8 957 fClusterMultGammas++;
9ef6f13f 958 written = kTRUE;
959 if(nottrackmatched){//Not removed and not supposed to be removed - signal
960 fEtNonRemovedGammas += clEt;
961 fMultNonRemovedGammas++;
962 fNeutralNotRemoved--;
963 fEnergyNeutralNotRemoved -= clEt;
964 fHistGammasAccepted->Fill(fDepositedEt);
0861cc1f 965 if(fCalcTrackMatchVsMult){
966 fHistGammasAcceptedMult->Fill(fDepositedEt,fClusterMult);
967 if(fClusterMult<25){fHistGammasAcceptedPeripheral->Fill(fDepositedEt);}
968 }
9ef6f13f 969 }
970 else{//removed but shouldn't have been
971 Int_t trackindex = (caloCluster->GetLabelsArray())->At(1);
972 if(caloCluster->GetLabel()!=trackindex){
973 fHistBadTrackMatches->Fill(part->Pt(),fReconstructedE);
974 fHistBadTrackMatchesdPhidEta->Fill(caloCluster->GetTrackDx(),caloCluster->GetTrackDz());
975// cout<<"Track matched, label cluster "<<caloCluster->GetLabel()<<" track "<<trackindex<<endl;
976// PrintFamilyTree(trackindex, stack);
977// cout<<"Cluster"<<endl;
978 }
979 fGammaRemoved++;
980 fGammaRemovedEt+=clEt;
981 fHistGammasCut->Fill(fDepositedEt);
0861cc1f 982 if(fCalcTrackMatchVsMult){
983 fHistGammasCutMult->Fill(fDepositedEt,fClusterMult);
984 if(fClusterMult<25){fHistGammasCutPeripheral->Fill(fDepositedEt);}
985 }
6a152780 986 fHistMatchedTracksEvspTSignal->Fill(matchedTrackp,fReconstructedE);
0861cc1f 987 if(fCalcTrackMatchVsMult){
6a152780 988 if(fClusterMult<25){fHistMatchedTracksEvspTSignalPeripheral->Fill(matchedTrackp,fReconstructedEt);}
ac610b08 989 fHistMatchedTracksEvspTSignalvsCent->Fill(matchedTrackp,fReconstructedEt, fCentClass);
990 fHistMatchedTracksEvspTSignalvsCentEffCorr->Fill(matchedTrackp,clEt, fCentClass);
0861cc1f 991 }
9ef6f13f 992 }
993 }
994 //all other cases - neutron, anti-neutron, not aware of other cases
995 if(!written){
996 fNeutralNotRemoved++;
997 fEnergyNeutralNotRemoved += clEt;
998 fHistRemovedOrNot->Fill(3.0, fCentClass);
999 fHistDxDzNonRemovedNeutral->Fill(caloCluster->GetTrackDx(), caloCluster->GetTrackDz());
1000 }
1001 }
f61cec2f 1002 fPrimaryTree->Fill();
9ef6f13f 1003 } // end of loop over clusters
e2ee5727 1004
b2c10007 1005
1006 std::sort(foundGammas.begin(), foundGammas.end());
1007 for (Int_t iPart = 0; iPart < stack->GetNtrack(); iPart++)
1008 {
1009
1010 if(!stack->IsPhysicalPrimary(iPart)) continue;
1011
1012 TParticle *part = stack->Particle(iPart);
1013
1014 if (!part)
1015 {
1016 Printf("ERROR: Could not get particle %d", iPart);
1017 continue;
1018 }
1019 TParticlePDG *pdg = part->GetPDG(0);
1020
1021 if (!pdg)
1022 {
1023 Printf("ERROR: Could not get particle PDG %d", iPart);
1024 continue;
1025 }
1026
1027 if(pdg->PdgCode()==fgGammaCode && fSelector->CutGeometricalAcceptance(*part))// TMath::Abs(part->Eta()) < 0.12)
1028 {
1029 fHistGammasGenerated->Fill(part->Energy());
5881f036 1030 fHistGammasGeneratedCent->Fill(part->Energy(),fCentClass);
b2c10007 1031 if(std::binary_search(foundGammas.begin(),foundGammas.end(),iPart))
1032 {
1033 fHistGammasFound->Fill(part->Energy());
5881f036 1034 fHistGammasFoundCent->Fill(part->Energy(),fCentClass);
b2c10007 1035 }
b2c10007 1036 }
6a152780 1037 if(pdg->PdgCode()==fgPiPlusCode || pdg->PdgCode()==fgPiMinusCode || pdg->PdgCode()==fgProtonCode || pdg->PdgCode()==fgAntiProtonCode){//section here for all hadrons generated
ac610b08 1038 fHistHadronsAllCent->Fill(part->Pt(), fCentClass);
6a152780 1039 }
1040
b2c10007 1041
1042 }
0861cc1f 1043 if(fCalcForKaonCorrection){
1044 Float_t etCuts[11] = {0.0,0.05,0.1,0.15,0.2,0.25, 0.3,0.35,0.4,0.45,0.5};
1045 Int_t nEtCuts = 11;
1046 //loop over simulated particles in order to find K0S
1047 for (Int_t iPart = 0; iPart < stack->GetNtrack(); iPart++){
1048 TParticle *part = stack->Particle(iPart);
1049 if (!part){
1050 //Printf("ERROR: Could not get particle %d", iPart);
1051 continue;
1052 }
1053 TParticlePDG *pdg = part->GetPDG(0);
1054 if (!pdg){
1055 //Printf("ERROR: Could not get particle PDG %d", iPart);
1056 continue;
1057 }
1058 //if(stack->IsPhysicalPrimary(iPart)){//if it is a K0 it might have decayed into four pions
2aab9269 1059 //fgK0SCode, fgGammaCode, fgPi0Code
1060 Int_t code = pdg->PdgCode();
0861cc1f 1061 if(code == fgK0SCode || code==fgKPlusCode || code==fgKMinusCode ||code==fgK0LCode || code==fgK0Code){//this is a kaon
1062 //cout<<"I am a kaon too! "<<stack->Particle(iPart)->GetName()<<" "<<code<<endl;
1063 Float_t pTk = stack->Particle(iPart)->Pt();
1064 if(TMath::Abs(stack->Particle(iPart)->Y())<0.5 && stack->IsPhysicalPrimary(iPart)){//these are particles which would be included in our spectra measurements
1065 fHistSimKaonsInAcceptance->Fill(pTk);
1066 if(code == fgK0SCode){fHistSimK0SInAcceptance->Fill(pTk);}
1067 if(code == fgK0LCode){fHistSimK0LInAcceptance->Fill(pTk);}
1068 if(code == fgKPlusCode){fHistSimKPlusInAcceptance->Fill(pTk);}
1069 if(code == fgKMinusCode){fHistSimKMinusInAcceptance->Fill(pTk);}
1070 if(code == fgK0Code){//Split K0's between the two
1071 fHistSimK0SInAcceptance->Fill(pTk,0.5);
1072 fHistSimK0LInAcceptance->Fill(pTk,0.5);
1073 }
1074 }
1075 else{
1076 fHistSimKaonsOutOfAcceptance->Fill(pTk);
1077 // if(!stack->IsPhysicalPrimary(iPart)){
1078 // PrintFamilyTree(iPart, stack);
1079 // }
1080 }
1081 Float_t totalGammaEts[11] = {0.0,0.0,0.0,0.0,0.0, 0.0,0.0,0.0,0.0,0.0, 0.0};
1082 Float_t totalClusterEts[11] = {0.0,0.0,0.0,0.0,0.0, 0.0,0.0,0.0,0.0,0.0, 0.0};
1083 for (int iCluster = 0; iCluster < nCluster; iCluster++ ){//if this cluster is from any of the decay daughters of any kaon... but there is no easy way to look at this so we loop over clusters...
1084 AliESDCaloCluster* caloCluster = ( AliESDCaloCluster* )caloClusters->At( iCluster );
5881f036 1085 if (!fSelector->CutGeometricalAcceptance(*caloCluster)) continue;
1086 const Int_t myPart = TMath::Abs(caloCluster->GetLabel());
0861cc1f 1087 //identify the primary particle which created this cluster
1088 int primIdx = myPart;
1089 if (!stack->IsPhysicalPrimary(myPart)){
1090 primIdx = GetPrimMother(iPart, stack);
1091 } // end of primary particle check
29b7ba7f 1092 TParticle *hitPart = stack->Particle(myPart);
1093 Bool_t hitsAsChargedKaon = kFALSE;
1094 if(hitPart->GetPdgCode()== fgKPlusCode || hitPart->GetPdgCode()== fgKPlusCode){
1095 if(myPart==primIdx){
1096 //The particle hits as a charged kaon and that kaon is a primary kaon - do not count because this is counted in the hadronic correction!
1097 hitsAsChargedKaon = kTRUE;
1098 //cout<<"Found primary charged kaon cluster!"<<endl;
1099 }
1100 }
1101 if(primIdx==iPart && primIdx>0 && !hitsAsChargedKaon){//This cluster is from our primary particle and our primary particle is a kaon
0861cc1f 1102 //cout<<"I have a particle match! prim code"<<code<<" id "<<primIdx <<endl;
1103 Float_t pos[3];
1104 caloCluster->GetPosition(pos);
1105 TVector3 cp(pos);
1106 Double_t clEt = caloCluster->E()*TMath::Sin(cp.Theta());
5881f036 1107 Double_t clEtCorr = CorrectForReconstructionEfficiency(*caloCluster,fCentClass);
2aab9269 1108 for(int l=0;l<nEtCuts;l++){//loop over cut values
0861cc1f 1109 if(clEt>=etCuts[l]){
1110 //cout<<", "<<clEt<<">="<<etCuts[l];
1111 totalClusterEts[l] += clEtCorr;//if cluster et is above the cut off energy add it
1112 totalGammaEts[l] += clEt;//if cluster et is above the cut off energy add it
2aab9269 1113 }
1114 }
2aab9269 1115 }
1116 }
0861cc1f 1117 // cout<<"Deposits: pT: "<<pTk;
1118 // for(int l=0;l<nEtCuts;l++){//loop over cut values
1119 // cout<<" "<<totalClusterEts[l];
1120 // }
1121 // cout<<endl;
1122 if(TMath::Abs(stack->Particle(iPart)->Y())<0.5 && stack->IsPhysicalPrimary(iPart)){//within the acceptance of our spectra and is a primary particle
1123 if(totalClusterEts[0]>0.0){fHistSimKaonsInAcceptanceWithDepositsPrimaries->Fill(pTk);}
1124 //cout<<"I have a particle match! prim code"<<code<<" id "<<iPart <<endl;
1125 for(int l=0;l<nEtCuts;l++){
d3ce32b8 1126 fHistK0EDepositsVsPtInAcceptance->Fill(pTk,totalClusterEts[l],etCuts[l]+0.001);
1127 fHistK0EGammaVsPtInAcceptance->Fill(pTk,totalGammaEts[l],etCuts[l]+0.001);
0861cc1f 1128 }
1129 }
1130 else{//outside the acceptance of our spectra
1131 if(totalClusterEts[0]>0.0){
1132 if(stack->IsPhysicalPrimary(iPart)){fHistSimKaonsOutOfAcceptanceWithDepositsPrimaries->Fill(pTk);}
1133 else{fHistSimKaonsOutOfAcceptanceWithDepositsSecondaries->Fill(pTk);}
1134 }
1135 for(int l=0;l<nEtCuts;l++){
d3ce32b8 1136 fHistK0EDepositsVsPtOutOfAcceptance->Fill(pTk,totalClusterEts[l],etCuts[l]+0.001);
1137 fHistK0EGammaVsPtOutOfAcceptance->Fill(pTk,totalGammaEts[l],etCuts[l]+0.001);
0861cc1f 1138 }
1139 }
1140
2aab9269 1141 }
1142 }
1143 }
d3ce32b8 1144 fHistMultChVsSignalVsMult->Fill(fClusterMultChargedTracks,fClusterMultGammas,fClusterMult);
e2ee5727 1145 FillHistograms();
1146 return 0;
1147}
2fbf38ac 1148
1149void AliAnalysisEtMonteCarlo::Init()
f61cec2f 1150{ // init
e2ee5727 1151 AliAnalysisEt::Init();
2fbf38ac 1152}
13b0d3c1 1153
ce546038 1154void AliAnalysisEtMonteCarlo::ResetEventValues()
f61cec2f 1155{ // reset event values
e2ee5727 1156 AliAnalysisEt::ResetEventValues();
c79c36cc 1157 if(!fIsMC) return;
e2ee5727 1158
f61cec2f 1159 fTotEtSecondary = 0;
1160 fTotEtSecondaryFromEmEtPrimary = 0;
1161 fTotEtWithSecondaryRemoved = 0;
1162
e2ee5727 1163 // collision geometry defaults for p+p:
1164 fImpactParameter = 0;
1165 fNcoll = 1;
1166 fNpart = 2;
1167
1168 fEtNonRemovedProtons = 0;
1169 fEtNonRemovedAntiProtons = 0;
1170 fEtNonRemovedPiPlus = 0;
1171 fEtNonRemovedPiMinus = 0;
1172 fEtNonRemovedKaonPlus = 0;
1173 fEtNonRemovedKaonMinus = 0;
ef647350 1174 fEtNonRemovedK0S = 0;
1175 fEtNonRemovedK0L = 0;
e2ee5727 1176 fEtNonRemovedLambdas = 0;
1177 fEtNonRemovedElectrons = 0;
1178 fEtNonRemovedPositrons = 0;
1179 fEtNonRemovedMuPlus = 0;
1180 fEtNonRemovedMuMinus = 0;
1181 fEtNonRemovedNeutrons = 0;
1182 fEtNonRemovedAntiNeutrons = 0;
1183 fEtNonRemovedGammas = 0;
1184 fEtNonRemovedGammasFromPi0 = 0;
1185
ef647350 1186 fEtRemovedProtons = 0;
1187 fEtRemovedAntiProtons = 0;
1188 fEtRemovedPiPlus = 0;
1189 fEtRemovedPiMinus = 0;
1190 fEtRemovedKaonPlus = 0;
1191 fEtRemovedKaonMinus = 0;
1192 fEtRemovedK0s = 0;
1193 fEtRemovedK0L = 0;
1194 fEtRemovedLambdas = 0;
1195 fEtRemovedElectrons = 0;
1196 fEtRemovedPositrons = 0;
1197 fEtRemovedMuPlus = 0;
1198 fEtRemovedMuMinus = 0;
1199 fEtRemovedNeutrons = 0;
1200
1201 fEtRemovedGammasFromPi0 = 0;
e2ee5727 1202 fEtRemovedGammas = 0;
1203 fEtRemovedNeutrons = 0;
1204 fEtRemovedAntiNeutrons = 0;
1205
1206 fMultNonRemovedProtons = 0;
1207 fMultNonRemovedAntiProtons = 0;
1208 fMultNonRemovedPiPlus = 0;
1209 fMultNonRemovedPiMinus = 0;
1210 fMultNonRemovedKaonPlus = 0;
1211 fMultNonRemovedKaonMinus = 0;
1212 fMultNonRemovedK0s = 0;
ef647350 1213 fMultNonRemovedK0L = 0;
e2ee5727 1214 fMultNonRemovedLambdas = 0;
1215 fMultNonRemovedElectrons = 0;
1216 fMultNonRemovedPositrons = 0;
1217 fMultNonRemovedMuPlus = 0;
1218 fMultNonRemovedMuMinus = 0;
1219 fMultNonRemovedNeutrons = 0;
1220 fMultNonRemovedAntiNeutrons = 0;
1221 fMultNonRemovedGammas = 0;
1222
ef647350 1223 fMultRemovedProtons = 0;
1224 fMultRemovedAntiProtons = 0;
1225 fMultRemovedPiPlus = 0;
1226 fMultRemovedPiMinus = 0;
1227 fMultRemovedKaonPlus = 0;
1228 fMultRemovedKaonMinus = 0;
1229 fMultRemovedK0s = 0;
1230 fMultRemovedK0L = 0;
1231 fMultRemovedLambdas = 0;
1232 fMultRemovedElectrons = 0;
1233 fMultRemovedPositrons = 0;
1234 fMultRemovedMuPlus = 0;
1235 fMultRemovedMuMinus = 0;
1236
e2ee5727 1237 fMultRemovedGammas = 0;
1238 fMultRemovedNeutrons = 0;
1239 fMultRemovedAntiNeutrons = 0;
1240
ef647350 1241 fEnergyChargedNotRemoved = 0;
1242 fEnergyChargedRemoved = 0;
1243 fEnergyNeutralNotRemoved = 0;
1244 fEnergyNeutralRemoved = 0;
f61cec2f 1245
ef647350 1246 fChargedNotRemoved = 0;
1247 fChargedRemoved = 0;
1248 fNeutralNotRemoved = 0;
1249 fNeutralRemoved = 0;
b2c10007 1250 fGammaRemoved = 0;
1251 fSecondaryNotRemoved = 0;
f61cec2f 1252
e2ee5727 1253 fTrackMultInAcc = 0;
1254
f61cec2f 1255 fTotNeutralEtAfterMinEnergyCut = 0;
b2c10007 1256
1257 fSecondaryNotRemoved = 0;
1258
1259 fTotPx = 0;
1260 fTotPy = 0;
1261 fTotPz = 0;
1262
1263
ce546038 1264}
1265
1266void AliAnalysisEtMonteCarlo::CreateHistograms()
f61cec2f 1267{ // histogram related additions
e2ee5727 1268 AliAnalysisEt::CreateHistograms();
c79c36cc 1269 if(!fIsMC) return;
f61cec2f 1270 if (fEventSummaryTree) {
1271 fEventSummaryTree->Branch("fImpactParameter",&fImpactParameter,"fImpactParameter/D");
1272 fEventSummaryTree->Branch("fNcoll",&fNcoll,"fNcoll/I");
1273 fEventSummaryTree->Branch("fNpart",&fNpart,"fNpart/I");
1274 fEventSummaryTree->Branch("fTotEtWithSecondaryRemoved", &fTotEtWithSecondaryRemoved, "fTotEtWithSecondaryRemoved/D");
1275 fEventSummaryTree->Branch("fTotEtSecondaryFromEmEtPrimary", &fTotEtSecondaryFromEmEtPrimary, "fTotEtSecondaryFromEmEtPrimary/D");
1276 fEventSummaryTree->Branch("fTotEtSecondary", &fTotEtSecondary, "fTotEtSecondary/D");
1277 fEventSummaryTree->Branch("fTotNeutralEtAfterMinEnergyCut", &fTotNeutralEtAfterMinEnergyCut, "fTotNeutralEtAfterMinEnergyCut/D");
b2c10007 1278 fEventSummaryTree->Branch("fSecondaryNotRemoved", &fSecondaryNotRemoved, "fSecondaryNotRemoved/I");
1279 fEventSummaryTree->Branch("fChargedNotRemoved", &fChargedNotRemoved, "fChargedNotRemoved/I");
1280 fEventSummaryTree->Branch("fNeutralNotRemoved", &fNeutralNotRemoved, "fNeutralNotRemoved/I");
1281 fEventSummaryTree->Branch("fChargedRemoved", &fChargedRemoved, "fChargedRemoved/I");
1282 fEventSummaryTree->Branch("fNeutralRemoved", &fNeutralRemoved, "fNeutralRemoved/I");
1283 fEventSummaryTree->Branch("fGammaRemoved", &fGammaRemoved, "fGammaRemoved/I");
1284 fEventSummaryTree->Branch("fTotPx", &fTotPx, "fTotPx/D");
1285 fEventSummaryTree->Branch("fTotPy", &fTotPy, "fTotPy/D");
1286 fEventSummaryTree->Branch("fTotPz", &fTotPz, "fTotPz/D");
1287// fEventSummaryTree->Branch("f
e2ee5727 1288 }
f61cec2f 1289
e2ee5727 1290 //fHistDecayVertexNonRemovedCharged = new TH3F("fHistDecayVertexNonRemovedCharged","fHistDecayVertexNonRemovedCharged", 500, -470, 30, 500, -300, 300, 40, -20, 20);
1291 //fHistDecayVertexRemovedCharged = new TH3F("fHistDecayVertexRemovedCharged","fHistDecayVertexRemovedCharged", 500, -470, 30, 500, -300, 300, 40, -20, 20);
1292 //fHistDecayVertexNonRemovedNeutral = new TH3F("fHistDecayVertexNonRemovedNeutral","fHistDecayVertexNonRemovedNeutral", 500, -470, 30, 500, -300, 300, 40, -20, 20);
1293 //fHistDecayVertexRemovedNeutral = new TH3F("fHistDecayVertexRemovedNeutral","fHistDecayVertexRemovedNeutral", 500, -470, 30, 500, -300, 300, 40, -20, 20);
1294
1295 fHistRemovedOrNot = new TH2F("fHistRemovedOrNot", "fHistRemovedOrNot", 4, -0.5, 3.5, 10, -0.5, 9.5);
1296
1297 fHistEtNonRemovedProtons = new TH2F("fHistEtNonRemovedProtons", "fHistEtNonRemovedProtons", 1500, 0, 30, 10, -0.5, 9.5);
1298 fHistEtNonRemovedAntiProtons = new TH2F("fHistEtNonRemovedAntiProtons", "fHistEtNonRemovedAntiProtons", 1500, 0, 30, 10, -0.5, 9.5);
1299 fHistEtNonRemovedPiPlus = new TH2F("fHistEtNonRemovedPiPlus", "fHistEtNonRemovedPiPlus", 1500, 0, 30, 10, -0.5, 9.5);
1300 fHistEtNonRemovedPiMinus = new TH2F("fHistEtNonRemovedPiMinus", "fHistEtNonRemovedPiMinus", 1500, 0, 30, 10, -0.5, 9.5);
1301 fHistEtNonRemovedKaonPlus = new TH2F("fHistEtNonRemovedKaonPlus", "fHistEtNonRemovedKaonPlus", 1500, 0, 30, 10, -0.5, 9.5);
1302 fHistEtNonRemovedKaonMinus = new TH2F("fHistEtNonRemovedKaonMinus", "fHistEtNonRemovedKaonMinus", 1500, 0, 30, 10, -0.5, 9.5);
1303 fHistEtNonRemovedK0s = new TH2F("fHistEtNonRemovedK0s", "fHistEtNonRemovedK0s", 1500, 0, 30, 10, -0.5, 9.5);
ef647350 1304 fHistEtNonRemovedK0L = new TH2F("fHistEtNonRemovedK0L", "fHistEtNonRemovedK0L", 1500, 0, 30, 10, -0.5, 9.5);
e2ee5727 1305 fHistEtNonRemovedLambdas = new TH2F("fHistEtNonRemovedLambdas", "fHistEtNonRemovedLambdas", 1500, 0, 30, 10, -0.5, 9.5);
1306 fHistEtNonRemovedElectrons = new TH2F("fHistEtNonRemovedElectrons", "fHistEtNonRemovedElectrons", 1500, 0, 30, 10, -0.5, 9.5);
1307 fHistEtNonRemovedPositrons = new TH2F("fHistEtNonRemovedPositrons", "fHistEtNonRemovedPositrons", 1500, 0, 30, 10, -0.5, 9.5);
1308 fHistEtNonRemovedMuPlus = new TH2F("fHistEtNonRemovedMuPlus", "fHistEtNonRemovedMuPlus", 1500, 0, 30, 10, -0.5, 9.5);
1309 fHistEtNonRemovedMuMinus = new TH2F("fHistEtNonRemovedMuMinus", "fHistEtNonRemovedMuMinus", 1500, 0, 30, 10, -0.5, 9.5);
1310 fHistEtNonRemovedNeutrons = new TH2F("fHistEtNonRemovedNeutrons", "fHistEtNonRemovedNeutrons", 1500, 0, 30, 10, -0.5, 9.5);
1311 fHistEtNonRemovedAntiNeutrons = new TH2F("fHistEtNonRemovedAntiNeutrons", "fHistEtNonRemovedAntiNeutrons", 1500, 0, 30, 10, -0.5, 9.5);
1312
1313 fHistEtNonRemovedGammas = new TH2F("fHistEtNonRemovedGammas", "fHistEtNonRemovedGammas", 1500, 0, 30, 10, -0.5, 9.5);
1314 fHistEtNonRemovedGammasFromPi0 = new TH2F("fHistEtNonRemovedGammasFromPi0", "fHistEtNonRemovedGammasFromPi0", 1500, 0, 30, 10, -0.5, 9.5);
1315
1316 fHistEtRemovedGammas = new TH2F("fHistEtRemovedGammas", "fHistEtRemovedGammas", 1500, 0, 30, 10, -0.5, 9.5);
1317 fHistEtRemovedNeutrons = new TH2F("fHistEtRemovedNeutrons", "fHistEtRemovedNeutrons", 1500, 0, 30, 10, -0.5, 9.5);
1318 fHistEtRemovedAntiNeutrons = new TH2F("fHistEtRemovedAntiNeutrons", "fHistEtRemovedAntiNeutrons", 1500, 0, 30, 10, -0.5, 9.5);
f61cec2f 1319
ef647350 1320 fHistEtRemovedCharged = new TH2F("fHistEtRemovedCharged", "fHistEtRemovedCharged", 1500, 0, 30, 10, -0.5, 9.5);
1321 fHistEtRemovedNeutrals = new TH2F("fHistEtRemovedNeutrals", "fHistEtRemovedNeutrals", 1500, 0, 30, 10, -0.5, 9.5);
e2ee5727 1322
ef647350 1323 fHistEtNonRemovedCharged = new TH2F("fHistEtNonRemovedCharged", "fHistEtNonRemovedCharged", 1500, 0, 30, 10, -0.5, 9.5);
1324 fHistEtNonRemovedNeutrals = new TH2F("fHistEtNonRemovedNeutrals", "fHistEtNonRemovedNeutrals", 1500, 0, 30, 10, -0.5, 9.5);
f61cec2f 1325
e2ee5727 1326 fHistMultNonRemovedProtons = new TH2F("fHistMultNonRemovedProtons", "fHistMultNonRemovedProtons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1327 fHistMultNonRemovedAntiProtons = new TH2F("fHistMultNonRemovedAntiProtons", "fHistMultNonRemovedAntiProtons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1328 fHistMultNonRemovedPiPlus = new TH2F("fHistMultNonRemovedPiPlus", "fHistMultNonRemovedPiPlus", 100, -0.5, 99.5, 10, -0.5, 9.5);
1329 fHistMultNonRemovedPiMinus = new TH2F("fHistMultNonRemovedPiMinus", "fHistMultNonRemovedPiMinus", 100, -0.5, 99.5, 10, -0.5, 9.5);
1330 fHistMultNonRemovedKaonPlus = new TH2F("fHistMultNonRemovedKaonPlus", "fHistMultNonRemovedKaonPlus", 100, -0.5, 99.5, 10, -0.5, 9.5);
1331 fHistMultNonRemovedKaonMinus = new TH2F("fHistMultNonRemovedKaonMinus", "fHistMultNonRemovedKaonMinus", 100, -0.5, 99.5, 10, -0.5, 9.5);
1332 fHistMultNonRemovedK0s = new TH2F("fHistMultNonRemovedK0s", "fHistMultNonRemovedK0s", 100, -0.5, 99.5, 10, -0.5, 9.5);
ef647350 1333 fHistMultNonRemovedK0L = new TH2F("fHistMultNonRemovedK0L", "fHistMultNonRemovedK0L", 100, -0.5, 99.5, 10, -0.5, 9.5);
e2ee5727 1334 fHistMultNonRemovedLambdas = new TH2F("fHistMultNonRemovedLambdas", "fHistMultNonRemovedLambdas", 100, -0.5, 99.5, 10, -0.5, 9.5);
1335 fHistMultNonRemovedElectrons = new TH2F("fHistMultNonRemovedElectrons", "fHistMultNonRemovedElectrons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1336 fHistMultNonRemovedPositrons = new TH2F("fHistMultNonRemovedPositrons", "fHistMultNonRemovedPositrons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1337 fHistMultNonRemovedMuPlus = new TH2F("fHistMultNonRemovedMuPlus", "fHistMultNonRemovedMuPlus", 100, -0.5, 99.5, 10, -0.5, 9.5);
1338 fHistMultNonRemovedMuMinus = new TH2F("fHistMultNonRemovedMuMinus", "fHistMultNonRemovedMuMinus", 100, -0.5, 99.5, 10, -0.5, 9.5);
1339 fHistMultNonRemovedNeutrons = new TH2F("fHistMultNonRemovedNeutrons", "fHistMultNonRemovedNeutrons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1340 fHistMultNonRemovedAntiNeutrons = new TH2F("fHistMultNonRemovedAntiNeutrons", "fHistMultNonRemovedAntiNeutrons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1341
ef647350 1342 fHistMultNonRemovedGammas = new TH2F("fHistMultNonRemovedGammas", "fHistMultNonRemovedGammas", 100, -0.5, 99.5, 100, -0.5, 99.5);
e2ee5727 1343
ef647350 1344 fHistMultRemovedGammas = new TH2F("fHistMultRemovedGammas", "fHistMultRemovedGammas", 100, -0.5, 99.5, 100, -0.5, 99.5);
e2ee5727 1345 fHistMultRemovedNeutrons = new TH2F("fHistMultRemovedNeutrons", "fHistMultRemovedNeutrons", 100, -0.5, 99.5, 10, -0.5, 9.5);
1346 fHistMultRemovedAntiNeutrons = new TH2F("fHistMultRemovedAntiNeutrons", "fHistMultRemovedAntiNeutrons", 100, -0.5, 99.5, 10, -0.5, 9.5);
f61cec2f 1347 /*
1348 fHistMultRemovedCharged = new TH2F("fHistMultRemovedCharged", "fHistMultRemovedCharged", 1500, 0, 30, 10, -0.5, 9.5);
1349 fHistMultRemovedNeutrals = new TH2F("fHistMultRemovedNeutrals", "fHistMultRemovedNeutrals", 1500, 0, 30, 10, -0.5, 9.5);
1350
1351 fHistMultNonRemovedCharged = new TH2F("fHistMultNonRemovedCharged", "fHistMultNonRemovedCharged", 1500, 0, 30, 10, -0.5, 9.5);
1352 fHistMultNonRemovedNeutrals = new TH2F("fHistMultNonRemovedNeutrals", "fHistMultNonRemovedNeutrals", 1500, 0, 30, 10, -0.5, 9.5);*/
ef647350 1353
ef647350 1354
ef647350 1355 fHistMultRemovedCharged = new TH2F("fHistMultRemovedCharged", "fHistMultRemovedCharged", 100, -0.5, 99.5, 100, -0.5, 99.5);
1356 fHistMultRemovedNeutrals = new TH2F("fHistMultRemovedNeutrals", "fHistMultRemovedNeutrals", 100, -0.5, 99.5, 100, -0.5, 99.5);
1357
1358 fHistMultNonRemovedCharged = new TH2F("fHistMultNonRemovedCharged", "fHistMultNonRemovedCharged", 100, -0.5, 99.5, 100, -0.5, 99.5);
1359 fHistMultNonRemovedNeutrals = new TH2F("fHistMultNonRemovedNeutrals", "fHistMultNonRemovedNeutrals", 100, -0.5, 99.5, 100, -0.5, 99.5);
e2ee5727 1360
1361 fHistTrackMultvsNonRemovedCharged = new TH2F("fHistTrackMultvsNonRemovedCharged", "fHistTrackMultvsNonRemovedCharged", 1000, -0.5, 999.5, 100, -0.5, 99.5);
1362 fHistTrackMultvsNonRemovedNeutral = new TH2F("fHistTrackMultvsNonRemovedNeutral", "fHistTrackMultvsNonRemovedNeutral", 1000, -0.5, 999.5, 100, -0.5, 99.5);
1363 fHistTrackMultvsRemovedGamma = new TH2F("fHistTrackMultvsRemovedGamma", "fHistTrackMultvsRemovedGamma", 1000, -0.5, 999.5, 100, -0.5, 99.5);
1364
1365 fHistClusterMultvsNonRemovedCharged = new TH2F("fHistClusterMultvsNonRemovedCharged", "fHistClusterMultvsNonRemovedCharged", 1000, -0.5, 999.5, 100, -0.5, 99.5);
1366 fHistClusterMultvsNonRemovedNeutral = new TH2F("fHistClusterMultvsNonRemovedNeutral", "fHistClusterMultvsNonRemovedNeutral", 1000, -0.5, 999.5, 100, -0.5, 99.5);
1367 fHistClusterMultvsRemovedGamma = new TH2F("fHistClusterMultvsRemovedGamma", "fHistClusterMultvsRemovedGamma", 1000, -0.5, 999.5, 100, -0.5, 99.5);
1368
1369 fHistDxDzNonRemovedCharged = new TH2F("fHistDxDzNonRemovedCharged", "fHistDxDzNonRemovedCharged", 800, -200, 200, 800, -200, 200);
1370 fHistDxDzRemovedCharged = new TH2F("fHistDxDzRemovedCharged", "fHistDxDzRemovedCharged", 800, -200, 200, 800, -200, 200);
1371 fHistDxDzNonRemovedNeutral = new TH2F("fHistDxDzNonRemovedNeutral", "fHistDxDzNonRemovedNeutral", 800, -200, 200, 800, -200, 200);
1372 fHistDxDzRemovedNeutral = new TH2F("fHistDxDzRemovedNeutral", "fHistDxDzRemovedNeutral", 800, -200, 200, 800, -200, 200);
1373
0861cc1f 1374 if(fCalcForKaonCorrection){
1375 Int_t nEtCut = 11;
1376 Float_t etCutAxis[12] = {0.00,0.05,0.10,0.15,0.20, 0.25,0.30,0.35,0.40,0.45, 0.50,0.55};
1377 fHistK0EDepositsVsPtInAcceptance = new TH3F("fHistK0EDepositsVsPtInAcceptance","Kaon deposits with corrections for kaons with y<0.5",fgNumOfPtBins,fgPtAxis,fgNumOfPtBins,fgPtAxis,nEtCut,etCutAxis);
1378 fHistK0EGammaVsPtInAcceptance = new TH3F("fHistK0EGammaVsPtInAcceptance","Kaon deposits without corrections for kaons with y<0.5",fgNumOfPtBins,fgPtAxis,fgNumOfPtBins,fgPtAxis,nEtCut,etCutAxis);
1379 fHistK0EDepositsVsPtOutOfAcceptance = new TH3F("fHistK0EDepositsVsPtOutOfAcceptance","Kaon deposits with corrections for kaons with y>0.5",fgNumOfPtBins,fgPtAxis,fgNumOfPtBins,fgPtAxis,nEtCut,etCutAxis);
1380 fHistK0EGammaVsPtOutOfAcceptance = new TH3F("fHistK0EGammaVsPtOutOfAcceptance","Kaon deposits without corrections for kaons with y>0.5",fgNumOfPtBins,fgPtAxis,fgNumOfPtBins,fgPtAxis,nEtCut,etCutAxis);
1381 fHistK0EDepositsVsPtInAcceptance->GetXaxis()->SetTitle("Kaon p_{T}");
1382 fHistK0EGammaVsPtInAcceptance->GetXaxis()->SetTitle("Kaon p_{T}");
1383 fHistK0EDepositsVsPtOutOfAcceptance->GetXaxis()->SetTitle("Kaon p_{T}");
1384 fHistK0EGammaVsPtOutOfAcceptance->GetXaxis()->SetTitle("Kaon p_{T}");
1385 fHistK0EDepositsVsPtInAcceptance->GetYaxis()->SetTitle("Deposited E_{T}");
1386 fHistK0EGammaVsPtInAcceptance->GetYaxis()->SetTitle("Deposited E_{T}");
1387 fHistK0EDepositsVsPtOutOfAcceptance->GetYaxis()->SetTitle("Deposited E_{T}");
1388 fHistK0EGammaVsPtOutOfAcceptance->GetYaxis()->SetTitle("Deposited E_{T}");
1389 fHistK0EDepositsVsPtInAcceptance->GetZaxis()->SetTitle("E_{T} cut");
1390 fHistK0EGammaVsPtInAcceptance->GetZaxis()->SetTitle("E_{T} cut");
1391 fHistK0EDepositsVsPtOutOfAcceptance->GetZaxis()->SetTitle("E_{T} cut");
1392 fHistK0EGammaVsPtOutOfAcceptance->GetZaxis()->SetTitle("E_{T} cut");
1393
1394 fHistSimKaonsInAcceptance = new TH1F("fHistSimKaonsInAcceptance","Kaons with y<0.5",fgNumOfPtBins,fgPtAxis);
1395 fHistSimK0SInAcceptance = new TH1F("fHistSimK0SInAcceptance","Kaons with y<0.5",fgNumOfPtBins,fgPtAxis);
1396 fHistSimK0LInAcceptance = new TH1F("fHistSimK0LInAcceptance","Kaons with y<0.5",fgNumOfPtBins,fgPtAxis);
1397 fHistSimKPlusInAcceptance = new TH1F("fHistSimKPlusInAcceptance","Kaons with y<0.5",fgNumOfPtBins,fgPtAxis);
1398 fHistSimKMinusInAcceptance = new TH1F("fHistSimKMinusInAcceptance","Kaons with y<0.5",fgNumOfPtBins,fgPtAxis);
1399 fHistSimKaonsOutOfAcceptance = new TH1F("fHistSimKaonsOutOfAcceptance","Kaons with y>0.5",fgNumOfPtBins,fgPtAxis);
1400 fHistSimKaonsInAcceptanceWithDepositsPrimaries = new TH1F("fHistSimKaonsInAcceptanceWithDepositsPrimaries","Primary Kaons which deposited energy in calorimeter with y>0.5",fgNumOfPtBins,fgPtAxis);
1401 fHistSimKaonsOutOfAcceptanceWithDepositsSecondaries = new TH1F("fHistSimKaonsOutOfAcceptanceWithDepositsSecondaries","Secondary Kaons which deposited energy in calorimeter with y>0.5",fgNumOfPtBins,fgPtAxis);
1402 fHistSimKaonsOutOfAcceptanceWithDepositsPrimaries = new TH1F("fHistSimKaonsOutOfAcceptanceWithDepositsPrimaries","Primary Kaons which deposited energy in calorimeter with y>0.5",fgNumOfPtBins,fgPtAxis);
1403 }
2aab9269 1404
e2ee5727 1405 fHistPiPlusMult = new TH1F("fHistPiPlusMult", "fHistPiPlusMult", 2000, -0.5, 1999.5);
1406 fHistPiMinusMult = new TH1F("fHistPiMinusMult", "fHistPiMinusMult", 2000, -0.5, 1999.5);
1407 fHistPiZeroMult = new TH1F("fHistPiZeroMult", "fHistPiZeroMult", 2000, -0.5, 1999.5);
1408
1409 fHistPiPlusMultAcc = new TH1F("fHistPiPlusMultAcc", "fHistPiPlusMultAcc", 2000, -0.5, 1999.5);
1410 fHistPiMinusMultAcc = new TH1F("fHistPiMinusMultAcc", "fHistPiMinusMultAcc", 2000, -0.5, 1999.5);
1411 fHistPiZeroMultAcc = new TH1F("fHistPiZeroMultAcc", "fHistPiZeroMultAcc", 2000, -0.5, 1999.5);
1412
f61cec2f 1413 if(fCuts->GetHistMakeTree())
1414 {
1415 TString treename = "fPrimaryTree" + fHistogramNameSuffix;
1416 fPrimaryTree = new TTree(treename, treename);
1417
1418 fPrimaryTree->Branch("fTotEt",&fTotEt,"fTotEt/D");
1419 fPrimaryTree->Branch("fNeutralMultiplicity", &fNeutralMultiplicity, "fNeutralMultiplicity/I");
1420 fPrimaryTree->Branch("fCentClass",&fCentClass,"fCentClass/I");
1421
1422 fPrimaryTree->Branch("fPrimaryCode", &fPrimaryCode, "fPrimaryCode/I");
1423 fPrimaryTree->Branch("fPrimaryCharge", &fPrimaryCharge, "fPrimaryCharge/I");
1424
1425 fPrimaryTree->Branch("fPrimaryE", &fPrimaryE, "fPrimaryE/D");
1426 fPrimaryTree->Branch("fPrimaryEt", &fPrimaryEt, "fPrimaryEt/D");
1427
1428 fPrimaryTree->Branch("fPrimaryPx", &fPrimaryPx, "fPrimaryPx/D");
1429 fPrimaryTree->Branch("fPrimaryPy", &fPrimaryPy, "fPrimaryPy/D");
1430 fPrimaryTree->Branch("fPrimaryPz", &fPrimaryPz, "fPrimaryPz/D");
1431
1432 fPrimaryTree->Branch("fPrimaryVx", &fPrimaryVx, "fPrimaryVx/D");
1433 fPrimaryTree->Branch("fPrimaryVy", &fPrimaryVy, "fPrimaryVy/D");
1434 fPrimaryTree->Branch("fPrimaryVz", &fPrimaryVz, "fPrimaryVz/D");
1435
1436 fPrimaryTree->Branch("fPrimaryAccepted", &fPrimaryAccepted, "fPrimaryAccepted/B");
b2c10007 1437 fPrimaryTree->Branch("fPrimaryMatched", &fPrimaryMatched, "fPrimaryMatched/B");
f61cec2f 1438
1439
1440 fPrimaryTree->Branch("fDepositedCode", &fDepositedCode, "fDepositedCode/I");
1441 fPrimaryTree->Branch("fDepositedCharge", &fDepositedCharge, "fDepositedCharge/I");
b2c10007 1442 fPrimaryTree->Branch("fDepositedE", &fDepositedE, "fDepositedE/D");
f61cec2f 1443 fPrimaryTree->Branch("fDepositedEt", &fDepositedEt, "fDepositedEt/D");
1444
1445 fPrimaryTree->Branch("fDepositedVx", &fDepositedVx, "fDepositedVx/D");
1446 fPrimaryTree->Branch("fDepositedVy", &fDepositedVy, "fDepositedVy/D");
1447 fPrimaryTree->Branch("fDepositedVz", &fDepositedVz, "fDepositedVz/D");
1448
b2c10007 1449 fPrimaryTree->Branch("fSecondary", &fSecondary, "fSecondary/I");
1450
1451
1452 fPrimaryTree->Branch("fReconstructedE", &fReconstructedE, "fReconstructedE/D");
1453 fPrimaryTree->Branch("fReconstructedEt", &fReconstructedEt, "fReconstructedEt/D");
1454
1455 fPrimaryTree->Branch("fClusterMult", &fClusterMult, "fClusterMult/I");
1456
1457
f61cec2f 1458 }
1459
9ef6f13f 1460 fHistGammasFound = new TH1F("fHistGammasFound", "fHistGammasFound",200, 0, 10);
1461 fHistGammasGenerated = new TH1F("fHistGammasGenerated", "fHistGammasGenerated",200, 0, 10);
5881f036 1462 fHistGammasFoundCent = new TH2F("fHistGammasFoundCent", "fHistGammasFoundCent",200, 0, 10,20,-0.5,19.5);
1463 fHistGammasGeneratedCent = new TH2F("fHistGammasGeneratedCent", "fHistGammasGeneratedCent",200, 0, 10,20,-0.5,19.5);
9ef6f13f 1464 fHistChargedTracksCut = new TH1F("fHistChargedTracksCut", "fHistChargedTracksCut",100, 0, 5);
1465 fHistChargedTracksAccepted = new TH1F("fHistChargedTracksAccepted", "fHistChargedTracksAccepted",100, 0, 5);
1466 fHistGammasCut = new TH1F("fHistGammasTracksCut", "fHistGammasTracksCut",100, 0, 5);
1467 fHistGammasAccepted = new TH1F("fHistGammasTracksAccepted", "fHistGammasTracksAccepted",100, 0, 5);
0861cc1f 1468
1469 if(fCalcTrackMatchVsMult){
1470 fHistChargedTracksCutMult = new TH2F("fHistChargedTracksCutMult", "fHistChargedTracksCutMult",100, 0, 5,10,0,100);
1471 fHistChargedTracksAcceptedMult = new TH2F("fHistChargedTracksAcceptedMult", "fHistChargedTracksAcceptedMult",100, 0, 5,10,0,100);
5881f036 1472 fHistChargedTracksAcceptedLowPtCent = new TH2F("fHistChargedTracksAcceptedLowPtCent", "fHistChargedTracksAcceptedLowPtCent",100, 0, 5,20,-0.5,19.5);
1473 fHistChargedTracksAcceptedLowPtCent500MeV = new TH2F("fHistChargedTracksAcceptedLowPtCent500MeV", "fHistChargedTracksAcceptedLowPtCent500MeV",100, 0, 5,20,-0.5,19.5);
1474 fHistChargedTracksAcceptedLowPtCentNoAntiProtons = new TH2F("fHistChargedTracksAcceptedLowPtCentNoAntiProtons", "fHistChargedTracksAcceptedLowPtCentNoAntiProtons",100, 0, 5,20,-0.5,19.5);
1475 fHistChargedTracksAcceptedLowPtCentAntiProtons = new TH2F("fHistChargedTracksAcceptedLowPtCentAntiProtons", "fHistChargedTracksAcceptedLowPtCentAntiProtons",100, 0, 5,20,-0.5,19.5);
0861cc1f 1476 fHistGammasCutMult = new TH2F("fHistGammasTracksCutMult", "fHistGammasTracksCutMult",100, 0, 5,10,0,100);
1477 fHistGammasAcceptedMult = new TH2F("fHistGammasTracksAcceptedMult", "fHistGammasTracksAcceptedMult",100, 0, 5,10,0,100);
1478 }
1479
9ef6f13f 1480 fHistBadTrackMatches = new TH1F("fHistBadTrackMatches", "fHistBadTrackMatches",100, 0, 5);
1481 fHistMatchedTracksEvspTBkgd = new TH2F("fHistMatchedTracksEvspTBkgd", "fHistMatchedTracksEvspTBkgd",100, 0, 3,100,0,3);
1482 fHistMatchedTracksEvspTSignal = new TH2F("fHistMatchedTracksEvspTSignal", "fHistMatchedTracksEvspTSignal",100, 0, 3,100,0,3);
0861cc1f 1483 if(fCalcTrackMatchVsMult){
1484 fHistMatchedTracksEvspTBkgdPeripheral = new TH2F("fHistMatchedTracksEvspTBkgdPeripheral", "fHistMatchedTracksEvspTBkgd",100, 0, 3,100,0,3);
1485 fHistMatchedTracksEvspTSignalPeripheral = new TH2F("fHistMatchedTracksEvspTSignalPeripheral", "fHistMatchedTracksEvspTSignal",100, 0, 3,100,0,3);
1486
5881f036 1487 fHistMatchedTracksEvspTBkgdvsCent = new TH3F("fHistMatchedTracksEvspTBkgdvsCent", "fHistMatchedTracksEvspTBkgdvsCent",100, 0, 3,100,0,3,20,-0.5,19.5);
1488 fHistMatchedTracksEvspTSignalvsCent = new TH3F("fHistMatchedTracksEvspTSignalvsCent", "fHistMatchedTracksEvspTSignalvsCent",100, 0, 3,100,0,3,20,-0.5,19.5);
1489 fHistMatchedTracksEvspTBkgdvsCentEffCorr = new TH3F("fHistMatchedTracksEvspTBkgdvsCentEffCorr", "fHistMatchedTracksEvspTBkgdvsCent",100, 0, 3,100,0,3,20,-0.5,19.5);
1490 fHistMatchedTracksEvspTSignalvsCentEffCorr = new TH3F("fHistMatchedTracksEvspTSignalvsCentEffCorr", "fHistMatchedTracksEvspTSignalvsCent",100, 0, 3,100,0,3,20,-0.5,19.5);
0861cc1f 1491
1492
1493 fHistChargedTracksCutPeripheral = new TH1F("fHistChargedTracksCutPeripheral", "fHistChargedTracksCut",100, 0, 5);
1494 fHistChargedTracksAcceptedPeripheral = new TH1F("fHistChargedTracksAcceptedPeripheral", "fHistChargedTracksAccepted",100, 0, 5);
1495 fHistGammasCutPeripheral = new TH1F("fHistGammasTracksCutPeripheral", "fHistGammasTracksCut",100, 0, 5);
1496 fHistGammasAcceptedPeripheral = new TH1F("fHistGammasTracksAcceptedPeripheral", "fHistGammasTracksAccepted",100, 0, 5);
1497 }
9ef6f13f 1498 fHistBadTrackMatchesdPhidEta = new TH2F("fHistBadTrackMatchesdPhidEta", "fHistBadTrackMatchesdPhidEta",20, -0.1, 0.1,20,-.1,0.1);
1499 fHistGoodTrackMatchesdPhidEta = new TH2F("fHistGoodTrackMatchesdPhidEta", "fHistGoodTrackMatchesdPhidEta",20, -0.1, 0.1,20,-.1,0.1);
c79c36cc 1500
1501 fHistHadronDepositsAll = new TH1F("fHistHadronDepositsAll","All Hadrons which deposited energy in calorimeter",fgNumOfPtBins,fgPtAxis);
1502 fHistHadronDepositsReco = new TH1F("fHistHadronDepositsReco","Reconstructed Hadrons which deposited energy in calorimeter",fgNumOfPtBins,fgPtAxis);
1503 //,10,0,100
1504 Int_t nMult = 20;
1505 Float_t nMultCuts[21] = { 0, 5,10,15,20, 25,30,35,40,45,
1506 50,55,60,65,70, 75,80,85,90,95,
1507 100};
ac610b08 1508 Int_t nCent = 20;
1509 Float_t nCentCuts[21] = { 0, 1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20};
c79c36cc 1510
ac610b08 1511 fHistHadronDepositsAllCent = new TH2F("fHistHadronDepositsAllCent","All Hadrons which deposited energy in calorimeter",fgNumOfPtBins,fgPtAxis,nCent,nCentCuts);
5881f036 1512 fHistHadronDepositsAllCent500MeV = new TH2F("fHistHadronDepositsAllCent500MeV","All Hadrons which deposited energy in calorimeter pT>500MeV",fgNumOfPtBins,fgPtAxis,nCent,nCentCuts);
ac610b08 1513 fHistHadronDepositsRecoCent = new TH2F("fHistHadronDepositsRecoCent","Reconstructed Hadrons which deposited energy in calorimeter",fgNumOfPtBins,fgPtAxis,nCent,nCentCuts);
d3ce32b8 1514
ac610b08 1515 fHistHadronsAllCent = new TH2F("fHistHadronsAllCent","All Hadrons vs cluster mult",fgNumOfPtBins,fgPtAxis,nCent,nCentCuts);
d3ce32b8 1516
1517 fHistMultChVsSignalVsMult = new TH3F("fHistMultChVsSignalVsMult","Charged particle Multiplicity vs Signal particle multiplicity vs Cluster Mult",nMult,nMultCuts,nMult,nMultCuts,nMult,nMultCuts);
ce546038 1518}
1519
0651f6b4 1520void AliAnalysisEtMonteCarlo::FillOutputList(TList *list)
f61cec2f 1521{ //fill the output list
e2ee5727 1522 AliAnalysisEt::FillOutputList(list);
1523
c79c36cc 1524 if(!fIsMC) return;
f61cec2f 1525 if(fCuts->GetHistMakeTree())
1526 {
1527 list->Add(fPrimaryTree);
1528 }
1529
e2ee5727 1530 list->Add(fHistRemovedOrNot);
1531
1532 list->Add(fHistEtNonRemovedProtons);
1533 list->Add(fHistEtNonRemovedAntiProtons);
1534 list->Add(fHistEtNonRemovedPiPlus);
1535 list->Add(fHistEtNonRemovedPiMinus);
1536 list->Add(fHistEtNonRemovedKaonPlus);
1537 list->Add(fHistEtNonRemovedKaonMinus);
1538 list->Add(fHistEtNonRemovedK0s);
ef647350 1539 list->Add(fHistEtNonRemovedK0L);
e2ee5727 1540 list->Add(fHistEtNonRemovedLambdas);
1541 list->Add(fHistEtNonRemovedElectrons);
1542 list->Add(fHistEtNonRemovedPositrons);
1543 list->Add(fHistEtNonRemovedMuPlus);
1544 list->Add(fHistEtNonRemovedMuMinus);
1545 list->Add(fHistEtNonRemovedNeutrons);
1546 list->Add(fHistEtNonRemovedAntiNeutrons);
1547 list->Add(fHistEtNonRemovedGammas);
1548 list->Add(fHistEtNonRemovedGammasFromPi0);
1549
1550 list->Add(fHistEtRemovedGammas);
1551 list->Add(fHistEtRemovedNeutrons);
1552 list->Add(fHistEtRemovedAntiNeutrons);
1553
ef647350 1554 list->Add(fHistEtRemovedCharged);
1555 list->Add(fHistEtRemovedNeutrals);
1556
1557 list->Add(fHistEtNonRemovedCharged);
1558 list->Add(fHistEtNonRemovedNeutrals);
e2ee5727 1559
1560 list->Add(fHistMultNonRemovedProtons);
1561 list->Add(fHistMultNonRemovedAntiProtons);
1562 list->Add(fHistMultNonRemovedPiPlus);
1563 list->Add(fHistMultNonRemovedPiMinus);
1564 list->Add(fHistMultNonRemovedKaonPlus);
1565 list->Add(fHistMultNonRemovedKaonMinus);
1566 list->Add(fHistMultNonRemovedK0s);
ef647350 1567 list->Add(fHistMultNonRemovedK0L);
e2ee5727 1568 list->Add(fHistMultNonRemovedLambdas);
1569 list->Add(fHistMultNonRemovedElectrons);
1570 list->Add(fHistMultNonRemovedPositrons);
1571 list->Add(fHistMultNonRemovedMuPlus);
1572 list->Add(fHistMultNonRemovedMuMinus);
1573 list->Add(fHistMultNonRemovedNeutrons);
1574 list->Add(fHistMultNonRemovedAntiNeutrons);
1575 list->Add(fHistMultNonRemovedGammas);
1576
1577 list->Add(fHistMultRemovedGammas);
1578 list->Add(fHistMultRemovedNeutrons);
1579 list->Add(fHistMultRemovedAntiNeutrons);
1580
ef647350 1581 list->Add(fHistMultRemovedCharged);
1582 list->Add(fHistMultRemovedNeutrals);
1583
1584 list->Add(fHistMultNonRemovedCharged);
1585 list->Add(fHistMultNonRemovedNeutrals);
1586
e2ee5727 1587 list->Add(fHistTrackMultvsNonRemovedCharged);
1588 list->Add(fHistTrackMultvsNonRemovedNeutral);
1589 list->Add(fHistTrackMultvsRemovedGamma);
1590
1591 list->Add(fHistClusterMultvsNonRemovedCharged);
1592 list->Add(fHistClusterMultvsNonRemovedNeutral);
1593 list->Add(fHistClusterMultvsRemovedGamma);
1594
1595 //list->Add(fHistDecayVertexNonRemovedCharged);
1596 //list->Add(fHistDecayVertexNonRemovedNeutral);
1597 //list->Add(fHistDecayVertexRemovedCharged);
1598 //list->Add(fHistDecayVertexRemovedNeutral);
1599
1600 list->Add(fHistDxDzNonRemovedCharged);
1601 list->Add(fHistDxDzRemovedCharged);
1602 list->Add(fHistDxDzNonRemovedNeutral);
1603 list->Add(fHistDxDzRemovedNeutral);
1604
0861cc1f 1605 if(fCalcForKaonCorrection){
1606 list->Add(fHistK0EDepositsVsPtInAcceptance);
1607 list->Add(fHistK0EGammaVsPtInAcceptance);
1608 list->Add(fHistK0EDepositsVsPtOutOfAcceptance);
1609 list->Add(fHistK0EGammaVsPtOutOfAcceptance);
1610 list->Add(fHistSimKaonsInAcceptance);
1611 list->Add(fHistSimK0SInAcceptance);
1612 list->Add(fHistSimK0LInAcceptance);
1613 list->Add(fHistSimKPlusInAcceptance);
1614 list->Add(fHistSimKMinusInAcceptance);
1615 list->Add(fHistSimKaonsOutOfAcceptance);
1616 list->Add(fHistSimKaonsInAcceptanceWithDepositsPrimaries);
1617 list->Add(fHistSimKaonsOutOfAcceptanceWithDepositsSecondaries);
1618 list->Add(fHistSimKaonsOutOfAcceptanceWithDepositsPrimaries);
1619 }
2aab9269 1620
e2ee5727 1621 list->Add(fHistPiPlusMult);
1622 list->Add(fHistPiMinusMult);
1623 list->Add(fHistPiZeroMult);
1624 list->Add(fHistPiPlusMultAcc);
1625 list->Add(fHistPiMinusMultAcc);
1626 list->Add(fHistPiZeroMultAcc);
b2c10007 1627
1628 list->Add(fHistGammasFound);
1629 list->Add(fHistGammasGenerated);
5881f036 1630 list->Add(fHistGammasFoundCent);
1631 list->Add(fHistGammasGeneratedCent);
9ef6f13f 1632 list->Add(fHistChargedTracksCut);
1633 list->Add(fHistChargedTracksAccepted);
1634 list->Add(fHistGammasCut);
1635 list->Add(fHistGammasAccepted);
0861cc1f 1636 if(fCalcTrackMatchVsMult){
1637 list->Add(fHistChargedTracksCutMult);
1638 list->Add(fHistChargedTracksAcceptedMult);
ac610b08 1639 list->Add(fHistChargedTracksAcceptedLowPtCent);
5881f036 1640 list->Add(fHistChargedTracksAcceptedLowPtCent500MeV);
ac610b08 1641 list->Add(fHistChargedTracksAcceptedLowPtCentNoAntiProtons);
1642 list->Add(fHistChargedTracksAcceptedLowPtCentAntiProtons);
0861cc1f 1643 list->Add(fHistGammasCutMult);
1644 list->Add(fHistGammasAcceptedMult);
1645 }
9ef6f13f 1646 list->Add(fHistBadTrackMatches);
1647 list->Add(fHistMatchedTracksEvspTBkgd);
1648 list->Add(fHistMatchedTracksEvspTSignal);
0861cc1f 1649 if(fCalcTrackMatchVsMult){
1650 list->Add(fHistMatchedTracksEvspTBkgdPeripheral);
1651 list->Add(fHistMatchedTracksEvspTSignalPeripheral);
ac610b08 1652 list->Add(fHistMatchedTracksEvspTBkgdvsCent);
1653 list->Add(fHistMatchedTracksEvspTSignalvsCent);
1654 list->Add(fHistMatchedTracksEvspTBkgdvsCentEffCorr);
1655 list->Add(fHistMatchedTracksEvspTSignalvsCentEffCorr);
0861cc1f 1656 list->Add(fHistChargedTracksCutPeripheral);
1657 list->Add(fHistChargedTracksAcceptedPeripheral);
1658 list->Add(fHistGammasCutPeripheral);
1659 list->Add(fHistGammasAcceptedPeripheral);
1660 }
9ef6f13f 1661 list->Add(fHistBadTrackMatchesdPhidEta);
1662 list->Add(fHistGoodTrackMatchesdPhidEta);
c79c36cc 1663 list->Add(fHistHadronDepositsAll);
1664 list->Add(fHistHadronDepositsReco);
ac610b08 1665 list->Add(fHistHadronDepositsAllCent);
5881f036 1666 list->Add(fHistHadronDepositsAllCent500MeV);
ac610b08 1667 list->Add(fHistHadronDepositsRecoCent);
1668 list->Add(fHistHadronsAllCent);
d3ce32b8 1669 list->Add(fHistMultChVsSignalVsMult);
e2ee5727 1670
0651f6b4 1671}
1672
1673
13b0d3c1 1674bool AliAnalysisEtMonteCarlo::TrackHitsCalorimeter(TParticle* part, Double_t magField)
1675{
e2ee5727 1676 // printf(" TrackHitsCalorimeter - magField %f\n", magField);
1677 AliESDtrack *esdTrack = new AliESDtrack(part);
1678 // Printf("MC Propagating track: eta: %f, phi: %f, pt: %f", esdTrack->Eta(), esdTrack->Phi(), esdTrack->Pt());
1679
1680 Bool_t prop = esdTrack->PropagateTo(fDetectorRadius, magField);
1681
1682 // if(prop) Printf("Track propagated, eta: %f, phi: %f, pt: %f", esdTrack->Eta(), esdTrack->Phi(), esdTrack->Pt());
1683
f61cec2f 1684 bool status = prop && fSelector->CutGeometricalAcceptance(*esdTrack);
e2ee5727 1685 delete esdTrack;
1686
1687 return status;
1688}
1689
1690void AliAnalysisEtMonteCarlo::FillHistograms()
f61cec2f 1691{ // let base class fill its histograms, and us fill the local ones
e2ee5727 1692 AliAnalysisEt::FillHistograms();
c79c36cc 1693 if(!fIsMC) return;
e2ee5727 1694 //std::cout << fEtNonRemovedPiPlus << " " << fCentClass << std::endl;
1695
1696 fHistEtNonRemovedProtons->Fill(fEtNonRemovedProtons, fCentClass);
1697 fHistEtNonRemovedAntiProtons->Fill(fEtNonRemovedAntiProtons, fCentClass);
1698 fHistEtNonRemovedKaonPlus->Fill(fEtNonRemovedKaonPlus, fCentClass);
1699 fHistEtNonRemovedKaonMinus->Fill(fEtNonRemovedKaonMinus, fCentClass);
ef647350 1700 fHistEtNonRemovedK0s->Fill(fEtNonRemovedK0S, fCentClass);
1701 fHistEtNonRemovedK0L->Fill(fEtNonRemovedK0L, fCentClass);
e2ee5727 1702 fHistEtNonRemovedLambdas->Fill(fEtNonRemovedLambdas, fCentClass);
1703 fHistEtNonRemovedPiPlus->Fill(fEtNonRemovedPiPlus, fCentClass);
1704 fHistEtNonRemovedPiMinus->Fill(fEtNonRemovedPiMinus, fCentClass);
1705 fHistEtNonRemovedElectrons->Fill(fEtNonRemovedElectrons, fCentClass);
1706 fHistEtNonRemovedPositrons->Fill(fEtNonRemovedPositrons, fCentClass);
1707 fHistEtNonRemovedMuPlus->Fill(fEtNonRemovedMuPlus, fCentClass);
1708 fHistEtNonRemovedMuMinus->Fill(fEtNonRemovedMuMinus, fCentClass);
1709 fHistEtNonRemovedNeutrons->Fill(fEtNonRemovedNeutrons, fCentClass);
1710 fHistEtNonRemovedAntiNeutrons->Fill(fEtNonRemovedAntiNeutrons, fCentClass);
1711 fHistEtNonRemovedGammas->Fill(fEtNonRemovedGammas, fCentClass);
1712 fHistEtNonRemovedGammasFromPi0->Fill(fEtNonRemovedGammasFromPi0, fCentClass);
1713
ef647350 1714 fHistEtRemovedGammas->Fill(fEtRemovedGammas, fNClusters);
e2ee5727 1715 fHistEtRemovedNeutrons->Fill(fEtRemovedNeutrons, fCentClass);
1716 fHistEtRemovedAntiNeutrons->Fill(fEtRemovedAntiNeutrons, fCentClass);
1717
ef647350 1718// fHistEtRemovedCharged->Fill(fEtRemovedAntiProtons+fEtRemovedElectrons+fEtRemovedKaonMinus+fEtRemovedKaonPlus
1719// +fEtRemovedMuMinus+fEtRemovedMuPlus+fEtRemovedPiMinus+fEtRemovedPiPlus+fEtRemovedPositrons
1720// +fEtRemovedProtons.
1721// fCentClass);
1722// fHistEtRemovedNeutrals->Fill(fEtRemovedNeutrons+fEtRemovedAntiNeutrons, fCentClass);
f61cec2f 1723//
ef647350 1724// fHistEtNonRemovedCharged->Fill(fEtNonRemovedAntiProtons+fEtNonRemovedElectrons+fEtNonRemovedKaonMinus+fEtNonRemovedKaonPlus
1725// +fEtNonRemovedMuMinus+fEtNonRemovedMuPlus+fEtNonRemovedPiMinus+fEtNonRemovedPiPlus+fEtNonRemovedPositrons
1726// +fEtNonRemovedProtons,
1727// fCentClass);
1728// fHistEtRemovedNeutrals->Fill(fEtNonRemovedNeutrons+fEtNonRemovedAntiNeutrons, fCentClass);
1729
1730 fHistEtRemovedCharged->Fill(fEnergyChargedRemoved, fNClusters);
1731 fHistEtRemovedNeutrals->Fill(fEnergyNeutralRemoved, fNClusters);
1732 fHistEtNonRemovedCharged->Fill(fEnergyChargedNotRemoved, fNClusters);
1733 fHistEtNonRemovedNeutrals->Fill(fEnergyNeutralNotRemoved, fNClusters);
f61cec2f 1734
ef647350 1735 fHistMultRemovedCharged->Fill(fChargedRemoved, fNClusters);
1736 fHistMultRemovedNeutrals->Fill(fNeutralRemoved, fNClusters);
1737 fHistMultNonRemovedCharged->Fill(fChargedNotRemoved, fNClusters);
1738 fHistMultNonRemovedNeutrals->Fill(fNeutralNotRemoved, fNClusters);
f61cec2f 1739
1740
e2ee5727 1741 fHistMultNonRemovedProtons->Fill(fMultNonRemovedProtons, fCentClass);
1742 fHistMultNonRemovedAntiProtons->Fill(fMultNonRemovedAntiProtons, fCentClass);
1743 fHistMultNonRemovedKaonPlus->Fill(fMultNonRemovedKaonPlus, fCentClass);
1744 fHistMultNonRemovedKaonMinus->Fill(fMultNonRemovedKaonMinus, fCentClass);
1745 fHistMultNonRemovedK0s->Fill(fMultNonRemovedK0s, fCentClass);
ef647350 1746 fHistMultNonRemovedK0L->Fill(fMultNonRemovedK0L, fCentClass);
e2ee5727 1747 fHistMultNonRemovedLambdas->Fill(fMultNonRemovedLambdas, fCentClass);
1748 fHistMultNonRemovedPiPlus->Fill(fMultNonRemovedPiPlus, fCentClass);
1749 fHistMultNonRemovedPiMinus->Fill(fMultNonRemovedPiMinus, fCentClass);
1750 fHistMultNonRemovedElectrons->Fill(fMultNonRemovedElectrons, fCentClass);
1751 fHistMultNonRemovedPositrons->Fill(fMultNonRemovedPositrons, fCentClass);
1752 fHistMultNonRemovedMuPlus->Fill(fMultNonRemovedMuPlus, fCentClass);
1753 fHistMultNonRemovedMuMinus->Fill(fMultNonRemovedMuMinus, fCentClass);
1754 fHistMultNonRemovedNeutrons->Fill(fMultNonRemovedNeutrons, fCentClass);
1755 fHistMultNonRemovedAntiNeutrons->Fill(fMultNonRemovedAntiNeutrons, fCentClass);
1756 fHistMultNonRemovedGammas->Fill(fMultNonRemovedGammas, fCentClass);
1757
1758 fHistMultRemovedGammas->Fill(fMultRemovedGammas, fCentClass);
1759 fHistMultRemovedNeutrons->Fill(fMultRemovedNeutrons, fCentClass);
1760 fHistMultRemovedAntiNeutrons->Fill(fMultRemovedAntiNeutrons, fCentClass);
1761
1762 fHistTrackMultvsNonRemovedCharged->Fill(fTrackMultInAcc,
1763 fMultNonRemovedAntiProtons+fMultNonRemovedElectrons+fMultNonRemovedKaonMinus+fMultNonRemovedKaonPlus
1764 +fMultNonRemovedMuMinus+fMultNonRemovedMuPlus+fMultNonRemovedPiMinus+fMultNonRemovedPiPlus+fMultNonRemovedPositrons
1765 +fMultNonRemovedProtons);
1766
1767 fHistTrackMultvsNonRemovedNeutral->Fill(fTrackMultInAcc,
b2c10007 1768 fMultNonRemovedNeutrons+fMultNonRemovedAntiNeutrons+fMultNonRemovedK0s+fMultNonRemovedK0L+fMultNonRemovedLambdas+fK0sMult);
e2ee5727 1769
1770 fHistTrackMultvsRemovedGamma->Fill(fTrackMultInAcc,
1771 fMultRemovedGammas);
1772
ef647350 1773 fHistClusterMultvsNonRemovedCharged->Fill(fNClusters,
e2ee5727 1774 fMultNonRemovedAntiProtons+fMultNonRemovedElectrons+fMultNonRemovedKaonMinus
1775 +fMultNonRemovedKaonPlus+fMultNonRemovedMuMinus+fMultNonRemovedMuPlus
1776 +fMultNonRemovedPiMinus+fMultNonRemovedPiPlus+fMultNonRemovedPositrons+fMultNonRemovedProtons);
1777
ef647350 1778 fHistClusterMultvsNonRemovedNeutral->Fill(fNClusters,
b2c10007 1779 fMultNonRemovedNeutrons+fMultNonRemovedAntiNeutrons+fMultNonRemovedK0s+fMultNonRemovedK0L+fMultNonRemovedLambdas+fK0sMult);
e2ee5727 1780
ef647350 1781 fHistClusterMultvsRemovedGamma->Fill(fNClusters,
e2ee5727 1782 fMultRemovedGammas);
1783
13b0d3c1 1784}
1785
e2ee5727 1786
1787
1788
ef647350 1789Int_t AliAnalysisEtMonteCarlo::PrintFamilyTree(Int_t partIdx, AliStack* stack)
4d376d01 1790{ // print family tree
f61cec2f 1791 TParticle *part = stack->Particle(partIdx);
b2c10007 1792// if(part->GetPdgCode() == fgK0SCode)
f61cec2f 1793 {
1794 std::cout << "This is index: " << partIdx << " (" << stack->Particle(partIdx)->GetName() <<") , is it primary: " << stack->IsPhysicalPrimary(partIdx)<< std::endl;
1795 std::cout << "PID: " << part->GetPdgCode() << "/" << part->GetName() << std::endl;
1796 std::cout << "Energy: " << part->Energy() << std::endl;
1797 std::cout << "Vertex: " << part->Vx() << ", " << part->Vy() << ", " << part->Vz() << std::endl;
1798 }
1799 return PrintMothers(partIdx, stack, 1);
ef647350 1800}
1801
1802Int_t AliAnalysisEtMonteCarlo::PrintMothers(Int_t partIdx, AliStack* stack, Int_t gen)
4d376d01 1803{ // print mothers
f61cec2f 1804 char *tabs = new char[gen+1];
1805 for(Int_t i = 0; i < gen; ++i)
1806 {
1807 //std::cout << i << std::endl;
1808 tabs[i] = '\t';
1809 }
1810 tabs[gen] = '\0';
1811 Int_t mothIdx = stack->Particle(partIdx)->GetMother(0);
1812 if(mothIdx < 0)
1813 {
1814 delete [] tabs;
1815 return 0;
1816 }
1817 TParticle *mother = stack->Particle(mothIdx);
b2c10007 1818// if(mother->GetPdgCode() == fgK0SCode)
f61cec2f 1819 {
b2c10007 1820 //std::cout << tabs << "Mother of index: " << partIdx << " (" << stack->Particle(partIdx)->GetName() <<") is: " << mothIdx << ", is it primary: " << stack->IsPhysicalPrimary(mothIdx)<< std::endl;
1821 std::cout << tabs << "Index: " << mothIdx << std::endl;
1822 std::cout << tabs << "Primary: " << stack->IsPhysicalPrimary(mothIdx) << std::endl;
f61cec2f 1823 std::cout << tabs << "PID: " << mother->GetPdgCode() << "/" << mother->GetName() << std::endl;
1824 std::cout << tabs << "Energy: " << mother->Energy() << std::endl;
b2c10007 1825 if(mother->GetFirstMother() >= 0)
1826 {
1827 std::cout << tabs << "Mother(s): " << stack->Particle(mother->GetFirstMother())->GetPdgCode();
1828 if(mother->GetSecondMother() >= 0) std::cout << ", " << stack->Particle(mother->GetSecondMother())->GetPdgCode();
1829 std::cout << std::endl;
1830 }
f61cec2f 1831 std::cout << tabs << "Vertex: " << mother->Vx() << ", " << mother->Vy() << ", " << mother->Vz() << std::endl;
1832 }
1833 if(mother->GetPdgCode() == fgK0SCode)
1834 {
ef647350 1835// std::cout << "K0S!!!!!!!!!!!!!11111!!!!!" << std::endl;
f61cec2f 1836 }
ef647350 1837// std::cout << "Mother of index: " << partIdx << " (" << stack->Particle(partIdx)->GetName() <<") is: " << mothIdx << std::endl;
1838// std::cout << "PID: " << mother->GetPdgCode() << "/" << mother->GetName() << std::endl;
1839// std::cout << "Energy: " << mother->Energy() << std::endl;
1840// std::cout << "Vertex: " << mother->Vx() << ", " << mother->Vy() << ", " << mother->Vz() << std::endl;
f61cec2f 1841
1842 delete [] tabs;
1843 return PrintMothers(mothIdx, stack, gen+1) + 1;
ef647350 1844}
1845
1846Int_t AliAnalysisEtMonteCarlo::GetPrimMother(Int_t partIdx, AliStack *stack)
4d376d01 1847{ // get primary mother
f61cec2f 1848 if(partIdx >= 0)
ef647350 1849 {
b2c10007 1850 //return stack->GetPrimary(partIdx);
1851
f61cec2f 1852 Int_t mothIdx = stack->Particle(partIdx)->GetMother(0);
1853 if(mothIdx < 0) return -1;
1854 TParticle *mother = stack->Particle(mothIdx);
1855 if(mother)
1856 {
f61cec2f 1857 if(stack->IsPhysicalPrimary(mothIdx)) return mothIdx;
1858 else return GetPrimMother(mothIdx, stack);
1859 }
1860 else
1861 {
1862 return -1;
1863 }
ef647350 1864 }
f61cec2f 1865 return -1;
ef647350 1866}
1867
1868Int_t AliAnalysisEtMonteCarlo::GetK0InFamily(Int_t partIdx, AliStack* stack)
4d376d01 1869{ // get K0 in family
f61cec2f 1870 if(partIdx >= 0)
ef647350 1871 {
2aab9269 1872 if(stack->Particle(partIdx)->GetPdgCode() == fgK0SCode) return partIdx;
f61cec2f 1873 Int_t mothIdx = stack->Particle(partIdx)->GetMother(0);
1874 if(mothIdx < 0) return -1;
1875 TParticle *mother = stack->Particle(mothIdx);
1876 if(mother)
ef647350 1877 {
2aab9269 1878 if(mother->GetPdgCode() == fgK0SCode)
ef647350 1879 {
f61cec2f 1880 return mothIdx;
ef647350 1881 }
f61cec2f 1882 return GetK0InFamily(mothIdx, stack);
1883 }
1884 else
1885 {
1886 return -1;
ef647350 1887 }
1888 }
f61cec2f 1889 return -1;
ef647350 1890}
e2ee5727 1891