Updated gamma jet task
[u/mrichter/AliRoot.git] / PWG4 / GammaConv / AliAnalysisTaskGammaConversion.cxx
CommitLineData
d7d7e825 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *
5 * Version 1.1 *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16////////////////////////////////////////////////
17//---------------------------------------------
18// Class used to do analysis on conversion pairs
19//---------------------------------------------
20////////////////////////////////////////////////
21
22// root
23#include <TChain.h>
24
25// analysis
26#include "AliAnalysisTaskGammaConversion.h"
27#include "AliStack.h"
28#include "AliLog.h"
29#include "AliESDtrackCuts.h"
30#include "TNtuple.h"
4a6157dc 31//#include "AliCFManager.h" // for CF
32//#include "AliCFContainer.h" // for CF
9c1cb6f7 33#include "AliESDInputHandler.h"
34#include "AliAnalysisManager.h"
4a6157dc 35#include "AliGammaConversionAODObject.h"
5e55d806 36#include "AliGammaConversionBGHandler.h"
6272370b 37#include "AliESDCaloCluster.h" // for combining PHOS and GammaConv
9c1cb6f7 38#include "AliKFVertex.h"
d707e3cf 39class AliESDTrackCuts;
4a6157dc 40class AliCFContainer;
41class AliCFManager;
d7d7e825 42class AliKFVertex;
43class AliAODHandler;
44class AliAODEvent;
45class ALiESDEvent;
46class AliMCEvent;
47class AliMCEventHandler;
48class AliESDInputHandler;
49class AliAnalysisManager;
50class Riostream;
51class TFile;
52class TInterpreter;
53class TSystem;
54class TROOT;
55
56ClassImp(AliAnalysisTaskGammaConversion)
57
58
59AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion():
60AliAnalysisTaskSE(),
61 fV0Reader(NULL),
62 fStack(NULL),
a0b94e5c 63 fMCTruth(NULL), // for CF
4a6157dc 64 fGCMCEvent(NULL), // for CF
d7d7e825 65 fESDEvent(NULL),
66 fOutputContainer(NULL),
a0b94e5c 67 fCFManager(0x0), // for CF
d7d7e825 68 fHistograms(NULL),
b5832f95 69 fTriggerCINT1B(kFALSE),
d7d7e825 70 fDoMCTruth(kFALSE),
71 fDoNeutralMeson(kFALSE),
6272370b 72 fDoOmegaMeson(kFALSE),
d7d7e825 73 fDoJet(kFALSE),
74 fDoChic(kFALSE),
9c1cb6f7 75 fRecalculateV0ForGamma(kFALSE),
6c84d371 76 fKFReconstructedGammasTClone(NULL),
6272370b 77 fKFReconstructedPi0sTClone(NULL),
9c1cb6f7 78 fKFRecalculatedGammasTClone(NULL),
6c84d371 79 fCurrentEventPosElectronTClone(NULL),
80 fCurrentEventNegElectronTClone(NULL),
81 fKFReconstructedGammasCutTClone(NULL),
82 fPreviousEventTLVNegElectronTClone(NULL),
83 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 84 fElectronv1(),
85 fElectronv2(),
6272370b 86 fGammav1(),
87 fGammav2(),
9c1cb6f7 88 fElectronRecalculatedv1(),
89 fElectronRecalculatedv2(),
d7d7e825 90 fElectronMass(-1),
91 fGammaMass(-1),
92 fPi0Mass(-1),
93 fEtaMass(-1),
94 fGammaWidth(-1),
95 fPi0Width(-1),
96 fEtaWidth(-1),
97 fMinOpeningAngleGhostCut(0.),
9640a3d1 98 fEsdTrackCuts(NULL),
d7d7e825 99 fCalculateBackground(kFALSE),
100 fWriteNtuple(kFALSE),
101 fGammaNtuple(NULL),
102 fNeutralMesonNtuple(NULL),
103 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 104 fChargedParticles(NULL),
d7d7e825 105 fChargedParticlesId(),
106 fGammaPtHighest(0.),
107 fMinPtForGammaJet(1.),
108 fMinIsoConeSize(0.2),
109 fMinPtIsoCone(0.7),
110 fMinPtGamChargedCorr(0.5),
111 fMinPtJetCone(0.5),
112 fLeadingChargedIndex(-1),
9640a3d1 113 fLowPtMapping(1.),
114 fHighPtMapping(3.),
1e7846f4 115 fDoCF(kFALSE),
d7d7e825 116 fAODBranch(NULL),
037dc2db 117 fAODBranchName("GammaConv"),
d765d400 118 fKFForceAOD(kFALSE),
332f1f44 119 fKFDeltaAODFileName(""),
037dc2db 120 fDoNeutralMesonV0MCCheck(kFALSE),
121 fKFReconstructedGammasV0Index()
d7d7e825 122{
123 // Default constructor
5a34881d 124
125 /* Kenneth: the default constructor should not have any define input/output or the call to SetESDtrackCuts
d7d7e825 126 // Common I/O in slot 0
127 DefineInput (0, TChain::Class());
128 DefineOutput(0, TTree::Class());
129
130 // Your private output
131 DefineOutput(1, TList::Class());
a0b94e5c 132
d7d7e825 133 // Define standard ESD track cuts for Gamma-hadron correlation
134 SetESDtrackCuts();
5a34881d 135 */
d7d7e825 136}
137
138AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name):
139 AliAnalysisTaskSE(name),
140 fV0Reader(NULL),
141 fStack(NULL),
a0b94e5c 142 fMCTruth(NULL), // for CF
4a6157dc 143 fGCMCEvent(NULL), // for CF
d7d7e825 144 fESDEvent(NULL),
145 fOutputContainer(0x0),
a0b94e5c 146 fCFManager(0x0), // for CF
d7d7e825 147 fHistograms(NULL),
b5832f95 148 fTriggerCINT1B(kFALSE),
d7d7e825 149 fDoMCTruth(kFALSE),
150 fDoNeutralMeson(kFALSE),
6272370b 151 fDoOmegaMeson(kFALSE),
d7d7e825 152 fDoJet(kFALSE),
153 fDoChic(kFALSE),
9c1cb6f7 154 fRecalculateV0ForGamma(kFALSE),
6c84d371 155 fKFReconstructedGammasTClone(NULL),
6272370b 156 fKFReconstructedPi0sTClone(NULL),
9c1cb6f7 157 fKFRecalculatedGammasTClone(NULL),
6c84d371 158 fCurrentEventPosElectronTClone(NULL),
159 fCurrentEventNegElectronTClone(NULL),
160 fKFReconstructedGammasCutTClone(NULL),
161 fPreviousEventTLVNegElectronTClone(NULL),
162 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 163 fElectronv1(),
164 fElectronv2(),
6272370b 165 fGammav1(),
166 fGammav2(),
9c1cb6f7 167 fElectronRecalculatedv1(),
168 fElectronRecalculatedv2(),
d7d7e825 169 fElectronMass(-1),
170 fGammaMass(-1),
171 fPi0Mass(-1),
172 fEtaMass(-1),
173 fGammaWidth(-1),
174 fPi0Width(-1),
175 fEtaWidth(-1),
176 fMinOpeningAngleGhostCut(0.),
9640a3d1 177 fEsdTrackCuts(NULL),
d7d7e825 178 fCalculateBackground(kFALSE),
179 fWriteNtuple(kFALSE),
180 fGammaNtuple(NULL),
181 fNeutralMesonNtuple(NULL),
182 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 183 fChargedParticles(NULL),
d7d7e825 184 fChargedParticlesId(),
185 fGammaPtHighest(0.),
186 fMinPtForGammaJet(1.),
187 fMinIsoConeSize(0.2),
188 fMinPtIsoCone(0.7),
189 fMinPtGamChargedCorr(0.5),
190 fMinPtJetCone(0.5),
191 fLeadingChargedIndex(-1),
9640a3d1 192 fLowPtMapping(1.),
193 fHighPtMapping(3.),
1e7846f4 194 fDoCF(kFALSE),
d7d7e825 195 fAODBranch(NULL),
d765d400 196 fAODBranchName(name),
197 fKFForceAOD(kFALSE),
332f1f44 198 fKFDeltaAODFileName(""),
037dc2db 199 fDoNeutralMesonV0MCCheck(kFALSE),
200 fKFReconstructedGammasV0Index()
d7d7e825 201{
202 // Common I/O in slot 0
203 DefineInput (0, TChain::Class());
204 DefineOutput(0, TTree::Class());
205
206 // Your private output
207 DefineOutput(1, TList::Class());
a0b94e5c 208 DefineOutput(2, AliCFContainer::Class()); // for CF
209
210
d7d7e825 211 // Define standard ESD track cuts for Gamma-hadron correlation
212 SetESDtrackCuts();
9c1cb6f7 213
d7d7e825 214}
215
216AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion()
217{
218 // Remove all pointers
219
220 if(fOutputContainer){
221 fOutputContainer->Clear() ;
222 delete fOutputContainer ;
223 }
224 if(fHistograms){
225 delete fHistograms;
226 }
227 if(fV0Reader){
228 delete fV0Reader;
229 }
a0b94e5c 230
231 // for CF
232 if(fCFManager){
233 delete fCFManager;
234 }
9640a3d1 235
236 if(fEsdTrackCuts){
237 delete fEsdTrackCuts;
238 }
239
d7d7e825 240 if (fAODBranch) {
241 fAODBranch->Clear();
242 delete fAODBranch ;
243 }
244}
245
246
247void AliAnalysisTaskGammaConversion::Init()
248{
249 // Initialization
4a6157dc 250 // AliLog::SetGlobalLogLevel(AliLog::kError);
d7d7e825 251}
252void AliAnalysisTaskGammaConversion::SetESDtrackCuts()
253{
254 // SetESDtrackCuts
9640a3d1 255 if (fEsdTrackCuts!=NULL){
256 delete fEsdTrackCuts;
257 }
d7d7e825 258 fEsdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
a0b94e5c 259 //standard cuts from:
260 //http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup
ebcfaa7e 261
262 // Cuts used up to 3rd of March
263
264 // fEsdTrackCuts->SetMinNClustersTPC(50);
265// fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
266// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
267// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
268// fEsdTrackCuts->SetMaxNsigmaToVertex(3);
269// fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE);
270
271 //------- To be tested-----------
d707e3cf 272 // Cuts used up to 26th of Agost
273// Int_t minNClustersTPC = 70;
274// Double_t maxChi2PerClusterTPC = 4.0;
275// Double_t maxDCAtoVertexXY = 2.4; // cm
276// Double_t maxDCAtoVertexZ = 3.2; // cm
277// fEsdTrackCuts->SetRequireSigmaToVertex(kFALSE);
278// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
279// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
280// // fEsdTrackCuts->SetRequireTPCStandAlone(kTRUE);
281// fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
282// fEsdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
283// fEsdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
284// fEsdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
285// fEsdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
286// fEsdTrackCuts->SetDCAToVertex2D(kTRUE);
287// fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
288// fEsdTrackCuts->SetPtRange(0.15);
289
290// fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
291
292
293// Using standard function for setting Cuts
294 Bool_t selectPrimaries=kTRUE;
295 fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selectPrimaries);
296 fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
297 fEsdTrackCuts->SetPtRange(0.15);
ebcfaa7e 298
299 //----- From Jacek 10.03.03 ------------------/
300// minNClustersTPC = 70;
301// maxChi2PerClusterTPC = 4.0;
302// maxDCAtoVertexXY = 2.4; // cm
303// maxDCAtoVertexZ = 3.2; // cm
304
305// esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
306// esdTrackCuts->SetRequireTPCRefit(kFALSE);
307// esdTrackCuts->SetRequireTPCStandAlone(kTRUE);
308// esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
309// esdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
310// esdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
311// esdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
312// esdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
313// esdTrackCuts->SetDCAToVertex2D(kTRUE);
314
315
316
d7d7e825 317 // fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
037dc2db 318 // fV0Reader->SetESDtrackCuts(fEsdTrackCuts);
d7d7e825 319}
320
c00009fb 321void AliAnalysisTaskGammaConversion::UserExec(Option_t */*option*/)
d7d7e825 322{
323 // Execute analysis for current event
6272370b 324
9c1cb6f7 325 // Load the esdpid from the esdhandler if exists (tender was applied) otherwise set the Bethe Bloch parameters
326
327 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
328 AliESDInputHandler *esdHandler=0x0;
329 if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){
330 AliV0Reader::SetESDpid(esdHandler->GetESDpid());
331 } else {
332 //load esd pid bethe bloch parameters depending on the existance of the MC handler
333 // yes: MC parameters
334 // no: data parameters
335 if (!AliV0Reader::GetESDpid()){
336 if (fMCEvent ) {
337 AliV0Reader::InitESDpid();
338 } else {
339 AliV0Reader::InitESDpid(1);
340 }
341 }
342 }
343
d765d400 344 //Must set fForceAOD to true for the AOD to get filled. Should only be done when running independent chain / train.
345 if(fKFForceAOD) {
346 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
347 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
348 }
48682642 349
350 if(fV0Reader == NULL){
351 // Write warning here cuts and so on are default if this ever happens
6272370b 352 }
48682642 353
354 fV0Reader->SetInputAndMCEvent(InputEvent(), MCEvent());
355
356 fV0Reader->Initialize();
357 fDoMCTruth = fV0Reader->GetDoMCTruth();
358
359
360 // if(fDoMCTruth){
361 //ConnectInputData("");
362 //}
d7d7e825 363 //Each event needs an empty branch
1e7846f4 364 // fAODBranch->Clear();
365 fAODBranch->Delete();
a0b94e5c 366
6c84d371 367 if(fKFReconstructedGammasTClone == NULL){
368 fKFReconstructedGammasTClone = new TClonesArray("AliKFParticle",0);
369 }
370 if(fCurrentEventPosElectronTClone == NULL){
371 fCurrentEventPosElectronTClone = new TClonesArray("AliESDtrack",0);
372 }
373 if(fCurrentEventNegElectronTClone == NULL){
374 fCurrentEventNegElectronTClone = new TClonesArray("AliESDtrack",0);
375 }
376 if(fKFReconstructedGammasCutTClone == NULL){
377 fKFReconstructedGammasCutTClone = new TClonesArray("AliKFParticle",0);
378 }
379 if(fPreviousEventTLVNegElectronTClone == NULL){
380 fPreviousEventTLVNegElectronTClone = new TClonesArray("TLorentzVector",0);
381 }
382 if(fPreviousEventTLVPosElectronTClone == NULL){
383 fPreviousEventTLVPosElectronTClone = new TClonesArray("TLorentzVector",0);
384 }
385 if(fChargedParticles == NULL){
386 fChargedParticles = new TClonesArray("AliESDtrack",0);
387 }
6272370b 388
389 if(fKFReconstructedPi0sTClone == NULL){
390 fKFReconstructedPi0sTClone = new TClonesArray("AliKFParticle",0);
391 }
9c1cb6f7 392
393 if(fKFRecalculatedGammasTClone == NULL){
394 fKFRecalculatedGammasTClone = new TClonesArray("AliKFParticle",0);
395 }
396
a0b94e5c 397
6c84d371 398 //clear TClones
1e7846f4 399 fKFReconstructedGammasTClone->Delete();
400 fCurrentEventPosElectronTClone->Delete();
401 fCurrentEventNegElectronTClone->Delete();
402 fKFReconstructedGammasCutTClone->Delete();
403 fPreviousEventTLVNegElectronTClone->Delete();
404 fPreviousEventTLVPosElectronTClone->Delete();
6272370b 405 fKFReconstructedPi0sTClone->Delete();
9c1cb6f7 406 fKFRecalculatedGammasTClone->Delete();
5e55d806 407
d7d7e825 408 //clear vectors
d7d7e825 409 fElectronv1.clear();
410 fElectronv2.clear();
6272370b 411
412 fGammav1.clear();
413 fGammav2.clear();
414
9c1cb6f7 415 fElectronRecalculatedv1.clear();
416 fElectronRecalculatedv2.clear();
417
d7d7e825 418
1e7846f4 419 fChargedParticles->Delete();
5e55d806 420
d7d7e825 421 fChargedParticlesId.clear();
037dc2db 422
423 fKFReconstructedGammasV0Index.clear();
a0b94e5c 424
d7d7e825 425 //Clear the data in the v0Reader
5e55d806 426 // fV0Reader->UpdateEventByEventData();
b5832f95 427
428 //Take Only events with proper trigger
5e55d806 429 /*
b5832f95 430 if(fTriggerCINT1B){
431 if(!fV0Reader->GetESDEvent()->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL")) return;
432 }
5e55d806 433 */
10e3319b 434
c8206114 435 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
436 // cout<< "Event not taken"<< endl;
437 return; // aborts if the primary vertex does not have contributors.
438 }
439
440
885114d1 441 if(!fV0Reader->CheckForPrimaryVertexZ() ){
442 return;
443 }
444
445
446
d7d7e825 447 // Process the MC information
448 if(fDoMCTruth){
449 ProcessMCData();
450 }
a0b94e5c 451
d7d7e825 452 //Process the v0 information with no cuts
453 ProcessV0sNoCut();
cb90a330 454
d7d7e825 455 // Process the v0 information
456 ProcessV0s();
cb90a330 457
5e55d806 458
d7d7e825 459 //Fill Gamma AOD
460 FillAODWithConversionGammas() ;
a0b94e5c 461
d7d7e825 462 // Process reconstructed gammas
463 if(fDoNeutralMeson == kTRUE){
a0b94e5c 464 ProcessGammasForNeutralMesonAnalysis();
77880bd8 465
48682642 466 }
467
77880bd8 468 if(fDoMCTruth == kTRUE){
469 CheckV0Efficiency();
470 }
d7d7e825 471 //Process reconstructed gammas electrons for Chi_c Analysis
6c84d371 472 if(fDoChic == kTRUE){
d7d7e825 473 ProcessGammaElectronsForChicAnalysis();
474 }
475 // Process reconstructed gammas for gamma Jet/hadron correlations
476 if(fDoJet == kTRUE){
477 ProcessGammasForGammaJetAnalysis();
478 }
a0b94e5c 479
5e55d806 480 //calculate background if flag is set
481 if(fCalculateBackground){
482 CalculateBackground();
483 }
48682642 484
485 if(fDoNeutralMeson == kTRUE){
486// ProcessConvPHOSGammasForNeutralMesonAnalysis();
487 if(fDoOmegaMeson == kTRUE){
488 ProcessGammasForOmegaMesonAnalysis();
489 }
490 }
491
5e55d806 492 //Clear the data in the v0Reader
493 fV0Reader->UpdateEventByEventData();
9c1cb6f7 494 if(fRecalculateV0ForGamma==kTRUE){
495 RecalculateV0ForGamma();
496 }
d7d7e825 497 PostData(1, fOutputContainer);
a0b94e5c 498 PostData(2, fCFManager->GetParticleContainer()); // for CF
d7d7e825 499
500}
501
48682642 502// void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t *option){
503// // see header file for documentation
504// // printf(" ConnectInputData %s\n", GetName());
8a685cf3 505
48682642 506// AliAnalysisTaskSE::ConnectInputData(option);
8a685cf3 507
48682642 508// if(fV0Reader == NULL){
509// // Write warning here cuts and so on are default if this ever happens
510// }
511// fV0Reader->Initialize();
512// fDoMCTruth = fV0Reader->GetDoMCTruth();
513// }
d7d7e825 514
515
516
517void AliAnalysisTaskGammaConversion::ProcessMCData(){
518 // see header file for documentation
48682642 519 //InputEvent(), MCEvent());
520 /* TestAnaMarin
d7d7e825 521 fStack = fV0Reader->GetMCStack();
a0b94e5c 522 fMCTruth = fV0Reader->GetMCTruth(); // for CF
4a6157dc 523 fGCMCEvent = fV0Reader->GetMCEvent(); // for CF
48682642 524 */
525 fStack= MCEvent()->Stack();
526 fGCMCEvent=MCEvent();
a0b94e5c 527
528 // for CF
1e7846f4 529 Double_t containerInput[3];
530 if(fDoCF){
531 if(!fGCMCEvent) cout << "NO MC INFO FOUND" << endl;
532 fCFManager->SetEventInfo(fGCMCEvent);
533 }
a0b94e5c 534 // end for CF
535
536
d7d7e825 537 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
538 return; // aborts if the primary vertex does not have contributors.
539 }
a0b94e5c 540
d7d7e825 541 for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
542 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
543
544 if (!particle) {
545 //print warning here
546 continue;
547 }
a0b94e5c 548
d7d7e825 549 ///////////////////////Begin Chic Analysis/////////////////////////////
a0b94e5c 550
a68437fb 551 if(particle->GetPdgCode() == 443){//Is JPsi
d7d7e825 552 if(particle->GetNDaughters()==2){
553 if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 &&
554 TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){
63e16c52 555
d7d7e825 556 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
557 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
558 if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9)
559 fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance
a0b94e5c 560
d7d7e825 561 if( TMath::Abs(daug0->Eta()) < 0.9){
562 if(daug0->GetPdgCode() == -11)
563 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
564 else
565 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
a0b94e5c 566
d7d7e825 567 }
568 if(TMath::Abs(daug1->Eta()) < 0.9){
569 if(daug1->GetPdgCode() == -11)
570 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
571 else
572 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
573 }
574 }
575 }
576 }
577 // const int CHI_C0 = 10441;
578 // const int CHI_C1 = 20443;
579 // const int CHI_C2 = 445
580 if(particle->GetPdgCode() == 22){//gamma from JPsi
581 if(particle->GetMother(0) > -1){
582 if(fStack->Particle(particle->GetMother(0))->GetPdgCode() == 10441 ||
583 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 20443 ||
584 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 445){
585 if(TMath::Abs(particle->Eta()) < 1.2)
586 fHistograms->FillTable("Table_Electrons",17);// gamma from chic inside accptance
587 }
588 }
589 }
590 if(particle->GetPdgCode() == 10441 || particle->GetPdgCode() == 20443 || particle->GetPdgCode() == 445){
591 if( particle->GetNDaughters() == 2){
592 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
593 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
a0b94e5c 594
d7d7e825 595 if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){
596 if( daug0->GetPdgCode() == 443){
597 TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter());
598 TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter());
599 if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
600 fHistograms->FillTable("Table_Electrons",18);
a0b94e5c 601
d7d7e825 602 }//if
603 else if (daug1->GetPdgCode() == 443){
604 TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter());
605 TParticle* daugE1 = fStack->Particle(daug1->GetLastDaughter());
606 if( TMath::Abs(daug0->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
607 fHistograms->FillTable("Table_Electrons",18);
608 }//else if
609 }//gamma o Jpsi
610 }//GetNDaughters
611 }
a0b94e5c 612
613
d7d7e825 614 /////////////////////End Chic Analysis////////////////////////////
a0b94e5c 615
616
d7d7e825 617 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
a0b94e5c 618
d7d7e825 619 if(particle->R()>fV0Reader->GetMaxRCut()) continue; // cuts on distance from collision point
620
621 Double_t tmpPhi=particle->Phi();
622
623 if(particle->Phi()> TMath::Pi()){
624 tmpPhi = particle->Phi()-(2*TMath::Pi());
625 }
626
627 Double_t rapidity;
628 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
629 rapidity=0;
630 }
631 else{
632 rapidity = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())));
633 }
634
635 //process the gammas
636 if (particle->GetPdgCode() == 22){
a607c218 637
638
d7d7e825 639 if(particle->GetMother(0) >-1 && fStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
640 continue; // no photon as mothers!
641 }
a607c218 642
d7d7e825 643 if(particle->GetMother(0) >= fStack->GetNprimary()){
644 continue; // the gamma has a mother, and it is not a primary particle
645 }
a0b94e5c 646
d7d7e825 647 fHistograms->FillHistogram("MC_allGamma_Energy", particle->Energy());
648 fHistograms->FillHistogram("MC_allGamma_Pt", particle->Pt());
649 fHistograms->FillHistogram("MC_allGamma_Eta", particle->Eta());
650 fHistograms->FillHistogram("MC_allGamma_Phi", tmpPhi);
651 fHistograms->FillHistogram("MC_allGamma_Rapid", rapidity);
652
a0b94e5c 653 // for CF
1e7846f4 654 if(fDoCF){
655 containerInput[0] = particle->Pt();
656 containerInput[1] = particle->Eta();
657 if(particle->GetMother(0) >=0){
658 containerInput[2] = fStack->Particle(particle->GetMother(0))->GetMass();
659 }
660 else{
661 containerInput[2]=-1;
662 }
a607c218 663
1e7846f4 664 fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated); // generated gamma
665 }
d7d7e825 666 if(particle->GetMother(0) < 0){ // direct gamma
667 fHistograms->FillHistogram("MC_allDirectGamma_Energy",particle->Energy());
668 fHistograms->FillHistogram("MC_allDirectGamma_Pt", particle->Pt());
669 fHistograms->FillHistogram("MC_allDirectGamma_Eta", particle->Eta());
670 fHistograms->FillHistogram("MC_allDirectGamma_Phi", tmpPhi);
671 fHistograms->FillHistogram("MC_allDirectGamma_Rapid", rapidity);
672 }
673
d7d7e825 674 // looking for conversion (electron + positron from pairbuilding (= 5) )
675 TParticle* ePos = NULL;
676 TParticle* eNeg = NULL;
677
678 if(particle->GetNDaughters() >= 2){
679 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
680 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
681 if(tmpDaughter->GetUniqueID() == 5){
682 if(tmpDaughter->GetPdgCode() == 11){
683 eNeg = tmpDaughter;
684 }
685 else if(tmpDaughter->GetPdgCode() == -11){
686 ePos = tmpDaughter;
687 }
688 }
689 }
690 }
691
692
693 if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
694 continue;
695 }
696
697
698 Double_t ePosPhi = ePos->Phi();
699 if(ePos->Phi()> TMath::Pi()) ePosPhi = ePos->Phi()-(2*TMath::Pi());
700
701 Double_t eNegPhi = eNeg->Phi();
702 if(eNeg->Phi()> TMath::Pi()) eNegPhi = eNeg->Phi()-(2*TMath::Pi());
703
704
705 if(ePos->Pt()<fV0Reader->GetPtCut() || eNeg->Pt()<fV0Reader->GetPtCut()){
706 continue; // no reconstruction below the Pt cut
707 }
a0b94e5c 708
d7d7e825 709 if(TMath::Abs(ePos->Eta())> fV0Reader->GetEtaCut() || TMath::Abs(eNeg->Eta())> fV0Reader->GetEtaCut()){
710 continue;
711 }
a0b94e5c 712
d7d7e825 713 if(ePos->R()>fV0Reader->GetMaxRCut()){
714 continue; // cuts on distance from collision point
715 }
a0b94e5c 716
717 if(TMath::Abs(ePos->Vz()) > fV0Reader->GetMaxZCut()){
d7d7e825 718 continue; // outside material
719 }
720
a0b94e5c 721
d7d7e825 722 if((TMath::Abs(ePos->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() > ePos->R()){
723 continue; // line cut to exclude regions where we do not reconstruct
724 }
725
a0b94e5c 726
727 // for CF
1e7846f4 728 if(fDoCF){
729 fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructable); // reconstructable gamma
730 }
d7d7e825 731 fHistograms->FillHistogram("MC_ConvGamma_Energy", particle->Energy());
732 fHistograms->FillHistogram("MC_ConvGamma_Pt", particle->Pt());
733 fHistograms->FillHistogram("MC_ConvGamma_Eta", particle->Eta());
734 fHistograms->FillHistogram("MC_ConvGamma_Phi", tmpPhi);
735 fHistograms->FillHistogram("MC_ConvGamma_Rapid", rapidity);
736 fHistograms->FillHistogram("MC_ConvGamma_Pt_Eta", particle->Pt(),particle->Eta());
737
738 fHistograms->FillHistogram("MC_E_Energy", eNeg->Energy());
739 fHistograms->FillHistogram("MC_E_Pt", eNeg->Pt());
740 fHistograms->FillHistogram("MC_E_Eta", eNeg->Eta());
741 fHistograms->FillHistogram("MC_E_Phi", eNegPhi);
742
743 fHistograms->FillHistogram("MC_P_Energy", ePos->Energy());
744 fHistograms->FillHistogram("MC_P_Pt", ePos->Pt());
745 fHistograms->FillHistogram("MC_P_Eta", ePos->Eta());
746 fHistograms->FillHistogram("MC_P_Phi", ePosPhi);
747
748
d7d7e825 749 // begin Mapping
750 Int_t rBin = fHistograms->GetRBin(ePos->R());
9640a3d1 751 Int_t zBin = fHistograms->GetZBin(ePos->Vz());
d7d7e825 752 Int_t phiBin = fHistograms->GetPhiBin(particle->Phi());
9c1cb6f7 753 Double_t rFMD=30;
9640a3d1 754
755 TVector3 vtxPos(ePos->Vx(),ePos->Vy(),ePos->Vz());
756
d7d7e825 757 TString nameMCMappingPhiR="";
e158cbc3 758 nameMCMappingPhiR.Form("MC_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 759 // fHistograms->FillHistogram(nameMCMappingPhiR, ePos->Vz(), particle->Eta());
d7d7e825 760
761 TString nameMCMappingPhi="";
e158cbc3 762 nameMCMappingPhi.Form("MC_Conversion_Mapping_Phi%02d",phiBin);
9640a3d1 763 // fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());
ebcfaa7e 764 //fHistograms->FillHistogram(nameMCMappingPhi, ePos->Vz(), particle->Eta());
d7d7e825 765
766 TString nameMCMappingR="";
e158cbc3 767 nameMCMappingR.Form("MC_Conversion_Mapping_R%02d",rBin);
9640a3d1 768 // fHistograms->FillHistogram(nameMCMappingR, particle->Eta());
ebcfaa7e 769 //fHistograms->FillHistogram(nameMCMappingR,ePos->Vz(), particle->Eta());
d7d7e825 770
771 TString nameMCMappingPhiInR="";
e158cbc3 772 nameMCMappingPhiInR.Form("MC_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 773 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
774 fHistograms->FillHistogram(nameMCMappingPhiInR, vtxPos.Phi());
775
776 TString nameMCMappingZInR="";
777 nameMCMappingZInR.Form("MC_Conversion_Mapping_Z_in_R_%02d",rBin);
778 fHistograms->FillHistogram(nameMCMappingZInR,ePos->Vz() );
779
780
781 TString nameMCMappingPhiInZ="";
782 nameMCMappingPhiInZ.Form("MC_Conversion_Mapping_Phi_in_Z_%02d",zBin);
783 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
784 fHistograms->FillHistogram(nameMCMappingPhiInZ, vtxPos.Phi());
785
9c1cb6f7 786
787 if(ePos->R()<rFMD){
788 TString nameMCMappingFMDPhiInZ="";
789 nameMCMappingFMDPhiInZ.Form("MC_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
790 fHistograms->FillHistogram(nameMCMappingFMDPhiInZ, vtxPos.Phi());
791 }
792
9640a3d1 793 TString nameMCMappingRInZ="";
794 nameMCMappingRInZ.Form("MC_Conversion_Mapping_R_in_Z_%02d",zBin);
795 fHistograms->FillHistogram(nameMCMappingRInZ,ePos->R() );
796
797 if(particle->Pt() > fLowPtMapping && particle->Pt()< fHighPtMapping){
798 TString nameMCMappingMidPtPhiInR="";
799 nameMCMappingMidPtPhiInR.Form("MC_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
800 fHistograms->FillHistogram(nameMCMappingMidPtPhiInR, vtxPos.Phi());
801
802 TString nameMCMappingMidPtZInR="";
803 nameMCMappingMidPtZInR.Form("MC_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
804 fHistograms->FillHistogram(nameMCMappingMidPtZInR,ePos->Vz() );
805
806
807 TString nameMCMappingMidPtPhiInZ="";
808 nameMCMappingMidPtPhiInZ.Form("MC_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
809 fHistograms->FillHistogram(nameMCMappingMidPtPhiInZ, vtxPos.Phi());
9c1cb6f7 810
811
812 if(ePos->R()<rFMD){
813 TString nameMCMappingMidPtFMDPhiInZ="";
814 nameMCMappingMidPtFMDPhiInZ.Form("MC_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
815 fHistograms->FillHistogram(nameMCMappingMidPtFMDPhiInZ, vtxPos.Phi());
816 }
817
9640a3d1 818
819 TString nameMCMappingMidPtRInZ="";
820 nameMCMappingMidPtRInZ.Form("MC_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
821 fHistograms->FillHistogram(nameMCMappingMidPtRInZ,ePos->R() );
822
823 }
824
d7d7e825 825 //end mapping
826
827 fHistograms->FillHistogram("MC_Conversion_R",ePos->R());
828 fHistograms->FillHistogram("MC_Conversion_ZR",ePos->Vz(),ePos->R());
829 fHistograms->FillHistogram("MC_Conversion_XY",ePos->Vx(),ePos->Vy());
830 fHistograms->FillHistogram("MC_Conversion_OpeningAngle",GetMCOpeningAngle(ePos, eNeg));
9640a3d1 831 fHistograms->FillHistogram("MC_ConvGamma_E_AsymmetryP",particle->P(),eNeg->P()/particle->P());
832 fHistograms->FillHistogram("MC_ConvGamma_P_AsymmetryP",particle->P(),ePos->P()/particle->P());
833
834
d7d7e825 835 if(particle->GetMother(0) < 0){ // no mother = direct gamma, still inside converted
836 fHistograms->FillHistogram("MC_ConvDirectGamma_Energy",particle->Energy());
837 fHistograms->FillHistogram("MC_ConvDirectGamma_Pt", particle->Pt());
838 fHistograms->FillHistogram("MC_ConvDirectGamma_Eta", particle->Eta());
839 fHistograms->FillHistogram("MC_ConvDirectGamma_Phi", tmpPhi);
840 fHistograms->FillHistogram("MC_ConvDirectGamma_Rapid", rapidity);
841
842 } // end direct gamma
843 else{ // mother exits
6c84d371 844 /* if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0
a0b94e5c 845 fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S
846 fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2
847 ){
848 fMCGammaChic.push_back(particle);
849 }
6c84d371 850 */
d7d7e825 851 } // end if mother exits
852 } // end if particle is a photon
853
854
855
856 // process motherparticles (2 gammas as daughters)
857 // the motherparticle had already to pass the R and the eta cut, but no line cut.
858 // the line cut is just valid for the conversions!
859
860 if(particle->GetNDaughters() == 2){
861
862 TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());
863 TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());
864
865 if(daughter0->GetPdgCode() != 22 || daughter1->GetPdgCode() != 22) continue; //check for gamma gamma daughters
866
867 // Check the acceptance for both gammas
01b7fdcc 868 Bool_t gammaEtaCut = kTRUE;
869 if(TMath::Abs(daughter0->Eta()) > fV0Reader->GetEtaCut() || TMath::Abs(daughter1->Eta()) > fV0Reader->GetEtaCut() ) gammaEtaCut = kFALSE;
870 Bool_t gammaRCut = kTRUE;
871 if(daughter0->R() > fV0Reader->GetMaxRCut() || daughter1->R() > fV0Reader->GetMaxRCut() ) gammaRCut = kFALSE;
d7d7e825 872
873
874
875 // check for conversions now -> have to pass eta, R and line cut!
876 Bool_t daughter0Electron = kFALSE;
877 Bool_t daughter0Positron = kFALSE;
878 Bool_t daughter1Electron = kFALSE;
879 Bool_t daughter1Positron = kFALSE;
880
881 if(daughter0->GetNDaughters() >= 2){ // first gamma
882 for(Int_t TrackIndex=daughter0->GetFirstDaughter();TrackIndex<=daughter0->GetLastDaughter();TrackIndex++){
883 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
884 if(tmpDaughter->GetUniqueID() == 5){
885 if(tmpDaughter->GetPdgCode() == 11){
886 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
887 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
888 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
889 daughter0Electron = kTRUE;
890 }
891 }
892
893 }
894 }
895 else if(tmpDaughter->GetPdgCode() == -11){
896 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
897 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
898 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
899 daughter0Positron = kTRUE;
900 }
a68437fb 901 }
902 }
d7d7e825 903 }
904 }
905 }
906 }
907
908
909 if(daughter1->GetNDaughters() >= 2){ // second gamma
910 for(Int_t TrackIndex=daughter1->GetFirstDaughter();TrackIndex<=daughter1->GetLastDaughter();TrackIndex++){
911 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
912 if(tmpDaughter->GetUniqueID() == 5){
913 if(tmpDaughter->GetPdgCode() == 11){
914 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
915 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
916 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
917 daughter1Electron = kTRUE;
918 }
919 }
920
921 }
922 }
923 else if(tmpDaughter->GetPdgCode() == -11){
924 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
925 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
926 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
927 daughter1Positron = kTRUE;
928 }
929 }
930
931 }
932
933 }
934 }
935 }
936 }
937
a0b94e5c 938
d7d7e825 939 if(particle->GetPdgCode()==111){ //Pi0
940 if( iTracks >= fStack->GetNprimary()){
941 fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());
942 fHistograms->FillHistogram("MC_Pi0_Secondaries_Rapid", rapidity);
943 fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);
944 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt", particle->Pt());
945 fHistograms->FillHistogram("MC_Pi0_Secondaries_Energy", particle->Energy());
946 fHistograms->FillHistogram("MC_Pi0_Secondaries_R", particle->R());
947 fHistograms->FillHistogram("MC_Pi0_Secondaries_ZR", particle->Vz(),particle->R());
948 fHistograms->FillHistogram("MC_Pi0_Secondaries_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
949 fHistograms->FillHistogram("MC_Pi0_Secondaries_XY", particle->Vx(),particle->Vy());//only fill from one daughter to avoid multiple filling
950
01b7fdcc 951 if(gammaEtaCut && gammaRCut){
d7d7e825 952 //if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
953 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
954 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
955 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
956 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
957 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
958 }
959 }
960 }
961 else{
962 fHistograms->FillHistogram("MC_Pi0_Eta", particle->Eta());
963 fHistograms->FillHistogram("MC_Pi0_Rapid", rapidity);
964 fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);
965 fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());
966 fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());
967 fHistograms->FillHistogram("MC_Pi0_R", particle->R());
968 fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());
969 fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
970 fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
63e16c52 971 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_Fiducial", particle->Pt());
972
01b7fdcc 973 if(gammaEtaCut && gammaRCut){
d7d7e825 974 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
975 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
976 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
63e16c52 977 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_withinAcceptance_Fiducial", particle->Pt());
978
d7d7e825 979 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
980 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
981 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
982 fHistograms->FillHistogram("MC_Pi0_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
10e3319b 983 fHistograms->FillHistogram("MC_Pi0_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
984 fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
985 fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
986
48682642 987 Double_t alfa=0.;
988 if((daughter0->Energy()+daughter1->Energy())!= 0.){
989 alfa= TMath::Abs((daughter0->Energy()-daughter1->Energy())/(daughter0->Energy()+daughter1->Energy()));
990 }
991 fHistograms->FillHistogram("MC_Pi0_alpha",alfa);
992 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_ConvGamma_withinAcceptance_Fiducial", particle->Pt());
993
d7d7e825 994 }
995 }
996 }
997 }
998
999 if(particle->GetPdgCode()==221){ //Eta
1000 fHistograms->FillHistogram("MC_Eta_Eta", particle->Eta());
1001 fHistograms->FillHistogram("MC_Eta_Rapid", rapidity);
1002 fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);
1003 fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());
1004 fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());
1005 fHistograms->FillHistogram("MC_Eta_R", particle->R());
1006 fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());
1007 fHistograms->FillHistogram("MC_Eta_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1008 fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
1009
01b7fdcc 1010 if(gammaEtaCut && gammaRCut){
a0b94e5c 1011 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 1012 fHistograms->FillHistogram("MC_Eta_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1013 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1014 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1015 fHistograms->FillHistogram("MC_Eta_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1016 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1017 fHistograms->FillHistogram("MC_Eta_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
10e3319b 1018 fHistograms->FillHistogram("MC_Eta_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
1019 fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
1020 fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
1021
d7d7e825 1022 }
1023
1024 }
1025
1026 }
1027
1028 // all motherparticles with 2 gammas as daughters
1029 fHistograms->FillHistogram("MC_Mother_R", particle->R());
1030 fHistograms->FillHistogram("MC_Mother_ZR", particle->Vz(),particle->R());
1031 fHistograms->FillHistogram("MC_Mother_XY", particle->Vx(),particle->Vy());
1032 fHistograms->FillHistogram("MC_Mother_Mass", particle->GetCalcMass());
1033 fHistograms->FillHistogram("MC_Mother_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1034 fHistograms->FillHistogram("MC_Mother_Energy", particle->Energy());
1035 fHistograms->FillHistogram("MC_Mother_Pt", particle->Pt());
1036 fHistograms->FillHistogram("MC_Mother_Eta", particle->Eta());
1037 fHistograms->FillHistogram("MC_Mother_Rapid", rapidity);
1038 fHistograms->FillHistogram("MC_Mother_Phi",tmpPhi);
1039 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt",particle->GetMass(),particle->Pt());
a0b94e5c 1040
01b7fdcc 1041 if(gammaEtaCut && gammaRCut){
a0b94e5c 1042 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 1043 fHistograms->FillHistogram("MC_Mother_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1044 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1045 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_withinAcceptance",particle->GetMass(),particle->Pt());
1046 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1047 fHistograms->FillHistogram("MC_Mother_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1048 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1049 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_ConvGamma_withinAcceptance",particle->GetMass(),particle->Pt());
a0b94e5c 1050
d7d7e825 1051 }
1052
1053
1054 } // end passed R and eta cut
a0b94e5c 1055
d7d7e825 1056 } // end if(particle->GetNDaughters() == 2)
1057
1058 }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++)
a0b94e5c 1059
d7d7e825 1060} // end ProcessMCData
1061
1062
1063
1064void AliAnalysisTaskGammaConversion::FillNtuple(){
1065 //Fills the ntuple with the different values
a0b94e5c 1066
d7d7e825 1067 if(fGammaNtuple == NULL){
1068 return;
1069 }
1070 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1071 for(Int_t i=0;i<numberOfV0s;i++){
1072 Float_t values[27] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
1073 AliESDv0 * cV0 = fV0Reader->GetV0(i);
1074 Double_t negPID=0;
1075 Double_t posPID=0;
1076 fV0Reader->GetPIDProbability(negPID,posPID);
1077 values[0]=cV0->GetOnFlyStatus();
1078 values[1]=fV0Reader->CheckForPrimaryVertex();
1079 values[2]=negPID;
1080 values[3]=posPID;
1081 values[4]=fV0Reader->GetX();
1082 values[5]=fV0Reader->GetY();
1083 values[6]=fV0Reader->GetZ();
1084 values[7]=fV0Reader->GetXYRadius();
1085 values[8]=fV0Reader->GetMotherCandidateNDF();
1086 values[9]=fV0Reader->GetMotherCandidateChi2();
1087 values[10]=fV0Reader->GetMotherCandidateEnergy();
1088 values[11]=fV0Reader->GetMotherCandidateEta();
1089 values[12]=fV0Reader->GetMotherCandidatePt();
1090 values[13]=fV0Reader->GetMotherCandidateMass();
1091 values[14]=fV0Reader->GetMotherCandidateWidth();
1092 // values[15]=fV0Reader->GetMotherMCParticle()->Pt(); MOVED TO THE END, HAS TO BE CALLED AFTER HasSameMother NB: still has the same entry in the array
1093 values[16]=fV0Reader->GetOpeningAngle();
1094 values[17]=fV0Reader->GetNegativeTrackEnergy();
1095 values[18]=fV0Reader->GetNegativeTrackPt();
1096 values[19]=fV0Reader->GetNegativeTrackEta();
1097 values[20]=fV0Reader->GetNegativeTrackPhi();
1098 values[21]=fV0Reader->GetPositiveTrackEnergy();
1099 values[22]=fV0Reader->GetPositiveTrackPt();
1100 values[23]=fV0Reader->GetPositiveTrackEta();
1101 values[24]=fV0Reader->GetPositiveTrackPhi();
1102 values[25]=fV0Reader->HasSameMCMother();
1103 if(values[25] != 0){
1104 values[26]=fV0Reader->GetMotherMCParticlePDGCode();
1105 values[15]=fV0Reader->GetMotherMCParticle()->Pt();
1106 }
1107 fTotalNumberOfAddedNtupleEntries++;
1108 fGammaNtuple->Fill(values);
1109 }
1110 fV0Reader->ResetV0IndexNumber();
1111
1112}
1113
1114void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){
1115 // Process all the V0's without applying any cuts to it
a0b94e5c 1116
d7d7e825 1117 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1118 for(Int_t i=0;i<numberOfV0s;i++){
1119 /*AliESDv0 * cV0 = */fV0Reader->GetV0(i);
a0b94e5c 1120
d7d7e825 1121 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
cb90a330 1122 continue;
d7d7e825 1123 }
9640a3d1 1124
77880bd8 1125 // if( !fV0Reader->GetV0(i)->GetOnFlyStatus()){
1126 if( !fV0Reader->CheckV0FinderStatus(i)){
cb90a330 1127 continue;
9640a3d1 1128 }
1129
cb90a330 1130
1131 if( !((fV0Reader->GetNegativeESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ||
1132 !((fV0Reader->GetPositiveESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ){
1133 continue;
9640a3d1 1134 }
cb90a330 1135
ebcfaa7e 1136 if( fV0Reader->GetNegativeESDTrack()->GetSign()== fV0Reader->GetPositiveESDTrack()->GetSign()){
1137 continue;
1138 }
1139
1140 if( fV0Reader->GetNegativeESDTrack()->GetKinkIndex(0) > 0 ||
1141 fV0Reader->GetPositiveESDTrack()->GetKinkIndex(0) > 0) {
1142 continue;
1143 }
9c1cb6f7 1144 if(TMath::Abs(fV0Reader->GetMotherCandidateEta())> fV0Reader->GetEtaCut()){
1145 continue;
1146 }
1147 if(TMath::Abs(fV0Reader->GetPositiveTrackEta())> fV0Reader->GetEtaCut()){
1148 continue;
1149 }
1150 if(TMath::Abs(fV0Reader->GetNegativeTrackEta())> fV0Reader->GetEtaCut()){
1151 continue;
1152 }
1153 if((TMath::Abs(fV0Reader->GetZ())*fV0Reader->GetLineCutZRSlope())-fV0Reader->GetLineCutZValue() > fV0Reader->GetXYRadius() ){ // cuts out regions where we do not reconstruct
1154 continue;
1155 }
d7d7e825 1156 if(fDoMCTruth){
a0b94e5c 1157
d7d7e825 1158 if(fV0Reader->HasSameMCMother() == kFALSE){
1159 continue;
1160 }
a0b94e5c 1161
d7d7e825 1162 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1163 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 1164
d7d7e825 1165 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1166 continue;
1167 }
26923b22 1168 if(negativeMC->GetPdgCode() == positiveMC->GetPdgCode()){
d7d7e825 1169 continue;
1170 }
a68437fb 1171
26923b22 1172 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){ // id 5 is conversion
a68437fb 1173 continue;
1174 }
a0b94e5c 1175
d7d7e825 1176 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a0b94e5c 1177
d7d7e825 1178 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1179 fHistograms->FillHistogram("ESD_NoCutConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1180 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1181 fHistograms->FillHistogram("ESD_NoCutConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1182 fHistograms->FillHistogram("ESD_NoCutConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1183 fHistograms->FillHistogram("ESD_NoCutConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1184 fHistograms->FillHistogram("ESD_NoCutConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1185 fHistograms->FillHistogram("ESD_NoCutConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1186 fHistograms->FillHistogram("ESD_NoCutConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1187 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
a0b94e5c 1188
d7d7e825 1189 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1190 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
1191
1192 fHistograms->FillHistogram("ESD_NoCutConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1193 fHistograms->FillHistogram("ESD_NoCutConversion_R", fV0Reader->GetXYRadius());
1194 fHistograms->FillHistogram("ESD_NoCutConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1195 fHistograms->FillHistogram("ESD_NoCutConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1196 fHistograms->FillHistogram("ESD_NoCutConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1197 fHistograms->FillHistogram("ESD_NoCutConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1198 fHistograms->FillHistogram("ESD_NoCutConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1199 fHistograms->FillHistogram("ESD_NoCutConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1200
1201 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1202 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1203 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1204 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1205
d7d7e825 1206 //store MCTruth properties
1207 fHistograms->FillHistogram("ESD_NoCutConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1208 fHistograms->FillHistogram("ESD_NoCutConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1209 fHistograms->FillHistogram("ESD_NoCutConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
d7d7e825 1210 }
1211 }
1212 }
1213 fV0Reader->ResetV0IndexNumber();
1214}
1215
1216void AliAnalysisTaskGammaConversion::ProcessV0s(){
1217 // see header file for documentation
037dc2db 1218
1219
d7d7e825 1220 if(fWriteNtuple == kTRUE){
1221 FillNtuple();
1222 }
1223
1224 Int_t nSurvivingV0s=0;
1225 while(fV0Reader->NextV0()){
1226 nSurvivingV0s++;
1227
1228
1229 //-------------------------- filling v0 information -------------------------------------
1230 fHistograms->FillHistogram("ESD_Conversion_R", fV0Reader->GetXYRadius());
1231 fHistograms->FillHistogram("ESD_Conversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1232 fHistograms->FillHistogram("ESD_Conversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1233 fHistograms->FillHistogram("ESD_Conversion_OpeningAngle", fV0Reader->GetOpeningAngle());
1234
1235 fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());
1236 fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());
1237 fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());
1238 fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());
ebcfaa7e 1239 fHistograms->FillHistogram("ESD_E_nTPCClusters", fV0Reader->GetNegativeTracknTPCClusters());
1240 fHistograms->FillHistogram("ESD_E_nITSClusters", fV0Reader->GetNegativeTracknITSClusters());
d7d7e825 1241
1242 fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());
1243 fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());
1244 fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());
1245 fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());
ebcfaa7e 1246 fHistograms->FillHistogram("ESD_P_nTPCClusters", fV0Reader->GetPositiveTracknTPCClusters());
1247 fHistograms->FillHistogram("ESD_P_nITSClusters", fV0Reader->GetPositiveTracknITSClusters());
d7d7e825 1248
1249 fHistograms->FillHistogram("ESD_ConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1250 fHistograms->FillHistogram("ESD_ConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1251 fHistograms->FillHistogram("ESD_ConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1252 fHistograms->FillHistogram("ESD_ConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1253 fHistograms->FillHistogram("ESD_ConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1254 fHistograms->FillHistogram("ESD_ConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1255 fHistograms->FillHistogram("ESD_ConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1256 fHistograms->FillHistogram("ESD_ConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1257 fHistograms->FillHistogram("ESD_ConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1258 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1259
1260 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1261 fHistograms->FillHistogram("ESD_ConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
9640a3d1 1262
1263 fHistograms->FillHistogram("ESD_ConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1264 fHistograms->FillHistogram("ESD_ConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1265 fHistograms->FillHistogram("ESD_ConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1266 fHistograms->FillHistogram("ESD_ConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1267
1268 fHistograms->FillHistogram("ESD_ConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1269 fHistograms->FillHistogram("ESD_ConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1270 fHistograms->FillHistogram("ESD_ConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1271 fHistograms->FillHistogram("ESD_ConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1272
9c1cb6f7 1273 Double_t negPID=0;
1274 Double_t posPID=0;
1275 fV0Reader->GetPIDProbability(negPID,posPID);
1276 fHistograms->FillHistogram("ESD_ConvGamma_E_EProbP",fV0Reader->GetNegativeTrackP(),negPID);
1277 fHistograms->FillHistogram("ESD_ConvGamma_P_EProbP",fV0Reader->GetPositiveTrackP(),posPID);
1278
1279 Double_t negPIDmupi=0;
1280 Double_t posPIDmupi=0;
1281 fV0Reader->GetPIDProbabilityMuonPion(negPIDmupi,posPIDmupi);
1282 fHistograms->FillHistogram("ESD_ConvGamma_E_mupiProbP",fV0Reader->GetNegativeTrackP(),negPIDmupi);
1283 fHistograms->FillHistogram("ESD_ConvGamma_P_mupiProbP",fV0Reader->GetPositiveTrackP(),posPIDmupi);
1284
48682642 1285 Double_t armenterosQtAlfa[2];
1286 fV0Reader->GetArmenterosQtAlfa(fV0Reader-> GetNegativeKFParticle(),
1287 fV0Reader-> GetPositiveKFParticle(),
1288 fV0Reader->GetMotherCandidateKFCombination(),
1289 armenterosQtAlfa);
1290
1291 fHistograms->FillHistogram("ESD_ConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
1292
037dc2db 1293
d7d7e825 1294 // begin mapping
1295 Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());
9640a3d1 1296 Int_t zBin = fHistograms->GetZBin(fV0Reader->GetZ());
d7d7e825 1297 Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());
9c1cb6f7 1298 Double_t rFMD=30;
1299
9640a3d1 1300 TVector3 vtxConv(fV0Reader->GetX(),fV0Reader->GetY(), fV0Reader->GetZ());
1301
ebcfaa7e 1302 // Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();
d7d7e825 1303
1304 TString nameESDMappingPhiR="";
e158cbc3 1305 nameESDMappingPhiR.Form("ESD_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 1306 //fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1307
1308 TString nameESDMappingPhi="";
e158cbc3 1309 nameESDMappingPhi.Form("ESD_Conversion_Mapping_Phi%02d",phiBin);
ebcfaa7e 1310 //fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1311
1312 TString nameESDMappingR="";
e158cbc3 1313 nameESDMappingR.Form("ESD_Conversion_Mapping_R%02d",rBin);
ebcfaa7e 1314 //fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1315
1316 TString nameESDMappingPhiInR="";
e158cbc3 1317 nameESDMappingPhiInR.Form("ESD_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 1318 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1319 fHistograms->FillHistogram(nameESDMappingPhiInR, vtxConv.Phi());
1320
1321 TString nameESDMappingZInR="";
1322 nameESDMappingZInR.Form("ESD_Conversion_Mapping_Z_in_R_%02d",rBin);
1323 fHistograms->FillHistogram(nameESDMappingZInR, fV0Reader->GetZ());
1324
1325 TString nameESDMappingPhiInZ="";
1326 nameESDMappingPhiInZ.Form("ESD_Conversion_Mapping_Phi_in_Z_%02d",zBin);
1327 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1328 fHistograms->FillHistogram(nameESDMappingPhiInZ, vtxConv.Phi());
1329
9c1cb6f7 1330 if(fV0Reader->GetXYRadius()<rFMD){
1331 TString nameESDMappingFMDPhiInZ="";
1332 nameESDMappingFMDPhiInZ.Form("ESD_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
1333 fHistograms->FillHistogram(nameESDMappingFMDPhiInZ, vtxConv.Phi());
1334 }
1335
1336
9640a3d1 1337 TString nameESDMappingRInZ="";
1338 nameESDMappingRInZ.Form("ESD_Conversion_Mapping_R_in_Z_%02d",zBin);
1339 fHistograms->FillHistogram(nameESDMappingRInZ, fV0Reader->GetXYRadius());
1340
1341 if(fV0Reader->GetMotherCandidatePt() > fLowPtMapping && fV0Reader->GetMotherCandidatePt()< fHighPtMapping){
1342 TString nameESDMappingMidPtPhiInR="";
1343 nameESDMappingMidPtPhiInR.Form("ESD_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
1344 fHistograms->FillHistogram(nameESDMappingMidPtPhiInR, vtxConv.Phi());
1345
1346 TString nameESDMappingMidPtZInR="";
1347 nameESDMappingMidPtZInR.Form("ESD_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
1348 fHistograms->FillHistogram(nameESDMappingMidPtZInR, fV0Reader->GetZ());
1349
1350 TString nameESDMappingMidPtPhiInZ="";
1351 nameESDMappingMidPtPhiInZ.Form("ESD_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
1352 fHistograms->FillHistogram(nameESDMappingMidPtPhiInZ, vtxConv.Phi());
9c1cb6f7 1353 if(fV0Reader->GetXYRadius()<rFMD){
1354 TString nameESDMappingMidPtFMDPhiInZ="";
1355 nameESDMappingMidPtFMDPhiInZ.Form("ESD_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
1356 fHistograms->FillHistogram(nameESDMappingMidPtFMDPhiInZ, vtxConv.Phi());
1357 }
1358
9640a3d1 1359
1360 TString nameESDMappingMidPtRInZ="";
1361 nameESDMappingMidPtRInZ.Form("ESD_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
1362 fHistograms->FillHistogram(nameESDMappingMidPtRInZ, fV0Reader->GetXYRadius());
1363
1364 }
1365
1366
d7d7e825 1367 // end mapping
1368
6c84d371 1369 new((*fKFReconstructedGammasTClone)[fKFReconstructedGammasTClone->GetEntriesFast()]) AliKFParticle(*fV0Reader->GetMotherCandidateKFCombination());
037dc2db 1370 fKFReconstructedGammasV0Index.push_back(fV0Reader->GetCurrentV0IndexNumber()-1);
6c84d371 1371 // fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());
d7d7e825 1372 fElectronv1.push_back(fV0Reader->GetCurrentV0()->GetPindex());
1373 fElectronv2.push_back(fV0Reader->GetCurrentV0()->GetNindex());
a0b94e5c 1374
d7d7e825 1375 //----------------------------------- checking for "real" conversions (MC match) --------------------------------------
1376 if(fDoMCTruth){
1377
a68437fb 1378 if(fV0Reader->HasSameMCMother() == kFALSE){
1379 continue;
1380 }
d7d7e825 1381 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
a68437fb 1382 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1383
1384 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1385 continue;
1386 }
1387
1388 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
1389 continue;
1390 }
037dc2db 1391 if(negativeMC->GetUniqueID() == 4 && positiveMC->GetUniqueID() ==4){// fill r distribution for Dalitz decays
1392 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 111){ //pi0
1393 fHistograms->FillHistogram("ESD_TrueDalitzContamination_R", fV0Reader->GetXYRadius());
1394 }
1395 }
a68437fb 1396
1397 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){// check if the daughters come from a conversion
1398 continue;
1399 }
a0b94e5c 1400
d7d7e825 1401 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a68437fb 1402
26923b22 1403 if(fDoCF){
1404 Double_t containerInput[3];
1405 containerInput[0] = fV0Reader->GetMotherCandidatePt();
1406 containerInput[1] = fV0Reader->GetMotherCandidateEta();
1407 containerInput[2] = fV0Reader->GetMotherCandidateMass();
1408 fCFManager->GetParticleContainer()->Fill(containerInput,kStepTrueGamma); // for CF
1409 }
1410
d7d7e825 1411 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1412 fHistograms->FillHistogram("ESD_TrueConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1413 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1414 fHistograms->FillHistogram("ESD_TrueConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1415 fHistograms->FillHistogram("ESD_TrueConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1416 fHistograms->FillHistogram("ESD_TrueConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1417 fHistograms->FillHistogram("ESD_TrueConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1418 fHistograms->FillHistogram("ESD_TrueConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1419 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1420 fHistograms->FillHistogram("ESD_TrueConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1421 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters());
1422 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters());
1423 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters(),fV0Reader->GetMotherCandidateMass());
1424 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters(),fV0Reader->GetMotherCandidateMass());
a0b94e5c 1425
d7d7e825 1426 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1427 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
a0b94e5c 1428
d7d7e825 1429
1430 fHistograms->FillHistogram("ESD_TrueConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1431 fHistograms->FillHistogram("ESD_TrueConversion_R", fV0Reader->GetXYRadius());
1432 fHistograms->FillHistogram("ESD_TrueConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1433 fHistograms->FillHistogram("ESD_TrueConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1434
1435 fHistograms->FillHistogram("ESD_TrueConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1436 fHistograms->FillHistogram("ESD_TrueConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1437 fHistograms->FillHistogram("ESD_TrueConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1438 fHistograms->FillHistogram("ESD_TrueConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
96ade8ef 1439 if (fV0Reader->GetMotherCandidateP() != 0) {
1440 fHistograms->FillHistogram("ESD_TrueConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1441 fHistograms->FillHistogram("ESD_TrueConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1442 } else { cout << "Error::fV0Reader->GetNegativeTrackP() == 0 !!!" << endl; }
1443
9640a3d1 1444 fHistograms->FillHistogram("ESD_TrueConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1445 fHistograms->FillHistogram("ESD_TrueConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
70ef88b5 1446 fHistograms->FillHistogram("ESD_TrueConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
1447
9640a3d1 1448
a0b94e5c 1449
d7d7e825 1450 //store MCTruth properties
1451 fHistograms->FillHistogram("ESD_TrueConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1452 fHistograms->FillHistogram("ESD_TrueConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1453 fHistograms->FillHistogram("ESD_TrueConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
48682642 1454
d7d7e825 1455 //resolution
1456 Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();
1457 Double_t esdpt = fV0Reader->GetMotherCandidatePt();
9c1cb6f7 1458 Double_t resdPt = 0.;
4a6157dc 1459 if(mcpt > 0){
9c1cb6f7 1460 resdPt = ((esdpt - mcpt)/mcpt)*100.;
d7d7e825 1461 }
4a6157dc 1462 else if(mcpt < 0){
1463 cout<<"Pt of MC particle is negative, this will cause wrong calculation of resPt"<<endl;
1464 }
d7d7e825 1465
ca6d4600 1466 fHistograms->FillHistogram("Resolution_Gamma_dPt_Pt", mcpt, resdPt);
d7d7e825 1467 fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);
1468 fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);
ca6d4600 1469 fHistograms->FillHistogram("Resolution_Gamma_dPt_Phi", fV0Reader->GetMotherCandidatePhi(), resdPt);
d7d7e825 1470
9c1cb6f7 1471 Double_t resdZ = 0.;
d7d7e825 1472 if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){
9c1cb6f7 1473 resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100.;
d7d7e825 1474 }
9c1cb6f7 1475 Double_t resdZAbs = 0.;
037dc2db 1476 resdZAbs = (fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz());
1477
1478 fHistograms->FillHistogram("Resolution_dZAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdZAbs);
d7d7e825 1479 fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);
1480 fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());
1481 fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());
ca6d4600 1482
1483 // new for dPt_Pt-histograms for Electron and Positron
1484 Double_t mcEpt = fV0Reader->GetNegativeMCParticle()->Pt();
9c1cb6f7 1485 Double_t resEdPt = 0.;
ca6d4600 1486 if (mcEpt > 0){
9c1cb6f7 1487 resEdPt = ((fV0Reader->GetNegativeTrackPt()-mcEpt)/mcEpt)*100.;
ca6d4600 1488 }
1489 UInt_t statusN = fV0Reader->GetNegativeESDTrack()->GetStatus();
1490 // AliESDtrack * negTrk = fV0Reader->GetNegativeESDTrack();
1491 UInt_t kTRDoutN = (statusN & AliESDtrack::kTRDout);
1492
1493 Int_t ITSclsE= fV0Reader->GetNegativeTracknITSClusters();
1494 // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
1495 switch(ITSclsE){
1496 case 0: // 0 ITS clusters
1497 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS0", mcEpt, resEdPt);
1498 break;
1499 case 1: // 1 ITS cluster
1500 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS1", mcEpt, resEdPt);
1501 break;
1502 case 2: // 2 ITS clusters
1503 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS2", mcEpt, resEdPt);
1504 break;
1505 case 3: // 3 ITS clusters
1506 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS3", mcEpt, resEdPt);
1507 break;
1508 case 4: // 4 ITS clusters
1509 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS4", mcEpt, resEdPt);
1510 break;
1511 case 5: // 5 ITS clusters
1512 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS5", mcEpt, resEdPt);
1513 break;
1514 case 6: // 6 ITS clusters
1515 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS6", mcEpt, resEdPt);
1516 break;
1517 }
1518 //Filling histograms with respect to Electron resolution
1519 fHistograms->FillHistogram("Resolution_E_dPt_Pt", mcEpt, resEdPt);
1520 fHistograms->FillHistogram("Resolution_E_dPt_Phi", fV0Reader->GetNegativeTrackPhi(), resEdPt);
1521 if(kTRDoutN){
1522 fHistograms->FillHistogram("Resolution_E_nTRDtracklets_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
1523 fHistograms->FillHistogram("Resolution_E_nTRDtracklets_MCPt", mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
1524 fHistograms->FillHistogram("Resolution_E_nTRDclusters_ESDPt",fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDncls());
1525 fHistograms->FillHistogram("Resolution_E_nTRDclusters_MCPt",mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDncls());
1526 fHistograms->FillHistogram("Resolution_E_TRDsignal_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDsignal());
1527 }
1528
1529 Double_t mcPpt = fV0Reader->GetPositiveMCParticle()->Pt();
1530 Double_t resPdPt = 0;
1531 if (mcPpt > 0){
9c1cb6f7 1532 resPdPt = ((fV0Reader->GetPositiveTrackPt()-mcPpt)/mcPpt)*100.;
ca6d4600 1533 }
1534
1535 UInt_t statusP = fV0Reader->GetPositiveESDTrack()->GetStatus();
1536// AliESDtrack * posTr= fV0Reader->GetPositiveESDTrack();
1537 UInt_t kTRDoutP = (statusP & AliESDtrack::kTRDout);
1538
1539 Int_t ITSclsP = fV0Reader->GetPositiveTracknITSClusters();
1540 // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
1541 switch(ITSclsP){
1542 case 0: // 0 ITS clusters
1543 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS0", mcPpt, resPdPt);
1544 break;
1545 case 1: // 1 ITS cluster
1546 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS1", mcPpt, resPdPt);
1547 break;
1548 case 2: // 2 ITS clusters
1549 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS2", mcPpt, resPdPt);
1550 break;
1551 case 3: // 3 ITS clusters
1552 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS3", mcPpt, resPdPt);
1553 break;
1554 case 4: // 4 ITS clusters
1555 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS4", mcPpt, resPdPt);
1556 break;
1557 case 5: // 5 ITS clusters
1558 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS5", mcPpt, resPdPt);
1559 break;
1560 case 6: // 6 ITS clusters
1561 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS6", mcPpt, resPdPt);
1562 break;
1563 }
1564 //Filling histograms with respect to Positron resolution
1565 fHistograms->FillHistogram("Resolution_P_dPt_Pt", mcPpt, resPdPt);
1566 fHistograms->FillHistogram("Resolution_P_dPt_Phi", fV0Reader->GetPositiveTrackPhi(), resPdPt);
1567 if(kTRDoutP){
1568 fHistograms->FillHistogram("Resolution_P_nTRDtracklets_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
1569 fHistograms->FillHistogram("Resolution_P_nTRDtracklets_MCPt", mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
1570 fHistograms->FillHistogram("Resolution_P_nTRDclusters_ESDPt",fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDncls());
1571 fHistograms->FillHistogram("Resolution_P_nTRDclusters_MCPt",mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDncls());
1572 fHistograms->FillHistogram("Resolution_P_TRDsignal_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDsignal());
1573 }
1574
1575
9c1cb6f7 1576 Double_t resdR = 0.;
d7d7e825 1577 if(fV0Reader->GetNegativeMCParticle()->R() != 0){
9c1cb6f7 1578 resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100.;
d7d7e825 1579 }
9c1cb6f7 1580 Double_t resdRAbs = 0.;
037dc2db 1581 resdRAbs = (fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R());
1582
1583 fHistograms->FillHistogram("Resolution_dRAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdRAbs);
d7d7e825 1584 fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);
1585 fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());
1586 fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());
ca6d4600 1587 fHistograms->FillHistogram("Resolution_R_dPt", fV0Reader->GetNegativeMCParticle()->R(), resdPt);
1588
9c1cb6f7 1589 Double_t resdPhiAbs=0.;
1590 resdPhiAbs=0.;
037dc2db 1591 resdPhiAbs= (fV0Reader->GetMotherCandidatePhi()-fV0Reader->GetNegativeMCParticle()->Phi());
1592 fHistograms->FillHistogram("Resolution_dPhiAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdPhiAbs);
1593
d7d7e825 1594 }//if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22)
d7d7e825 1595 }//if(fDoMCTruth)
1596 }//while(fV0Reader->NextV0)
1597 fHistograms->FillHistogram("ESD_NumberOfSurvivingV0s", nSurvivingV0s);
1598 fHistograms->FillHistogram("ESD_NumberOfV0s", fV0Reader->GetNumberOfV0s());
b5832f95 1599 fHistograms->FillHistogram("ESD_NumberOfContributorsVtx", fV0Reader->GetNumberOfContributorsVtx());
cb90a330 1600
1601 fV0Reader->ResetV0IndexNumber();
d7d7e825 1602}
1603
1604void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){
1605 // Fill AOD with reconstructed Gamma
a0b94e5c 1606
6c84d371 1607 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
1608 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
d7d7e825 1609 //Create AOD particle object from AliKFParticle
d7d7e825 1610 //You could add directly AliKFParticle objects to the AOD, avoiding dependences with PartCorr
1611 //but this means that I have to work a little bit more in my side.
1612 //AODPWG4Particle objects are simpler and lighter, I think
332f1f44 1613 /*
1614 AliKFParticle * gammakf = dynamic_cast<AliKFParticle*>(fKFReconstructedGammasTClone->At(gammaIndex));
d7d7e825 1615 AliAODPWG4Particle gamma = AliAODPWG4Particle(gammakf->Px(),gammakf->Py(),gammakf->Pz(), gammakf->E());
332f1f44 1616 //gamma.SetLabel(-1);//How to get the MC label of the reconstructed gamma?
1617 gamma.SetTrackLabel( fElectronv1[gammaIndex], fElectronv2[gammaIndex] ); //How to get the MC label of the 2 electrons that form the gamma?
a0b94e5c 1618 gamma.SetDetector("CTS"); //tag the gamma as reconstructed in the central barrel
332f1f44 1619 gamma.SetPdg(AliPID::kEleCon); //photon id
a0b94e5c 1620 gamma.SetTag(-1); //Here I usually put a flag saying that montecarlo says it is prompt, decay fragmentation photon, or hadrons or whatever
332f1f44 1621 gamma.SetChi2(gammakf->Chi2());
a0b94e5c 1622 Int_t i = fAODBranch->GetEntriesFast();
1623 new((*fAODBranch)[i]) AliAODPWG4Particle(gamma);
d7d7e825 1624 */
332f1f44 1625
6c84d371 1626 AliKFParticle * gammakf = (AliKFParticle *)fKFReconstructedGammasTClone->At(gammaIndex);
d7d7e825 1627 AliGammaConversionAODObject aodObject;
1628 aodObject.SetPx(gammakf->GetPx());
1629 aodObject.SetPy(gammakf->GetPy());
1630 aodObject.SetPz(gammakf->GetPz());
1631 aodObject.SetLabel1(fElectronv1[gammaIndex]);
1632 aodObject.SetLabel2(fElectronv2[gammaIndex]);
332f1f44 1633 aodObject.SetChi2(gammakf->Chi2());
d7d7e825 1634 Int_t i = fAODBranch->GetEntriesFast();
6c84d371 1635 new((*fAODBranch)[i]) AliGammaConversionAODObject(aodObject);
a0b94e5c 1636 }
1637
d7d7e825 1638}
1639
6272370b 1640void AliAnalysisTaskGammaConversion::ProcessGammasForOmegaMesonAnalysis(){
1641 // omega meson analysis pi0+gamma decay
1642 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
48682642 1643 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
6272370b 1644 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
48682642 1645
6272370b 1646 AliKFParticle * omegaCandidateGammaDaughter = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
1647 if(fGammav1[firstPi0Index]==firstGammaIndex || fGammav2[firstPi0Index]==firstGammaIndex){
1648 continue;
1649 }
1650
48682642 1651 AliKFParticle omegaCandidate(*omegaCandidatePi0Daughter,*omegaCandidateGammaDaughter);
6272370b 1652 Double_t massOmegaCandidate = 0.;
1653 Double_t widthOmegaCandidate = 0.;
1654
48682642 1655 omegaCandidate.GetMass(massOmegaCandidate,widthOmegaCandidate);
6272370b 1656
48682642 1657 fHistograms->FillHistogram("ESD_Omega_InvMass_vs_Pt",massOmegaCandidate ,omegaCandidate.GetPt());
6272370b 1658 fHistograms->FillHistogram("ESD_Omega_InvMass",massOmegaCandidate);
48682642 1659
1660 //delete omegaCandidate;
1661
1662 }// end of omega reconstruction in pi0+gamma channel
1663
1664 if(fDoJet == kTRUE){
1665 AliKFParticle* negPiKF=NULL;
1666 AliKFParticle* posPiKF=NULL;
1667
1668 // look at the pi+pi+pi0 channel
1669 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
1670 AliESDtrack* posTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
1671 if (posTrack->GetSign()<0) continue;
c0d9051e 1672 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion))>2.) continue;
48682642 1673 if (posPiKF) delete posPiKF; posPiKF=NULL;
1674 posPiKF = new AliKFParticle( *(posTrack) ,211);
1675
1676 for(Int_t jCh=0;jCh<fChargedParticles->GetEntriesFast();jCh++){
1677 AliESDtrack* negTrack = (AliESDtrack*)(fChargedParticles->At(jCh));
1678 if( negTrack->GetSign()>0) continue;
c0d9051e 1679 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion))>2.) continue;
48682642 1680 if (negPiKF) delete negPiKF; negPiKF=NULL;
1681 negPiKF = new AliKFParticle( *(negTrack) ,-211);
1682 AliKFParticle omegaCandidatePipPinPi0(*omegaCandidatePi0Daughter,*posPiKF,*negPiKF);
1683 Double_t massOmegaCandidatePipPinPi0 = 0.;
1684 Double_t widthOmegaCandidatePipPinPi0 = 0.;
1685
1686 omegaCandidatePipPinPi0.GetMass(massOmegaCandidatePipPinPi0,widthOmegaCandidatePipPinPi0);
1687
1688 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass_vs_Pt",massOmegaCandidatePipPinPi0 ,omegaCandidatePipPinPi0.GetPt());
1689 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass",massOmegaCandidatePipPinPi0);
1690
1691 // delete omegaCandidatePipPinPi0;
1692 }
1693 }
1694 } // checking ig gammajet because in that case the chargedparticle list is created
1695
1696
6272370b 1697
6272370b 1698 }
48682642 1699
1700 if(fCalculateBackground){
1701 // Background calculation for the omega
1702 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
1703 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
1704 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
1705 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
1706 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
1707 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
1708 AliKFParticle * omegaBckCandidate = new AliKFParticle(*omegaCandidatePi0Daughter,previousGoodV0);
1709 Double_t massOmegaBckCandidate = 0.;
1710 Double_t widthOmegaBckCandidate = 0.;
1711
1712 omegaBckCandidate->GetMass(massOmegaBckCandidate,widthOmegaBckCandidate);
1713 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass_vs_Pt",massOmegaBckCandidate ,omegaBckCandidate->GetPt());
1714 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass",massOmegaBckCandidate);
1715
1716 delete omegaBckCandidate;
1717 }
1718 }
1719 }
1720 } // end of checking if background calculation is available
6272370b 1721}
d7d7e825 1722
1723void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
1724 // see header file for documentation
1725
6c84d371 1726 // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
1727 // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
037dc2db 1728
1729 if(fKFReconstructedGammasTClone->GetEntriesFast()>fV0Reader->GetNumberOfV0s()){
1730 cout<<"Warning, number of entries in the tclone is bigger than number of v0s"<<endl;
1731 }
1732
6c84d371 1733 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
1734 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
d7d7e825 1735
6c84d371 1736 // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
1737 // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
a0b94e5c 1738
6c84d371 1739 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
1740 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
a0b94e5c 1741
d7d7e825 1742 if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
1743 continue;
1744 }
1745 if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
1746 continue;
1747 }
a0b94e5c 1748
d7d7e825 1749 AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
1750
1751 Double_t massTwoGammaCandidate = 0.;
1752 Double_t widthTwoGammaCandidate = 0.;
1753 Double_t chi2TwoGammaCandidate =10000.;
1754 twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
d707e3cf 1755 // if(twoGammaCandidate->GetNDF()>0){
1756 // chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
1757 chi2TwoGammaCandidate = twoGammaCandidate->GetChi2();
d7d7e825 1758
70ef88b5 1759 fHistograms->FillHistogram("ESD_Mother_Chi2",chi2TwoGammaCandidate);
d707e3cf 1760 //if(chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()){
d7d7e825 1761
1762 TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
1763 TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
1764
1765 Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
1766 Double_t rapidity;
96ade8ef 1767 if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() <= 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() <= 0){
1768 cout << "Error: |Pz| > E !!!! " << endl;
d7d7e825 1769 rapidity=0;
1770 }
1771 else{
1772 rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
1773 }
48682642 1774 Double_t alfa=0.0;
1775 if( (twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()) != 0){
1776 alfa=TMath::Abs((twoGammaDecayCandidateDaughter0->GetE()-twoGammaDecayCandidateDaughter1->GetE())
1777 /(twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()));
1778 }
d7d7e825 1779
1e7846f4 1780 if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
1781 delete twoGammaCandidate;
1782 continue; // minimum opening angle to avoid using ghosttracks
1783 }
1784
d7d7e825 1785 fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
1786 fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
1787 fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
1788 fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
1789 fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
1790 fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
1791 fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
48682642 1792 fHistograms->FillHistogram("ESD_Mother_alfa", alfa);
d7d7e825 1793 fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
1794 fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
1795 fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
1796 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
67381a40 1797 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
9c1cb6f7 1798 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1799 }
d7d7e825 1800 fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
9640a3d1 1801
10e3319b 1802 /* Kenneth, just for testing*/
2c9fed75 1803 AliGammaConversionBGHandler * bgHandlerTest = fV0Reader->GetBGHandler();
10e3319b 1804
2c9fed75 1805 Int_t zbin= bgHandlerTest->GetZBinIndex(fV0Reader->GetVertexZ());
1806 Int_t mbintest= bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
10e3319b 1807
2c9fed75 1808 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass",zbin,mbintest),massTwoGammaCandidate);
10e3319b 1809
1810 /* end Kenneth, just for testing*/
1811
9c1cb6f7 1812 if(fCalculateBackground){
1813 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
1814 Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
1815 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1816 }
ebcfaa7e 1817 // if(fDoNeutralMesonV0MCCheck){
1818 if(fDoMCTruth){
037dc2db 1819 //Kenneth: Checking the eta of the gamma to check the difference between 0.9 and 1.2
1820 Int_t indexKF1 = fKFReconstructedGammasV0Index.at(firstGammaIndex);
1821 if(indexKF1<fV0Reader->GetNumberOfV0s()){
1822 fV0Reader->GetV0(indexKF1);//updates to the correct v0
1823 Double_t eta1 = fV0Reader->GetMotherCandidateEta();
1824 Bool_t isRealPi0=kFALSE;
10e3319b 1825 Bool_t isRealEta=kFALSE;
037dc2db 1826 Int_t gamma1MotherLabel=-1;
1827 if(fV0Reader->HasSameMCMother() == kTRUE){
1828 //cout<<"This v0 is a real v0!!!!"<<endl;
1829 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1830 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1831 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
1832 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
1833 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
1834 gamma1MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
1835 }
1836 }
1837 }
1838 }
1839 Int_t indexKF2 = fKFReconstructedGammasV0Index.at(secondGammaIndex);
1840 if(indexKF1 == indexKF2){
1841 cout<<"index of the two KF particles are the same.... should not happen"<<endl;
1842 }
1843 if(indexKF2<fV0Reader->GetNumberOfV0s()){
1844 fV0Reader->GetV0(indexKF2);
1845 Double_t eta2 = fV0Reader->GetMotherCandidateEta();
1846 Int_t gamma2MotherLabel=-1;
1847 if(fV0Reader->HasSameMCMother() == kTRUE){
1848 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1849 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1850 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
1851 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
1852 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
1853 gamma2MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
1854 }
1855 }
1856 }
1857 }
1858 if(gamma1MotherLabel>=0 && gamma1MotherLabel==gamma2MotherLabel){
1859 if(fV0Reader->CheckIfPi0IsMother(gamma1MotherLabel)){
1860 isRealPi0=kTRUE;
1861 }
10e3319b 1862 if(fV0Reader->CheckIfEtaIsMother(gamma1MotherLabel)){
1863 isRealEta=kTRUE;
1864 }
1865
037dc2db 1866 }
1867
1868 if(TMath::Abs(eta1)>0.9 && TMath::Abs(eta2)>0.9){
ebcfaa7e 1869 // fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
1870 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
10e3319b 1871 if(isRealPi0 || isRealEta){
037dc2db 1872 fHistograms->FillHistogram("ESD_TruePi0_InvMass_1212",massTwoGammaCandidate);
1873 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_1212",openingAngleTwoGammaCandidate);
1874 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 1875 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1876 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 1877 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
1878 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1879 }
037dc2db 1880 }
1881 }
1882 else if(TMath::Abs(eta1)>0.9 || TMath::Abs(eta2)>0.9){
ebcfaa7e 1883 // fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
1884 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
10e3319b 1885 if(isRealPi0 || isRealEta){
037dc2db 1886 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0912",massTwoGammaCandidate);
1887 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0912",openingAngleTwoGammaCandidate);
1888 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 1889 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1890 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 1891 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
1892 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1893 }
037dc2db 1894 }
1895 }
1896 else{
ebcfaa7e 1897 // fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
1898 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
10e3319b 1899 if(isRealPi0 || isRealEta){
037dc2db 1900 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0909",massTwoGammaCandidate);
1901 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0909",openingAngleTwoGammaCandidate);
1902 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 1903 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1904 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 1905 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
1906 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1907 }
037dc2db 1908 }
1909 }
1910 }
1911 }
1912 }
9640a3d1 1913 if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
1914 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1915 fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
1916 }
ebcfaa7e 1917
1918 if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
1919 fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
1920 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
1921 }
1922 else if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 || TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
1923 fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
1924 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
1925 }
1926 else{
1927 fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
1928 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
1929 }
6272370b 1930 Double_t lowMassPi0=0.1;
1931 Double_t highMassPi0=0.15;
1932 if (massTwoGammaCandidate > lowMassPi0 && massTwoGammaCandidate < highMassPi0 ){
1933 new((*fKFReconstructedPi0sTClone)[fKFReconstructedPi0sTClone->GetEntriesFast()]) AliKFParticle(*twoGammaCandidate);
1934 fGammav1.push_back(firstGammaIndex);
1935 fGammav2.push_back(secondGammaIndex);
1936 }
1937
d707e3cf 1938 // }
1939 //}
1940 delete twoGammaCandidate;
d7d7e825 1941 }
1942 }
1943}
1944
6272370b 1945void AliAnalysisTaskGammaConversion::ProcessConvPHOSGammasForNeutralMesonAnalysis(){
48682642 1946/*
6272370b 1947 // see header file for documentation
1948 // Analyse Pi0 with one photon from Phos and 1 photon from conversions
1949
1950
1951
1952 Double_t vtx[3];
1953 vtx[0] = fV0Reader->GetPrimaryVertex()->GetX();
1954 vtx[1] = fV0Reader->GetPrimaryVertex()->GetY();
1955 vtx[2] = fV0Reader->GetPrimaryVertex()->GetZ();
1956
1957
1958 // Loop over all CaloClusters and consider only the PHOS ones:
1959 AliESDCaloCluster *clu;
1960 TLorentzVector pPHOS;
1961 TLorentzVector gammaPHOS;
1962 TLorentzVector gammaGammaConv;
1963 TLorentzVector pi0GammaConvPHOS;
1964 TLorentzVector gammaGammaConvBck;
1965 TLorentzVector pi0GammaConvPHOSBck;
1966
1967
1968 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
1969 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
1970 if ( !clu->IsPHOS() || clu->E()<0.1 ) continue;
1971 clu ->GetMomentum(pPHOS ,vtx);
1972 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
1973 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
1974 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz(),0.);
1975 gammaPHOS.SetXYZM(pPHOS.Px(),pPHOS.Py(),pPHOS.Pz(),0.);
1976 pi0GammaConvPHOS=gammaGammaConv+gammaPHOS;
1977 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS",pi0GammaConvPHOS.M());
1978 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvPHOS",pi0GammaConvPHOS.M(),pi0GammaConvPHOS.Pt());
1979
1980 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz());
1981 TVector3 v3D1(gammaPHOS.Px(),gammaPHOS.Py(),gammaPHOS.Pz());
1982 Double_t opanConvPHOS= v3D0.Angle(v3D1);
1983 if ( opanConvPHOS < 0.35){
1984 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanLow",pi0GammaConvPHOS.M());
1985 }else{
1986 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanHigh",pi0GammaConvPHOS.M());
1987 }
1988
1989 }
1990
1991 // Now the LorentVector pPHOS is obtained and can be paired with the converted proton
1992 }
1993 //==== End of the PHOS cluster selection ============
1994 TLorentzVector pEMCAL;
1995 TLorentzVector gammaEMCAL;
1996 TLorentzVector pi0GammaConvEMCAL;
1997 TLorentzVector pi0GammaConvEMCALBck;
1998
1999 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
2000 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2001 if ( !clu->IsEMCAL() || clu->E()<0.1 ) continue;
2002 if (clu->GetNCells() <= 1) continue;
2003 if ( clu->GetTOF()*1e9 < 550 || clu->GetTOF()*1e9 > 750) continue;
2004
2005 clu ->GetMomentum(pEMCAL ,vtx);
2006 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2007 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2008 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),
2009 twoGammaDecayCandidateDaughter0->Py(),
2010 twoGammaDecayCandidateDaughter0->Pz(),0.);
2011 gammaEMCAL.SetXYZM(pEMCAL.Px(),pEMCAL.Py(),pEMCAL.Pz(),0.);
2012 pi0GammaConvEMCAL=gammaGammaConv+gammaEMCAL;
2013 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL",pi0GammaConvEMCAL.M());
2014 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL",pi0GammaConvEMCAL.M(),pi0GammaConvEMCAL.Pt());
2015 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),
2016 twoGammaDecayCandidateDaughter0->Py(),
2017 twoGammaDecayCandidateDaughter0->Pz());
2018 TVector3 v3D1(gammaEMCAL.Px(),gammaEMCAL.Py(),gammaEMCAL.Pz());
2019
2020
2021 Double_t opanConvEMCAL= v3D0.Angle(v3D1);
2022 if ( opanConvEMCAL < 0.35){
2023 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanLow",pi0GammaConvEMCAL.M());
2024 }else{
2025 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanHigh",pi0GammaConvEMCAL.M());
2026 }
2027
2028 }
48682642 2029 if(fCalculateBackground){
2030 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2031 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
2032 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2033 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2034 gammaGammaConvBck.SetXYZM(previousGoodV0.Px(),
2035 previousGoodV0.Py(),
2036 previousGoodV0.Pz(),0.);
2037 pi0GammaConvEMCALBck=gammaGammaConvBck+gammaEMCAL;
2038 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M());
2039 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M(),
2040 pi0GammaConvEMCALBck.Pt());
2041 }
6272370b 2042 }
48682642 2043
2044 // Now the LorentVector pEMCAL is obtained and can be paired with the converted proton
2045 } // end of checking if background photons are available
2046 }
6272370b 2047 //==== End of the PHOS cluster selection ============
48682642 2048*/
6272370b 2049}
2050
d7d7e825 2051void AliAnalysisTaskGammaConversion::CalculateBackground(){
2052 // see header file for documentation
5e55d806 2053
2054
2055 TClonesArray * currentEventV0s = fV0Reader->GetCurrentEventGoodV0s();
2056
2057 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2058 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
2059 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2060 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2061 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2062 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2063
2064 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
d7d7e825 2065
5e55d806 2066 Double_t massBG =0.;
2067 Double_t widthBG = 0.;
2068 Double_t chi2BG =10000.;
2069 backgroundCandidate->GetMass(massBG,widthBG);
d707e3cf 2070 // if(backgroundCandidate->GetNDF()>0){
2071 // chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
2072 chi2BG = backgroundCandidate->GetChi2();
2073 // if(chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()){
d7d7e825 2074
5e55d806 2075 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2076 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
d7d7e825 2077
5e55d806 2078 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
a0b94e5c 2079
5e55d806 2080 Double_t rapidity;
96ade8ef 2081
2082 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() <= 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() <= 0){
2083 cout << "Error: |Pz| > E !!!! " << endl;
2084 rapidity=0;
2085 } else {
2086 rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2087 }
d7d7e825 2088
9c1cb6f7 2089 Double_t alfa=0.0;
2090 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
2091 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
2092 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
2093 }
2094
d7d7e825 2095
5e55d806 2096 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2097 delete backgroundCandidate;
2098 continue; // minimum opening angle to avoid using ghosttracks
2099 }
037dc2db 2100
2101 // original
5e55d806 2102 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2103 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2104 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2105 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2106 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2107 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2108 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2109 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2110 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2111 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2112 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2113 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
2114
67381a40 2115 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
9c1cb6f7 2116 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
2117 }
5e55d806 2118 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2119 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
2120 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
2121 }
037dc2db 2122
2123
2124 // test
2125 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2126
2127 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
2128 Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2129
2130 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
2131 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
2132 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
2133 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
2134 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
2135 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
2136 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
2137 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2138 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2139 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
2140 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2141 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
9640a3d1 2142
037dc2db 2143 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2144 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2145 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
2146 }
d707e3cf 2147 // }
2148 // }
5e55d806 2149 delete backgroundCandidate;
d7d7e825 2150 }
d7d7e825 2151 }
2152 }
2153}
2154
2155
d7d7e825 2156void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){
2157 //ProcessGammasForGammaJetAnalysis
a0b94e5c 2158
d7d7e825 2159 Double_t distIsoMin;
a0b94e5c 2160
d7d7e825 2161 CreateListOfChargedParticles();
a0b94e5c 2162
2163
6c84d371 2164 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
2165 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
2166 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
01b7fdcc 2167 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 2168
01b7fdcc 2169 if( momentumVectorCurrentGamma.Pt()> fMinPtForGammaJet){
d7d7e825 2170 distIsoMin=GetMinimumDistanceToCharge(gammaIndex);
a0b94e5c 2171
d7d7e825 2172 if (distIsoMin > fMinIsoConeSize && fLeadingChargedIndex>=0){
01b7fdcc 2173 CalculateJetCone(gammaIndex);
d7d7e825 2174 }
2175 }
2176 }
2177}
2178
6272370b 2179//____________________________________________________________________
2180Bool_t AliAnalysisTaskGammaConversion::IsGoodImpPar(AliESDtrack *const track)
2181{
2182//
2183// check whether particle has good DCAr(Pt) impact
2184// parameter. Only for TPC+ITS tracks (7*sigma cut)
2185// Origin: Andrea Dainese
2186//
2187
2188Float_t d0z0[2],covd0z0[3];
2189track->GetImpactParameters(d0z0,covd0z0);
2190Float_t sigma= 0.0050+0.0060/TMath::Power(track->Pt(),0.9);
2191Float_t d0max = 7.*sigma;
2192if(TMath::Abs(d0z0[0]) < d0max) return kTRUE;
2193
2194return kFALSE;
2195}
2196
2197
d7d7e825 2198void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){
2199 // CreateListOfChargedParticles
a0b94e5c 2200
d7d7e825 2201 fESDEvent = fV0Reader->GetESDEvent();
037dc2db 2202 Int_t numberOfESDTracks=0;
d7d7e825 2203 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
2204 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 2205
d7d7e825 2206 if(!curTrack){
2207 continue;
2208 }
d707e3cf 2209 // Not needed if Standard function used.
2210// if(!IsGoodImpPar(curTrack)){
2211// continue;
2212// }
a0b94e5c 2213
d7d7e825 2214 if(fEsdTrackCuts->AcceptTrack(curTrack) ){
6c84d371 2215 new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack);
2216 // fChargedParticles.push_back(curTrack);
d7d7e825 2217 fChargedParticlesId.push_back(iTracks);
037dc2db 2218 numberOfESDTracks++;
d7d7e825 2219 }
2220 }
037dc2db 2221 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",numberOfESDTracks);
e40fd7e2 2222
2223 if (fV0Reader->GetNumberOfContributorsVtx()>=1){
2224 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",numberOfESDTracks);
2225 }
d7d7e825 2226}
9c1cb6f7 2227void AliAnalysisTaskGammaConversion::RecalculateV0ForGamma(){
2228
2229 Double_t massE=0.00051099892;
2230 TLorentzVector curElecPos;
2231 TLorentzVector curElecNeg;
2232 TLorentzVector curGamma;
2233
2234 TLorentzVector curGammaAt;
2235 TLorentzVector curElecPosAt;
2236 TLorentzVector curElecNegAt;
2237 AliKFVertex primVtxGamma(*(fESDEvent->GetPrimaryVertex()));
2238 AliKFVertex primVtxImprovedGamma = primVtxGamma;
2239
2240 const AliESDVertex *vtxT3D=fESDEvent->GetPrimaryVertex();
2241
2242 Double_t xPrimaryVertex=vtxT3D->GetXv();
2243 Double_t yPrimaryVertex=vtxT3D->GetYv();
2244 Double_t zPrimaryVertex=vtxT3D->GetZv();
70ef88b5 2245 // Float_t primvertex[3]={xPrimaryVertex,yPrimaryVertex,zPrimaryVertex};
9c1cb6f7 2246
2247 Float_t nsigmaTPCtrackPos;
2248 Float_t nsigmaTPCtrackNeg;
2249 Float_t nsigmaTPCtrackPosToPion;
2250 Float_t nsigmaTPCtrackNegToPion;
2251 AliKFParticle* negKF=NULL;
2252 AliKFParticle* posKF=NULL;
2253
2254 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
2255 AliESDtrack* posTrack = fESDEvent->GetTrack(iTracks);
2256 if(!posTrack){
2257 continue;
2258 }
2259 if (posKF) delete posKF; posKF=NULL;
2260 if(posTrack->GetSign()<0) continue;
2261 if(!(posTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
2262 if(posTrack->GetKinkIndex(0)>0 ) continue;
2263 if(posTrack->GetNcls(1)<50)continue;
2264 Double_t momPos[3];
2265 // posTrack->GetConstrainedPxPyPz(momPos);
2266 posTrack->GetPxPyPz(momPos);
2267 AliESDtrack *ptrk=fESDEvent->GetTrack(iTracks);
2268 curElecPos.SetXYZM(momPos[0],momPos[1],momPos[2],massE);
2269 if(TMath::Abs(curElecPos.Eta())<0.9) continue;
2270 posKF = new AliKFParticle( *(posTrack),-11);
2271
2272 nsigmaTPCtrackPos = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kElectron);
2273 nsigmaTPCtrackPosToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion);
2274
2275 if ( nsigmaTPCtrackPos>5.|| nsigmaTPCtrackPos<-2.){
2276 continue;
2277 }
2278
2279 if(pow((momPos[0]*momPos[0]+momPos[1]*momPos[1]+momPos[2]*momPos[2]),0.5)>0.5 && nsigmaTPCtrackPosToPion<1){
2280 continue;
2281 }
2282
2283
2284
2285 for(Int_t jTracks = 0; jTracks < fESDEvent->GetNumberOfTracks(); jTracks++){
2286 AliESDtrack* negTrack = fESDEvent->GetTrack(jTracks);
2287 if(!negTrack){
2288 continue;
2289 }
2290 if (negKF) delete negKF; negKF=NULL;
2291 if(negTrack->GetSign()>0) continue;
2292 if(!(negTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
2293 if(negTrack->GetKinkIndex(0)>0 ) continue;
2294 if(negTrack->GetNcls(1)<50)continue;
2295 Double_t momNeg[3];
2296 // negTrack->GetConstrainedPxPyPz(momNeg);
2297 negTrack->GetPxPyPz(momNeg);
2298
2299 nsigmaTPCtrackNeg = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kElectron);
2300 nsigmaTPCtrackNegToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion);
2301 if ( nsigmaTPCtrackNeg>5. || nsigmaTPCtrackNeg<-2.){
2302 continue;
2303 }
2304 if(pow((momNeg[0]*momNeg[0]+momNeg[1]*momNeg[1]+momNeg[2]*momNeg[2]),0.5)>0.5 && nsigmaTPCtrackNegToPion<1){
2305 continue;
2306 }
2307 AliESDtrack *ntrk=fESDEvent->GetTrack(jTracks);
2308 curElecNeg.SetXYZM(momNeg[0],momNeg[1],momNeg[2],massE);
2309 if(TMath::Abs(curElecNeg.Eta())<0.9) continue;
2310 negKF = new AliKFParticle( *(negTrack) ,11);
2311
2312 Double_t b=fESDEvent->GetMagneticField();
2313 Double_t xn, xp, dca=ntrk->GetDCA(ptrk,b,xn,xp);
2314 AliExternalTrackParam nt(*ntrk), pt(*ptrk);
2315 nt.PropagateTo(xn,b); pt.PropagateTo(xp,b);
2316
2317
2318 //--- Like in ITSV0Finder
2319 AliExternalTrackParam ntAt0(*ntrk), ptAt0(*ptrk);
2320 Double_t xxP,yyP,alphaP;
2321 Double_t rP[3];
2322
2323 // if (!ptAt0.GetGlobalXYZat(ptAt0->GetX(),xxP,yyP,zzP)) continue;
2324 if (!ptAt0.GetXYZAt(ptAt0.GetX(),b,rP)) continue;
2325 xxP=rP[0];
2326 yyP=rP[1];
2327 alphaP = TMath::ATan2(yyP,xxP);
2328
2329
2330 ptAt0.Propagate(alphaP,0,b);
2331 Float_t ptfacP = (1.+100.*TMath::Abs(ptAt0.GetC(b)));
2332
2333 // Double_t distP = ptAt0.GetY();
2334 Double_t normP = ptfacP*TMath::Sqrt(ptAt0.GetSigmaY2());
2335 Double_t normdist0P = TMath::Abs(ptAt0.GetY()/normP);
2336 Double_t normdist1P = TMath::Abs((ptAt0.GetZ()-zPrimaryVertex)/(ptfacP*TMath::Sqrt(ptAt0.GetSigmaZ2())));
2337 Double_t normdistP = TMath::Sqrt(normdist0P*normdist0P+normdist1P*normdist1P);
2338
2339
2340 Double_t xxN,yyN,alphaN;
2341 Double_t rN[3];
2342 // if (!ntAt0.GetGlobalXYZat(ntAt0->GetX(),xxN,yyN,zzN)) continue;
2343 if (!ntAt0.GetXYZAt(ntAt0.GetX(),b,rN)) continue;
2344 xxN=rN[0];
2345 yyN=rN[1];
2346
2347 alphaN = TMath::ATan2(yyN,xxN);
2348
2349 ntAt0.Propagate(alphaN,0,b);
2350
2351 Float_t ptfacN = (1.+100.*TMath::Abs(ntAt0.GetC(b)));
2352 // Double_t distN = ntAt0.GetY();
2353 Double_t normN = ptfacN*TMath::Sqrt(ntAt0.GetSigmaY2());
2354 Double_t normdist0N = TMath::Abs(ntAt0.GetY()/normN);
2355 Double_t normdist1N = TMath::Abs((ntAt0.GetZ()-zPrimaryVertex)/(ptfacN*TMath::Sqrt(ntAt0.GetSigmaZ2())));
2356 Double_t normdistN = TMath::Sqrt(normdist0N*normdist0N+normdist1N*normdist1N);
2357
2358 //-----------------------------
2359
2360 Double_t momNegAt[3];
2361 nt.GetPxPyPz(momNegAt);
2362 curElecNegAt.SetXYZM(momNegAt[0],momNegAt[1],momNegAt[2],massE);
2363
2364 Double_t momPosAt[3];
2365 pt.GetPxPyPz(momPosAt);
2366 curElecPosAt.SetXYZM(momPosAt[0],momPosAt[1],momPosAt[2],massE);
2367 if(dca>1){
2368 continue;
2369 }
2370
2371 // Double_t dneg= negTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
2372 // Double_t dpos= posTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
2373 AliESDv0 vertex(nt,jTracks,pt,iTracks);
2374
2375
2376 Float_t cpa=vertex.GetV0CosineOfPointingAngle(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex);
2377
9c1cb6f7 2378
70ef88b5 2379
9c1cb6f7 2380 // cout<< "v0Rr::"<< v0Rr<<endl;
2381 // if (pvertex.GetRr()<0.5){
2382 // continue;
2383 //}
2384 // cout<<"vertex.GetChi2V0()"<<vertex.GetChi2V0()<<endl;
2385 if(cpa<0.9)continue;
2386 // if (vertex.GetChi2V0() > 30) continue;
2387 // cout<<"xp+xn::"<<xp<<" "<<xn<<endl;
2388 if ((xn+xp) < 0.4) continue;
2389 if (TMath::Abs(ntrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05)
2390 if (TMath::Abs(ptrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05) continue;
2391
2392 //cout<<"pass"<<endl;
2393
2394 AliKFParticle v0GammaC;
2395 v0GammaC+=(*negKF);
2396 v0GammaC+=(*posKF);
2397 v0GammaC.SetMassConstraint(0,0.001);
2398 primVtxImprovedGamma+=v0GammaC;
2399 v0GammaC.SetProductionVertex(primVtxImprovedGamma);
2400
2401
2402 curGamma=curElecNeg+curElecPos;
2403 curGammaAt=curElecNegAt+curElecPosAt;
2404
2405 // invariant mass versus pt of K0short
2406
2407 Double_t chi2V0GammaC=100000.;
2408 if( v0GammaC.GetNDF() != 0) {
2409 chi2V0GammaC = v0GammaC.GetChi2()/v0GammaC.GetNDF();
2410 }else{
2411 cout<< "ERROR::v0K0C.GetNDF()" << endl;
2412 }
2413
2414 if(chi2V0GammaC<200 &&chi2V0GammaC>0 ){
2415 if(fHistograms != NULL){
2416 fHistograms->FillHistogram("ESD_RecalculateV0_InvMass",v0GammaC.GetMass());
2417 fHistograms->FillHistogram("ESD_RecalculateV0_Pt",v0GammaC.GetPt());
2418 fHistograms->FillHistogram("ESD_RecalculateV0_E_dEdxP",curElecNegAt.P(),negTrack->GetTPCsignal());
2419 fHistograms->FillHistogram("ESD_RecalculateV0_P_dEdxP",curElecPosAt.P(),posTrack->GetTPCsignal());
2420 fHistograms->FillHistogram("ESD_RecalculateV0_cpa",cpa);
2421 fHistograms->FillHistogram("ESD_RecalculateV0_dca",dca);
9c1cb6f7 2422 fHistograms->FillHistogram("ESD_RecalculateV0_normdistP",normdistP);
2423 fHistograms->FillHistogram("ESD_RecalculateV0_normdistN",normdistN);
2424
2425 new((*fKFRecalculatedGammasTClone)[fKFRecalculatedGammasTClone->GetEntriesFast()]) AliKFParticle(v0GammaC);
2426 fElectronRecalculatedv1.push_back(iTracks);
2427 fElectronRecalculatedv2.push_back(jTracks);
2428 }
2429 }
2430 }
2431 }
2432
2433 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();firstGammaIndex++){
2434 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();secondGammaIndex++){
2435 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(firstGammaIndex);
2436 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(secondGammaIndex);
2437
2438 if(fElectronRecalculatedv1[firstGammaIndex]==fElectronRecalculatedv1[secondGammaIndex]){
2439 continue;
2440 }
2441 if( fElectronRecalculatedv2[firstGammaIndex]==fElectronRecalculatedv2[secondGammaIndex]){
2442 continue;
2443 }
2444
2445 AliKFParticle twoGammaCandidate(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
2446 fHistograms->FillHistogram("ESD_RecalculateGG_InvMass",twoGammaCandidate.GetMass());
2447 fHistograms->FillHistogram("ESD_RecalculateGG_InvMass_vs_Pt",twoGammaCandidate.GetMass(),twoGammaCandidate.GetPt());
2448 }
2449 }
2450}
01b7fdcc 2451void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){
6c84d371 2452 // CaculateJetCone
a0b94e5c 2453
d7d7e825 2454 Double_t cone;
2455 Double_t coneSize=0.3;
2456 Double_t ptJet=0;
a0b94e5c 2457
6c84d371 2458 // AliKFParticle * currentGamma = &fKFReconstructedGammas[gammaIndex];
2459 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
98778c17 2460
01b7fdcc 2461 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 2462
6c84d371 2463 AliESDtrack* leadingCharged = (AliESDtrack*)(fChargedParticles->At(fLeadingChargedIndex));
98778c17 2464
d7d7e825 2465 Double_t momLeadingCharged[3];
2466 leadingCharged->GetConstrainedPxPyPz(momLeadingCharged);
a0b94e5c 2467
01b7fdcc 2468 TVector3 momentumVectorLeadingCharged(momLeadingCharged[0],momLeadingCharged[1],momLeadingCharged[2]);
a0b94e5c 2469
01b7fdcc 2470 Double_t phi1=momentumVectorLeadingCharged.Phi();
2471 Double_t eta1=momentumVectorLeadingCharged.Eta();
2472 Double_t phi3=momentumVectorCurrentGamma.Phi();
a0b94e5c 2473
6c84d371 2474 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
2475 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 2476 Int_t chId = fChargedParticlesId[iCh];
2477 if(fLeadingChargedIndex==chId || fLeadingChargedIndex==chId) continue;
2478 Double_t mom[3];
2479 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 2480 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
2481 Double_t phi2=momentumVectorChargedParticle.Phi();
2482 Double_t eta2=momentumVectorChargedParticle.Eta();
a0b94e5c 2483
2484
d7d7e825 2485 cone=100.;
2486 if( TMath::Abs(phi2 - phi1) <= ( TMath::TwoPi()-coneSize) ){
2487 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-phi1),2) );
2488 }else{
2489 if( (phi2 - phi1)> TMath::TwoPi()-coneSize ){
2490 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-TMath::TwoPi()-phi1),2) );
2491 }
2492 if( (phi2 - phi1)< -(TMath::TwoPi()-coneSize) ){
2493 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2+TMath::TwoPi()-phi1),2) );
2494 }
2495 }
a0b94e5c 2496
01b7fdcc 2497 if(cone <coneSize&& momentumVectorChargedParticle.Pt()>fMinPtJetCone ){
2498 ptJet+= momentumVectorChargedParticle.Pt();
2499 Double_t ffzHdrGam = momentumVectorChargedParticle.Pt()/momentumVectorCurrentGamma.Pt();
2500 Double_t imbalanceHdrGam=-momentumVectorChargedParticle.Dot(momentumVectorCurrentGamma)/momentumVectorCurrentGamma.Mag2();
d7d7e825 2501 fHistograms->FillHistogram("ESD_FFzHdrGam",ffzHdrGam);
2502 fHistograms->FillHistogram("ESD_ImbalanceHdrGam",imbalanceHdrGam);
a0b94e5c 2503
d7d7e825 2504 }
a0b94e5c 2505
d7d7e825 2506 Double_t dphiHdrGam=phi3-phi2;
2507 if ( dphiHdrGam < (-TMath::PiOver2())){
2508 dphiHdrGam+=(TMath::TwoPi());
2509 }
a0b94e5c 2510
d7d7e825 2511 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
2512 dphiHdrGam-=(TMath::TwoPi());
2513 }
a0b94e5c 2514
01b7fdcc 2515 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 2516 fHistograms->FillHistogram("ESD_dphiHdrGamIsolated",dphiHdrGam);
2517 }
2518 }//track loop
a0b94e5c 2519
2520
d7d7e825 2521}
2522
2523Double_t AliAnalysisTaskGammaConversion::GetMinimumDistanceToCharge(Int_t indexHighestPtGamma){
2524 // GetMinimumDistanceToCharge
a0b94e5c 2525
d7d7e825 2526 Double_t fIsoMin=100.;
2527 Double_t ptLeadingCharged=-1.;
98778c17 2528
2529 fLeadingChargedIndex=-1;
a0b94e5c 2530
6c84d371 2531 AliKFParticle * gammaHighestPt = (AliKFParticle*)fKFReconstructedGammasTClone->At(indexHighestPtGamma);
01b7fdcc 2532 TVector3 momentumVectorgammaHighestPt(gammaHighestPt->GetPx(),gammaHighestPt->GetPy(),gammaHighestPt->GetPz());
a0b94e5c 2533
01b7fdcc 2534 Double_t phi1=momentumVectorgammaHighestPt.Phi();
2535 Double_t eta1=momentumVectorgammaHighestPt.Eta();
a0b94e5c 2536
6c84d371 2537 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
2538 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 2539 Int_t chId = fChargedParticlesId[iCh];
2540 if(fElectronv1[indexHighestPtGamma]==chId || fElectronv2[indexHighestPtGamma]==chId) continue;
2541 Double_t mom[3];
2542 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 2543 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
2544 Double_t phi2=momentumVectorChargedParticle.Phi();
2545 Double_t eta2=momentumVectorChargedParticle.Eta();
d7d7e825 2546 Double_t iso=pow( (pow( (eta1-eta2),2)+ pow((phi1-phi2),2)),0.5 );
a0b94e5c 2547
01b7fdcc 2548 if(momentumVectorChargedParticle.Pt()>fMinPtIsoCone ){
d7d7e825 2549 if (iso<fIsoMin){
2550 fIsoMin=iso;
2551 }
2552 }
a0b94e5c 2553
d7d7e825 2554 Double_t dphiHdrGam=phi1-phi2;
2555 if ( dphiHdrGam < (-TMath::PiOver2())){
2556 dphiHdrGam+=(TMath::TwoPi());
2557 }
a0b94e5c 2558
d7d7e825 2559 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
2560 dphiHdrGam-=(TMath::TwoPi());
2561 }
01b7fdcc 2562 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 2563 fHistograms->FillHistogram("ESD_dphiHdrGam",dphiHdrGam);
2564 }
a0b94e5c 2565
d7d7e825 2566 if (dphiHdrGam>0.9*TMath::Pi() && dphiHdrGam<1.1*TMath::Pi()){
a0b94e5c 2567 if (momentumVectorChargedParticle.Pt()> ptLeadingCharged && momentumVectorChargedParticle.Pt()>0.1*momentumVectorgammaHighestPt.Pt()){
01b7fdcc 2568 ptLeadingCharged=momentumVectorChargedParticle.Pt();
d7d7e825 2569 fLeadingChargedIndex=iCh;
2570 }
2571 }
a0b94e5c 2572
d7d7e825 2573 }//track loop
2574 fHistograms->FillHistogram("ESD_MinimumIsoDistance",fIsoMin);
2575 return fIsoMin;
a0b94e5c 2576
d7d7e825 2577}
2578
a0b94e5c 2579Int_t AliAnalysisTaskGammaConversion::GetIndexHighestPtGamma(){
2580 //GetIndexHighestPtGamma
2581
d7d7e825 2582 Int_t indexHighestPtGamma=-1;
01b7fdcc 2583 //Double_t
2584 fGammaPtHighest = -100.;
a0b94e5c 2585
6c84d371 2586 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2587 AliKFParticle * gammaHighestPtCandidate = (AliKFParticle*)fKFReconstructedGammasTClone->At(firstGammaIndex);
a0b94e5c 2588 TVector3 momentumVectorgammaHighestPtCandidate(gammaHighestPtCandidate->GetPx(),gammaHighestPtCandidate->GetPy(),gammaHighestPtCandidate->GetPz());
2589 if (momentumVectorgammaHighestPtCandidate.Pt() > fGammaPtHighest){
2590 fGammaPtHighest=momentumVectorgammaHighestPtCandidate.Pt();
2591 //gammaHighestPt = gammaHighestPtCandidate;
2592 indexHighestPtGamma=firstGammaIndex;
2593 }
d7d7e825 2594 }
a0b94e5c 2595
d7d7e825 2596 return indexHighestPtGamma;
a0b94e5c 2597
d7d7e825 2598}
2599
2600
2601void AliAnalysisTaskGammaConversion::Terminate(Option_t */*option*/)
2602{
2603 // Terminate analysis
2604 //
2605 AliDebug(1,"Do nothing in Terminate");
2606}
2607
2608void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()
2609{
2610 //AOD
1e7846f4 2611 if(fAODBranch==NULL){
2612 fAODBranch = new TClonesArray("AliGammaConversionAODObject", 0);
2613 }
2614 fAODBranch->Delete();
d7d7e825 2615 fAODBranch->SetName(fAODBranchName);
332f1f44 2616
2617 //If delta AOD file name set, add in separate file. Else add in standard aod file.
2618 if(fKFDeltaAODFileName.Length() > 0) {
2619 AddAODBranch("TClonesArray", &fAODBranch, fKFDeltaAODFileName.Data());
2620 AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(fKFDeltaAODFileName.Data());
2621 } else {
2622 AddAODBranch("TClonesArray", &fAODBranch);
2623 }
d7d7e825 2624 // Create the output container
2625 if(fOutputContainer != NULL){
2626 delete fOutputContainer;
2627 fOutputContainer = NULL;
2628 }
2629 if(fOutputContainer == NULL){
2630 fOutputContainer = new TList();
b58d3c74 2631 fOutputContainer->SetOwner(kTRUE);
d7d7e825 2632 }
2633
2634 //Adding the histograms to the output container
2635 fHistograms->GetOutputContainer(fOutputContainer);
2636
2637
2638 if(fWriteNtuple){
2639 if(fGammaNtuple == NULL){
2640 fGammaNtuple = new TNtuple("V0ntuple","V0ntuple","OnTheFly:HasVertex:NegPIDProb:PosPIDProb:X:Y:Z:R:MotherCandidateNDF:MotherCandidateChi2:MotherCandidateEnergy:MotherCandidateEta:MotherCandidatePt:MotherCandidateMass:MotherCandidateWidth:MCMotherCandidatePT:EPOpeningAngle:ElectronEnergy:ElectronPt:ElectronEta:ElectronPhi:PositronEnergy:PositronPt:PositronEta:PositronPhi:HasSameMCMother:MotherMCParticlePIDCode",50000);
2641 }
2642 if(fNeutralMesonNtuple == NULL){
2643 fNeutralMesonNtuple = new TNtuple("NeutralMesonNtuple","NeutralMesonNtuple","test");
2644 }
2645 TList * ntupleTList = new TList();
b58d3c74 2646 ntupleTList->SetOwner(kTRUE);
d7d7e825 2647 ntupleTList->SetName("Ntuple");
2648 ntupleTList->Add((TNtuple*)fGammaNtuple);
2649 fOutputContainer->Add(ntupleTList);
2650 }
2651
2652 fOutputContainer->SetName(GetName());
2653}
2654
2655Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(TParticle* const daughter0, TParticle* const daughter1) const{
2656 //helper function
2657 TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());
2658 TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());
2659 return v3D0.Angle(v3D1);
2660}
2661
2662void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){
2663 // see header file for documentation
5e55d806 2664
d7d7e825 2665 vector<Int_t> indexOfGammaParticle;
a0b94e5c 2666
d7d7e825 2667 fStack = fV0Reader->GetMCStack();
a0b94e5c 2668
d7d7e825 2669 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
2670 return; // aborts if the primary vertex does not have contributors.
2671 }
a0b94e5c 2672
d7d7e825 2673 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
2674 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
2675 if(particle->GetPdgCode()==22){ //Gamma
2676 if(particle->GetNDaughters() >= 2){
2677 TParticle* electron=NULL;
2678 TParticle* positron=NULL;
2679 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
2680 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
2681 if(tmpDaughter->GetUniqueID() == 5){
2682 if(tmpDaughter->GetPdgCode() == 11){
2683 electron = tmpDaughter;
2684 }
2685 else if(tmpDaughter->GetPdgCode() == -11){
2686 positron = tmpDaughter;
2687 }
2688 }
2689 }
2690 if(electron!=NULL && positron!=0){
2691 if(electron->R()<160){
2692 indexOfGammaParticle.push_back(iTracks);
2693 }
2694 }
2695 }
2696 }
2697 }
a0b94e5c 2698
d7d7e825 2699 Int_t nFoundGammas=0;
2700 Int_t nNotFoundGammas=0;
a0b94e5c 2701
d7d7e825 2702 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
2703 for(Int_t i=0;i<numberOfV0s;i++){
2704 fV0Reader->GetV0(i);
a0b94e5c 2705
d7d7e825 2706 if(fV0Reader->HasSameMCMother() == kFALSE){
2707 continue;
2708 }
a0b94e5c 2709
d7d7e825 2710 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2711 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 2712
d7d7e825 2713 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
2714 continue;
2715 }
2716 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
2717 continue;
2718 }
a0b94e5c 2719
d7d7e825 2720 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2721 //TParticle * v0Gamma = fV0Reader->GetMotherMCParticle();
2722 for(UInt_t mcIndex=0;mcIndex<indexOfGammaParticle.size();mcIndex++){
2723 if(negativeMC->GetFirstMother()==indexOfGammaParticle[mcIndex]){
2724 nFoundGammas++;
2725 }
2726 else{
2727 nNotFoundGammas++;
2728 }
2729 }
2730 }
2731 }
d7d7e825 2732}
2733
2734
2735void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){
2736 // see header file for documantation
a0b94e5c 2737
d7d7e825 2738 fESDEvent = fV0Reader->GetESDEvent();
a0b94e5c 2739
2740
2741 TClonesArray * vESDeNegTemp = new TClonesArray("AliESDtrack",0);
6c84d371 2742 TClonesArray * vESDePosTemp = new TClonesArray("AliESDtrack",0);
2743 TClonesArray * vESDxNegTemp = new TClonesArray("AliESDtrack",0);
2744 TClonesArray * vESDxPosTemp = new TClonesArray("AliESDtrack",0);
2745 TClonesArray * vESDeNegNoJPsi = new TClonesArray("AliESDtrack",0);
2746 TClonesArray * vESDePosNoJPsi = new TClonesArray("AliESDtrack",0);
a0b94e5c 2747
6c84d371 2748 /*
2749 vector <AliESDtrack*> vESDeNegTemp(0);
2750 vector <AliESDtrack*> vESDePosTemp(0);
2751 vector <AliESDtrack*> vESDxNegTemp(0);
2752 vector <AliESDtrack*> vESDxPosTemp(0);
2753 vector <AliESDtrack*> vESDeNegNoJPsi(0);
2754 vector <AliESDtrack*> vESDePosNoJPsi(0);
2755 */
a0b94e5c 2756
2757
d7d7e825 2758 fHistograms->FillTable("Table_Electrons",0);//Count number of Events
a0b94e5c 2759
d7d7e825 2760 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
2761 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 2762
d7d7e825 2763 if(!curTrack){
2764 //print warning here
2765 continue;
2766 }
a0b94e5c 2767
d7d7e825 2768 double p[3];if(!curTrack->GetConstrainedPxPyPz(p))continue;
2769 double r[3];curTrack->GetConstrainedXYZ(r);
a0b94e5c 2770
d7d7e825 2771 TVector3 rXYZ(r);
a0b94e5c 2772
d7d7e825 2773 fHistograms->FillTable("Table_Electrons",4);//Count number of ESD tracks
a0b94e5c 2774
d7d7e825 2775 Bool_t flagKink = kTRUE;
2776 Bool_t flagTPCrefit = kTRUE;
2777 Bool_t flagTRDrefit = kTRUE;
2778 Bool_t flagITSrefit = kTRUE;
2779 Bool_t flagTRDout = kTRUE;
2780 Bool_t flagVertex = kTRUE;
a0b94e5c 2781
2782
d7d7e825 2783 //Cuts ---------------------------------------------------------------
a0b94e5c 2784
d7d7e825 2785 if(curTrack->GetKinkIndex(0) > 0){
2786 fHistograms->FillHistogram("Table_Electrons",5);//Count kink
2787 flagKink = kFALSE;
2788 }
a0b94e5c 2789
d7d7e825 2790 ULong_t trkStatus = curTrack->GetStatus();
a0b94e5c 2791
d7d7e825 2792 ULong_t tpcRefit = (trkStatus & AliESDtrack::kTPCrefit);
a0b94e5c 2793
d7d7e825 2794 if(!tpcRefit){
2795 fHistograms->FillHistogram("Table_Electrons",9);//Count not TPCrefit
2796 flagTPCrefit = kFALSE;
2797 }
a0b94e5c 2798
d7d7e825 2799 ULong_t itsRefit = (trkStatus & AliESDtrack::kITSrefit);
2800 if(!itsRefit){
2801 fHistograms->FillHistogram("Table_Electrons",10);//Count not ITSrefit
2802 flagITSrefit = kFALSE;
2803 }
a0b94e5c 2804
d7d7e825 2805 ULong_t trdRefit = (trkStatus & AliESDtrack::kTRDrefit);
a0b94e5c 2806
d7d7e825 2807 if(!trdRefit){
2808 fHistograms->FillHistogram("Table_Electrons",8); //Count not TRDrefit
2809 flagTRDrefit = kFALSE;
2810 }
a0b94e5c 2811
d7d7e825 2812 ULong_t trdOut = (trkStatus & AliESDtrack::kTRDout);
a0b94e5c 2813
d7d7e825 2814 if(!trdOut) {
2815 fHistograms->FillHistogram("Table_Electrons",7); //Count not TRDout
2816 flagTRDout = kFALSE;
2817 }
a0b94e5c 2818
d7d7e825 2819 double nsigmaToVxt = GetSigmaToVertex(curTrack);
a0b94e5c 2820
d7d7e825 2821 if(nsigmaToVxt > 3){
2822 fHistograms->FillHistogram("Table_Electrons",6); //Count Tracks with number of sigmas > 3
2823 flagVertex = kFALSE;
2824 }
a0b94e5c 2825
d7d7e825 2826 if(! (flagKink && flagTPCrefit && flagITSrefit && flagTRDrefit && flagTRDout && flagVertex ) ) continue;
2827 fHistograms->FillHistogram("Table_Electrons",11);//Count Tracks passed Cuts
a0b94e5c 2828
2829
d7d7e825 2830 Stat_t pid, weight;
2831 GetPID(curTrack, pid, weight);
a0b94e5c 2832
d7d7e825 2833 if(pid!=0){
2834 fHistograms->FillHistogram("Table_Electrons",12); //Count Tracks with pid != 0
2835 }
a0b94e5c 2836
d7d7e825 2837 if(pid == 0){
2838 fHistograms->FillHistogram("Table_Electrons",13); //Count Tracks with pid != 0
2839 }
a0b94e5c 2840
2841
2842
2843
a0b94e5c 2844
2845
d7d7e825 2846 TLorentzVector curElec;
2847 curElec.SetXYZM(p[0],p[1],p[2],fElectronMass);
a0b94e5c 2848
2849
113d8432 2850 if(fDoMCTruth){
2851 Int_t labelMC = TMath::Abs(curTrack->GetLabel());
2852 TParticle* curParticle = fStack->Particle(labelMC);
2853 if(curTrack->GetSign() > 0){
2854 if( pid == 0){
2855 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
2856 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
2857 }
2858 else{
2859 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
2860 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
2861 }
2862 }
2863 }
a0b94e5c 2864
2865
d7d7e825 2866 if(curTrack->GetSign() > 0){
a0b94e5c 2867
6c84d371 2868 // vESDxPosTemp.push_back(curTrack);
2869 new((*vESDxPosTemp)[vESDxPosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 2870
d7d7e825 2871 if( pid == 0){
a0b94e5c 2872
d7d7e825 2873 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
2874 fHistograms->FillHistogram("ESD_ElectronPosPt",curElec.Pt());
113d8432 2875 // fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
d7d7e825 2876 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
113d8432 2877 // fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
6c84d371 2878 // vESDePosTemp.push_back(curTrack);
2879 new((*vESDePosTemp)[vESDePosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 2880
d7d7e825 2881 }
a0b94e5c 2882
d7d7e825 2883 }
2884 else {