]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/GammaConv/AliAnalysisTaskGammaConversion.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Cleanup and changing of cuts (Ana)
[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
33#include "AliGammaConversionAODObject.h"
5e55d806 34#include "AliGammaConversionBGHandler.h"
d7d7e825 35
4a6157dc 36class AliCFContainer;
37class AliCFManager;
d7d7e825 38class AliKFVertex;
39class AliAODHandler;
40class AliAODEvent;
41class ALiESDEvent;
42class AliMCEvent;
43class AliMCEventHandler;
44class AliESDInputHandler;
45class AliAnalysisManager;
46class Riostream;
47class TFile;
48class TInterpreter;
49class TSystem;
50class TROOT;
51
52ClassImp(AliAnalysisTaskGammaConversion)
53
54
55AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion():
56AliAnalysisTaskSE(),
57 fV0Reader(NULL),
58 fStack(NULL),
a0b94e5c 59 fMCTruth(NULL), // for CF
4a6157dc 60 fGCMCEvent(NULL), // for CF
d7d7e825 61 fESDEvent(NULL),
62 fOutputContainer(NULL),
a0b94e5c 63 fCFManager(0x0), // for CF
d7d7e825 64 fHistograms(NULL),
b5832f95 65 fTriggerCINT1B(kFALSE),
d7d7e825 66 fDoMCTruth(kFALSE),
67 fDoNeutralMeson(kFALSE),
68 fDoJet(kFALSE),
69 fDoChic(kFALSE),
6c84d371 70 fKFReconstructedGammasTClone(NULL),
71 fCurrentEventPosElectronTClone(NULL),
72 fCurrentEventNegElectronTClone(NULL),
73 fKFReconstructedGammasCutTClone(NULL),
74 fPreviousEventTLVNegElectronTClone(NULL),
75 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 76 fElectronv1(),
77 fElectronv2(),
d7d7e825 78 fElectronMass(-1),
79 fGammaMass(-1),
80 fPi0Mass(-1),
81 fEtaMass(-1),
82 fGammaWidth(-1),
83 fPi0Width(-1),
84 fEtaWidth(-1),
85 fMinOpeningAngleGhostCut(0.),
9640a3d1 86 fEsdTrackCuts(NULL),
d7d7e825 87 fCalculateBackground(kFALSE),
88 fWriteNtuple(kFALSE),
89 fGammaNtuple(NULL),
90 fNeutralMesonNtuple(NULL),
91 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 92 fChargedParticles(NULL),
d7d7e825 93 fChargedParticlesId(),
94 fGammaPtHighest(0.),
95 fMinPtForGammaJet(1.),
96 fMinIsoConeSize(0.2),
97 fMinPtIsoCone(0.7),
98 fMinPtGamChargedCorr(0.5),
99 fMinPtJetCone(0.5),
100 fLeadingChargedIndex(-1),
9640a3d1 101 fLowPtMapping(1.),
102 fHighPtMapping(3.),
1e7846f4 103 fDoCF(kFALSE),
d7d7e825 104 fAODBranch(NULL),
037dc2db 105 fAODBranchName("GammaConv"),
106 fDoNeutralMesonV0MCCheck(kFALSE),
107 fKFReconstructedGammasV0Index()
d7d7e825 108{
109 // Default constructor
5a34881d 110
111 /* Kenneth: the default constructor should not have any define input/output or the call to SetESDtrackCuts
d7d7e825 112 // Common I/O in slot 0
113 DefineInput (0, TChain::Class());
114 DefineOutput(0, TTree::Class());
115
116 // Your private output
117 DefineOutput(1, TList::Class());
a0b94e5c 118
d7d7e825 119 // Define standard ESD track cuts for Gamma-hadron correlation
120 SetESDtrackCuts();
5a34881d 121 */
d7d7e825 122}
123
124AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name):
125 AliAnalysisTaskSE(name),
126 fV0Reader(NULL),
127 fStack(NULL),
a0b94e5c 128 fMCTruth(NULL), // for CF
4a6157dc 129 fGCMCEvent(NULL), // for CF
d7d7e825 130 fESDEvent(NULL),
131 fOutputContainer(0x0),
a0b94e5c 132 fCFManager(0x0), // for CF
d7d7e825 133 fHistograms(NULL),
b5832f95 134 fTriggerCINT1B(kFALSE),
d7d7e825 135 fDoMCTruth(kFALSE),
136 fDoNeutralMeson(kFALSE),
137 fDoJet(kFALSE),
138 fDoChic(kFALSE),
6c84d371 139 fKFReconstructedGammasTClone(NULL),
140 fCurrentEventPosElectronTClone(NULL),
141 fCurrentEventNegElectronTClone(NULL),
142 fKFReconstructedGammasCutTClone(NULL),
143 fPreviousEventTLVNegElectronTClone(NULL),
144 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 145 fElectronv1(),
146 fElectronv2(),
d7d7e825 147 fElectronMass(-1),
148 fGammaMass(-1),
149 fPi0Mass(-1),
150 fEtaMass(-1),
151 fGammaWidth(-1),
152 fPi0Width(-1),
153 fEtaWidth(-1),
154 fMinOpeningAngleGhostCut(0.),
9640a3d1 155 fEsdTrackCuts(NULL),
d7d7e825 156 fCalculateBackground(kFALSE),
157 fWriteNtuple(kFALSE),
158 fGammaNtuple(NULL),
159 fNeutralMesonNtuple(NULL),
160 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 161 fChargedParticles(NULL),
d7d7e825 162 fChargedParticlesId(),
163 fGammaPtHighest(0.),
164 fMinPtForGammaJet(1.),
165 fMinIsoConeSize(0.2),
166 fMinPtIsoCone(0.7),
167 fMinPtGamChargedCorr(0.5),
168 fMinPtJetCone(0.5),
169 fLeadingChargedIndex(-1),
9640a3d1 170 fLowPtMapping(1.),
171 fHighPtMapping(3.),
1e7846f4 172 fDoCF(kFALSE),
d7d7e825 173 fAODBranch(NULL),
037dc2db 174 fAODBranchName("GammaConv"),
175 fDoNeutralMesonV0MCCheck(kFALSE),
176 fKFReconstructedGammasV0Index()
d7d7e825 177{
178 // Common I/O in slot 0
179 DefineInput (0, TChain::Class());
180 DefineOutput(0, TTree::Class());
181
182 // Your private output
183 DefineOutput(1, TList::Class());
a0b94e5c 184 DefineOutput(2, AliCFContainer::Class()); // for CF
185
186
d7d7e825 187 // Define standard ESD track cuts for Gamma-hadron correlation
188 SetESDtrackCuts();
189}
190
191AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion()
192{
193 // Remove all pointers
194
195 if(fOutputContainer){
196 fOutputContainer->Clear() ;
197 delete fOutputContainer ;
198 }
199 if(fHistograms){
200 delete fHistograms;
201 }
202 if(fV0Reader){
203 delete fV0Reader;
204 }
a0b94e5c 205
206 // for CF
207 if(fCFManager){
208 delete fCFManager;
209 }
9640a3d1 210
211 if(fEsdTrackCuts){
212 delete fEsdTrackCuts;
213 }
214
d7d7e825 215 if (fAODBranch) {
216 fAODBranch->Clear();
217 delete fAODBranch ;
218 }
219}
220
221
222void AliAnalysisTaskGammaConversion::Init()
223{
224 // Initialization
4a6157dc 225 // AliLog::SetGlobalLogLevel(AliLog::kError);
d7d7e825 226}
227void AliAnalysisTaskGammaConversion::SetESDtrackCuts()
228{
229 // SetESDtrackCuts
9640a3d1 230 if (fEsdTrackCuts!=NULL){
231 delete fEsdTrackCuts;
232 }
d7d7e825 233 fEsdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
a0b94e5c 234 //standard cuts from:
235 //http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup
ebcfaa7e 236
237 // Cuts used up to 3rd of March
238
239 // fEsdTrackCuts->SetMinNClustersTPC(50);
240// fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
241// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
242// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
243// fEsdTrackCuts->SetMaxNsigmaToVertex(3);
244// fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE);
245
246 //------- To be tested-----------
247
248 Double_t minNClustersTPC = 50;
249 Double_t maxChi2PerClusterTPC = 4.0;
250 Double_t maxDCAtoVertexXY = 2.4; // cm
251 Double_t maxDCAtoVertexZ = 3.2; // cm
252 fEsdTrackCuts->SetRequireSigmaToVertex(kFALSE);
253 fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
254 fEsdTrackCuts->SetRequireITSRefit(kTRUE);
255 // fEsdTrackCuts->SetRequireTPCStandAlone(kTRUE);
256 fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
257 fEsdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
258 fEsdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
259 fEsdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
260 fEsdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
261 fEsdTrackCuts->SetDCAToVertex2D(kTRUE);
262 fEsdTrackCuts->SetEtaRange(-0.9, 0.9);
263 fEsdTrackCuts->SetPtRange(0.15);
264
265
266
267
268 //----- From Jacek 10.03.03 ------------------/
269// minNClustersTPC = 70;
270// maxChi2PerClusterTPC = 4.0;
271// maxDCAtoVertexXY = 2.4; // cm
272// maxDCAtoVertexZ = 3.2; // cm
273
274// esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
275// esdTrackCuts->SetRequireTPCRefit(kFALSE);
276// esdTrackCuts->SetRequireTPCStandAlone(kTRUE);
277// esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
278// esdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
279// esdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
280// esdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
281// esdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
282// esdTrackCuts->SetDCAToVertex2D(kTRUE);
283
284
285
d7d7e825 286 // fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
037dc2db 287 // fV0Reader->SetESDtrackCuts(fEsdTrackCuts);
d7d7e825 288}
289
c00009fb 290void AliAnalysisTaskGammaConversion::UserExec(Option_t */*option*/)
d7d7e825 291{
292 // Execute analysis for current event
d7d7e825 293 ConnectInputData("");
294
295 //Each event needs an empty branch
1e7846f4 296 // fAODBranch->Clear();
297 fAODBranch->Delete();
a0b94e5c 298
6c84d371 299 if(fKFReconstructedGammasTClone == NULL){
300 fKFReconstructedGammasTClone = new TClonesArray("AliKFParticle",0);
301 }
302 if(fCurrentEventPosElectronTClone == NULL){
303 fCurrentEventPosElectronTClone = new TClonesArray("AliESDtrack",0);
304 }
305 if(fCurrentEventNegElectronTClone == NULL){
306 fCurrentEventNegElectronTClone = new TClonesArray("AliESDtrack",0);
307 }
308 if(fKFReconstructedGammasCutTClone == NULL){
309 fKFReconstructedGammasCutTClone = new TClonesArray("AliKFParticle",0);
310 }
311 if(fPreviousEventTLVNegElectronTClone == NULL){
312 fPreviousEventTLVNegElectronTClone = new TClonesArray("TLorentzVector",0);
313 }
314 if(fPreviousEventTLVPosElectronTClone == NULL){
315 fPreviousEventTLVPosElectronTClone = new TClonesArray("TLorentzVector",0);
316 }
317 if(fChargedParticles == NULL){
318 fChargedParticles = new TClonesArray("AliESDtrack",0);
319 }
a0b94e5c 320
6c84d371 321 //clear TClones
1e7846f4 322 fKFReconstructedGammasTClone->Delete();
323 fCurrentEventPosElectronTClone->Delete();
324 fCurrentEventNegElectronTClone->Delete();
325 fKFReconstructedGammasCutTClone->Delete();
326 fPreviousEventTLVNegElectronTClone->Delete();
327 fPreviousEventTLVPosElectronTClone->Delete();
5e55d806 328
d7d7e825 329 //clear vectors
d7d7e825 330 fElectronv1.clear();
331 fElectronv2.clear();
d7d7e825 332
1e7846f4 333 fChargedParticles->Delete();
5e55d806 334
d7d7e825 335 fChargedParticlesId.clear();
037dc2db 336
337 fKFReconstructedGammasV0Index.clear();
a0b94e5c 338
d7d7e825 339 //Clear the data in the v0Reader
5e55d806 340 // fV0Reader->UpdateEventByEventData();
b5832f95 341
342 //Take Only events with proper trigger
5e55d806 343 /*
b5832f95 344 if(fTriggerCINT1B){
345 if(!fV0Reader->GetESDEvent()->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL")) return;
346 }
5e55d806 347 */
a0b94e5c 348
d7d7e825 349 // Process the MC information
350 if(fDoMCTruth){
351 ProcessMCData();
352 }
a0b94e5c 353
d7d7e825 354 //Process the v0 information with no cuts
355 ProcessV0sNoCut();
cb90a330 356
d7d7e825 357 // Process the v0 information
358 ProcessV0s();
cb90a330 359
5e55d806 360
d7d7e825 361 //Fill Gamma AOD
362 FillAODWithConversionGammas() ;
a0b94e5c 363
d7d7e825 364 // Process reconstructed gammas
365 if(fDoNeutralMeson == kTRUE){
a0b94e5c 366 ProcessGammasForNeutralMesonAnalysis();
d7d7e825 367 }
77880bd8 368
369 if(fDoMCTruth == kTRUE){
370 CheckV0Efficiency();
371 }
d7d7e825 372 //Process reconstructed gammas electrons for Chi_c Analysis
6c84d371 373 if(fDoChic == kTRUE){
d7d7e825 374 ProcessGammaElectronsForChicAnalysis();
375 }
376 // Process reconstructed gammas for gamma Jet/hadron correlations
377 if(fDoJet == kTRUE){
378 ProcessGammasForGammaJetAnalysis();
379 }
a0b94e5c 380
5e55d806 381 //calculate background if flag is set
382 if(fCalculateBackground){
383 CalculateBackground();
384 }
385
386 //Clear the data in the v0Reader
387 fV0Reader->UpdateEventByEventData();
388
d7d7e825 389 PostData(1, fOutputContainer);
a0b94e5c 390 PostData(2, fCFManager->GetParticleContainer()); // for CF
d7d7e825 391
392}
393
8a685cf3 394void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t *option){
d7d7e825 395 // see header file for documentation
8a685cf3 396
397 AliAnalysisTaskSE::ConnectInputData(option);
398
d7d7e825 399 if(fV0Reader == NULL){
400 // Write warning here cuts and so on are default if this ever happens
401 }
402 fV0Reader->Initialize();
61374d97 403 fDoMCTruth = fV0Reader->GetDoMCTruth();
d7d7e825 404}
405
406
407
408void AliAnalysisTaskGammaConversion::ProcessMCData(){
409 // see header file for documentation
410
411 fStack = fV0Reader->GetMCStack();
a0b94e5c 412 fMCTruth = fV0Reader->GetMCTruth(); // for CF
4a6157dc 413 fGCMCEvent = fV0Reader->GetMCEvent(); // for CF
a0b94e5c 414
415
416 // for CF
1e7846f4 417 Double_t containerInput[3];
418 if(fDoCF){
419 if(!fGCMCEvent) cout << "NO MC INFO FOUND" << endl;
420 fCFManager->SetEventInfo(fGCMCEvent);
421 }
a0b94e5c 422 // end for CF
423
424
d7d7e825 425 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
426 return; // aborts if the primary vertex does not have contributors.
427 }
a0b94e5c 428
d7d7e825 429 for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
430 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
431
432 if (!particle) {
433 //print warning here
434 continue;
435 }
a0b94e5c 436
d7d7e825 437 ///////////////////////Begin Chic Analysis/////////////////////////////
a0b94e5c 438
a68437fb 439 if(particle->GetPdgCode() == 443){//Is JPsi
d7d7e825 440 if(particle->GetNDaughters()==2){
441 if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 &&
442 TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){
443 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
444 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
445 if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9)
446 fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance
a0b94e5c 447
d7d7e825 448 if( TMath::Abs(daug0->Eta()) < 0.9){
449 if(daug0->GetPdgCode() == -11)
450 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
451 else
452 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
a0b94e5c 453
d7d7e825 454 }
455 if(TMath::Abs(daug1->Eta()) < 0.9){
456 if(daug1->GetPdgCode() == -11)
457 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
458 else
459 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
460 }
461 }
462 }
463 }
464 // const int CHI_C0 = 10441;
465 // const int CHI_C1 = 20443;
466 // const int CHI_C2 = 445
467 if(particle->GetPdgCode() == 22){//gamma from JPsi
468 if(particle->GetMother(0) > -1){
469 if(fStack->Particle(particle->GetMother(0))->GetPdgCode() == 10441 ||
470 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 20443 ||
471 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 445){
472 if(TMath::Abs(particle->Eta()) < 1.2)
473 fHistograms->FillTable("Table_Electrons",17);// gamma from chic inside accptance
474 }
475 }
476 }
477 if(particle->GetPdgCode() == 10441 || particle->GetPdgCode() == 20443 || particle->GetPdgCode() == 445){
478 if( particle->GetNDaughters() == 2){
479 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
480 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
a0b94e5c 481
d7d7e825 482 if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){
483 if( daug0->GetPdgCode() == 443){
484 TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter());
485 TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter());
486 if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
487 fHistograms->FillTable("Table_Electrons",18);
a0b94e5c 488
d7d7e825 489 }//if
490 else if (daug1->GetPdgCode() == 443){
491 TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter());
492 TParticle* daugE1 = fStack->Particle(daug1->GetLastDaughter());
493 if( TMath::Abs(daug0->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
494 fHistograms->FillTable("Table_Electrons",18);
495 }//else if
496 }//gamma o Jpsi
497 }//GetNDaughters
498 }
a0b94e5c 499
500
d7d7e825 501 /////////////////////End Chic Analysis////////////////////////////
a0b94e5c 502
503
d7d7e825 504 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
a0b94e5c 505
d7d7e825 506 if(particle->R()>fV0Reader->GetMaxRCut()) continue; // cuts on distance from collision point
507
508 Double_t tmpPhi=particle->Phi();
509
510 if(particle->Phi()> TMath::Pi()){
511 tmpPhi = particle->Phi()-(2*TMath::Pi());
512 }
513
514 Double_t rapidity;
515 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
516 rapidity=0;
517 }
518 else{
519 rapidity = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())));
520 }
521
522 //process the gammas
523 if (particle->GetPdgCode() == 22){
a607c218 524
525
d7d7e825 526 if(particle->GetMother(0) >-1 && fStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
527 continue; // no photon as mothers!
528 }
a607c218 529
d7d7e825 530 if(particle->GetMother(0) >= fStack->GetNprimary()){
531 continue; // the gamma has a mother, and it is not a primary particle
532 }
a0b94e5c 533
d7d7e825 534 fHistograms->FillHistogram("MC_allGamma_Energy", particle->Energy());
535 fHistograms->FillHistogram("MC_allGamma_Pt", particle->Pt());
536 fHistograms->FillHistogram("MC_allGamma_Eta", particle->Eta());
537 fHistograms->FillHistogram("MC_allGamma_Phi", tmpPhi);
538 fHistograms->FillHistogram("MC_allGamma_Rapid", rapidity);
539
a0b94e5c 540 // for CF
1e7846f4 541 if(fDoCF){
542 containerInput[0] = particle->Pt();
543 containerInput[1] = particle->Eta();
544 if(particle->GetMother(0) >=0){
545 containerInput[2] = fStack->Particle(particle->GetMother(0))->GetMass();
546 }
547 else{
548 containerInput[2]=-1;
549 }
a607c218 550
1e7846f4 551 fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated); // generated gamma
552 }
d7d7e825 553 if(particle->GetMother(0) < 0){ // direct gamma
554 fHistograms->FillHistogram("MC_allDirectGamma_Energy",particle->Energy());
555 fHistograms->FillHistogram("MC_allDirectGamma_Pt", particle->Pt());
556 fHistograms->FillHistogram("MC_allDirectGamma_Eta", particle->Eta());
557 fHistograms->FillHistogram("MC_allDirectGamma_Phi", tmpPhi);
558 fHistograms->FillHistogram("MC_allDirectGamma_Rapid", rapidity);
559 }
560
d7d7e825 561 // looking for conversion (electron + positron from pairbuilding (= 5) )
562 TParticle* ePos = NULL;
563 TParticle* eNeg = NULL;
564
565 if(particle->GetNDaughters() >= 2){
566 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
567 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
568 if(tmpDaughter->GetUniqueID() == 5){
569 if(tmpDaughter->GetPdgCode() == 11){
570 eNeg = tmpDaughter;
571 }
572 else if(tmpDaughter->GetPdgCode() == -11){
573 ePos = tmpDaughter;
574 }
575 }
576 }
577 }
578
579
580 if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
581 continue;
582 }
583
584
585 Double_t ePosPhi = ePos->Phi();
586 if(ePos->Phi()> TMath::Pi()) ePosPhi = ePos->Phi()-(2*TMath::Pi());
587
588 Double_t eNegPhi = eNeg->Phi();
589 if(eNeg->Phi()> TMath::Pi()) eNegPhi = eNeg->Phi()-(2*TMath::Pi());
590
591
592 if(ePos->Pt()<fV0Reader->GetPtCut() || eNeg->Pt()<fV0Reader->GetPtCut()){
593 continue; // no reconstruction below the Pt cut
594 }
a0b94e5c 595
d7d7e825 596 if(TMath::Abs(ePos->Eta())> fV0Reader->GetEtaCut() || TMath::Abs(eNeg->Eta())> fV0Reader->GetEtaCut()){
597 continue;
598 }
a0b94e5c 599
d7d7e825 600 if(ePos->R()>fV0Reader->GetMaxRCut()){
601 continue; // cuts on distance from collision point
602 }
a0b94e5c 603
604 if(TMath::Abs(ePos->Vz()) > fV0Reader->GetMaxZCut()){
d7d7e825 605 continue; // outside material
606 }
607
a0b94e5c 608
d7d7e825 609 if((TMath::Abs(ePos->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() > ePos->R()){
610 continue; // line cut to exclude regions where we do not reconstruct
611 }
612
a0b94e5c 613
614 // for CF
1e7846f4 615 if(fDoCF){
616 fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructable); // reconstructable gamma
617 }
d7d7e825 618 fHistograms->FillHistogram("MC_ConvGamma_Energy", particle->Energy());
619 fHistograms->FillHistogram("MC_ConvGamma_Pt", particle->Pt());
620 fHistograms->FillHistogram("MC_ConvGamma_Eta", particle->Eta());
621 fHistograms->FillHistogram("MC_ConvGamma_Phi", tmpPhi);
622 fHistograms->FillHistogram("MC_ConvGamma_Rapid", rapidity);
623 fHistograms->FillHistogram("MC_ConvGamma_Pt_Eta", particle->Pt(),particle->Eta());
624
625 fHistograms->FillHistogram("MC_E_Energy", eNeg->Energy());
626 fHistograms->FillHistogram("MC_E_Pt", eNeg->Pt());
627 fHistograms->FillHistogram("MC_E_Eta", eNeg->Eta());
628 fHistograms->FillHistogram("MC_E_Phi", eNegPhi);
629
630 fHistograms->FillHistogram("MC_P_Energy", ePos->Energy());
631 fHistograms->FillHistogram("MC_P_Pt", ePos->Pt());
632 fHistograms->FillHistogram("MC_P_Eta", ePos->Eta());
633 fHistograms->FillHistogram("MC_P_Phi", ePosPhi);
634
635
d7d7e825 636 // begin Mapping
637 Int_t rBin = fHistograms->GetRBin(ePos->R());
9640a3d1 638 Int_t zBin = fHistograms->GetZBin(ePos->Vz());
d7d7e825 639 Int_t phiBin = fHistograms->GetPhiBin(particle->Phi());
9640a3d1 640
641 TVector3 vtxPos(ePos->Vx(),ePos->Vy(),ePos->Vz());
642
d7d7e825 643 TString nameMCMappingPhiR="";
e158cbc3 644 nameMCMappingPhiR.Form("MC_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 645 // fHistograms->FillHistogram(nameMCMappingPhiR, ePos->Vz(), particle->Eta());
d7d7e825 646
647 TString nameMCMappingPhi="";
e158cbc3 648 nameMCMappingPhi.Form("MC_Conversion_Mapping_Phi%02d",phiBin);
9640a3d1 649 // fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());
ebcfaa7e 650 //fHistograms->FillHistogram(nameMCMappingPhi, ePos->Vz(), particle->Eta());
d7d7e825 651
652 TString nameMCMappingR="";
e158cbc3 653 nameMCMappingR.Form("MC_Conversion_Mapping_R%02d",rBin);
9640a3d1 654 // fHistograms->FillHistogram(nameMCMappingR, particle->Eta());
ebcfaa7e 655 //fHistograms->FillHistogram(nameMCMappingR,ePos->Vz(), particle->Eta());
d7d7e825 656
657 TString nameMCMappingPhiInR="";
e158cbc3 658 nameMCMappingPhiInR.Form("MC_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 659 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
660 fHistograms->FillHistogram(nameMCMappingPhiInR, vtxPos.Phi());
661
662 TString nameMCMappingZInR="";
663 nameMCMappingZInR.Form("MC_Conversion_Mapping_Z_in_R_%02d",rBin);
664 fHistograms->FillHistogram(nameMCMappingZInR,ePos->Vz() );
665
666
667 TString nameMCMappingPhiInZ="";
668 nameMCMappingPhiInZ.Form("MC_Conversion_Mapping_Phi_in_Z_%02d",zBin);
669 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
670 fHistograms->FillHistogram(nameMCMappingPhiInZ, vtxPos.Phi());
671
672 TString nameMCMappingRInZ="";
673 nameMCMappingRInZ.Form("MC_Conversion_Mapping_R_in_Z_%02d",zBin);
674 fHistograms->FillHistogram(nameMCMappingRInZ,ePos->R() );
675
676 if(particle->Pt() > fLowPtMapping && particle->Pt()< fHighPtMapping){
677 TString nameMCMappingMidPtPhiInR="";
678 nameMCMappingMidPtPhiInR.Form("MC_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
679 fHistograms->FillHistogram(nameMCMappingMidPtPhiInR, vtxPos.Phi());
680
681 TString nameMCMappingMidPtZInR="";
682 nameMCMappingMidPtZInR.Form("MC_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
683 fHistograms->FillHistogram(nameMCMappingMidPtZInR,ePos->Vz() );
684
685
686 TString nameMCMappingMidPtPhiInZ="";
687 nameMCMappingMidPtPhiInZ.Form("MC_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
688 fHistograms->FillHistogram(nameMCMappingMidPtPhiInZ, vtxPos.Phi());
689
690 TString nameMCMappingMidPtRInZ="";
691 nameMCMappingMidPtRInZ.Form("MC_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
692 fHistograms->FillHistogram(nameMCMappingMidPtRInZ,ePos->R() );
693
694 }
695
d7d7e825 696 //end mapping
697
698 fHistograms->FillHistogram("MC_Conversion_R",ePos->R());
699 fHistograms->FillHistogram("MC_Conversion_ZR",ePos->Vz(),ePos->R());
700 fHistograms->FillHistogram("MC_Conversion_XY",ePos->Vx(),ePos->Vy());
701 fHistograms->FillHistogram("MC_Conversion_OpeningAngle",GetMCOpeningAngle(ePos, eNeg));
9640a3d1 702 fHistograms->FillHistogram("MC_ConvGamma_E_AsymmetryP",particle->P(),eNeg->P()/particle->P());
703 fHistograms->FillHistogram("MC_ConvGamma_P_AsymmetryP",particle->P(),ePos->P()/particle->P());
704
705
d7d7e825 706 if(particle->GetMother(0) < 0){ // no mother = direct gamma, still inside converted
707 fHistograms->FillHistogram("MC_ConvDirectGamma_Energy",particle->Energy());
708 fHistograms->FillHistogram("MC_ConvDirectGamma_Pt", particle->Pt());
709 fHistograms->FillHistogram("MC_ConvDirectGamma_Eta", particle->Eta());
710 fHistograms->FillHistogram("MC_ConvDirectGamma_Phi", tmpPhi);
711 fHistograms->FillHistogram("MC_ConvDirectGamma_Rapid", rapidity);
712
713 } // end direct gamma
714 else{ // mother exits
6c84d371 715 /* if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0
a0b94e5c 716 fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S
717 fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2
718 ){
719 fMCGammaChic.push_back(particle);
720 }
6c84d371 721 */
d7d7e825 722 } // end if mother exits
723 } // end if particle is a photon
724
725
726
727 // process motherparticles (2 gammas as daughters)
728 // the motherparticle had already to pass the R and the eta cut, but no line cut.
729 // the line cut is just valid for the conversions!
730
731 if(particle->GetNDaughters() == 2){
732
733 TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());
734 TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());
735
736 if(daughter0->GetPdgCode() != 22 || daughter1->GetPdgCode() != 22) continue; //check for gamma gamma daughters
737
738 // Check the acceptance for both gammas
01b7fdcc 739 Bool_t gammaEtaCut = kTRUE;
740 if(TMath::Abs(daughter0->Eta()) > fV0Reader->GetEtaCut() || TMath::Abs(daughter1->Eta()) > fV0Reader->GetEtaCut() ) gammaEtaCut = kFALSE;
741 Bool_t gammaRCut = kTRUE;
742 if(daughter0->R() > fV0Reader->GetMaxRCut() || daughter1->R() > fV0Reader->GetMaxRCut() ) gammaRCut = kFALSE;
d7d7e825 743
744
745
746 // check for conversions now -> have to pass eta, R and line cut!
747 Bool_t daughter0Electron = kFALSE;
748 Bool_t daughter0Positron = kFALSE;
749 Bool_t daughter1Electron = kFALSE;
750 Bool_t daughter1Positron = kFALSE;
751
752 if(daughter0->GetNDaughters() >= 2){ // first gamma
753 for(Int_t TrackIndex=daughter0->GetFirstDaughter();TrackIndex<=daughter0->GetLastDaughter();TrackIndex++){
754 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
755 if(tmpDaughter->GetUniqueID() == 5){
756 if(tmpDaughter->GetPdgCode() == 11){
757 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
758 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
759 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
760 daughter0Electron = kTRUE;
761 }
762 }
763
764 }
765 }
766 else if(tmpDaughter->GetPdgCode() == -11){
767 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
768 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
769 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
770 daughter0Positron = kTRUE;
771 }
a68437fb 772 }
773 }
d7d7e825 774 }
775 }
776 }
777 }
778
779
780 if(daughter1->GetNDaughters() >= 2){ // second gamma
781 for(Int_t TrackIndex=daughter1->GetFirstDaughter();TrackIndex<=daughter1->GetLastDaughter();TrackIndex++){
782 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
783 if(tmpDaughter->GetUniqueID() == 5){
784 if(tmpDaughter->GetPdgCode() == 11){
785 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
786 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
787 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
788 daughter1Electron = kTRUE;
789 }
790 }
791
792 }
793 }
794 else if(tmpDaughter->GetPdgCode() == -11){
795 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
796 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
797 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
798 daughter1Positron = kTRUE;
799 }
800 }
801
802 }
803
804 }
805 }
806 }
807 }
808
a0b94e5c 809
d7d7e825 810 if(particle->GetPdgCode()==111){ //Pi0
811 if( iTracks >= fStack->GetNprimary()){
812 fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());
813 fHistograms->FillHistogram("MC_Pi0_Secondaries_Rapid", rapidity);
814 fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);
815 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt", particle->Pt());
816 fHistograms->FillHistogram("MC_Pi0_Secondaries_Energy", particle->Energy());
817 fHistograms->FillHistogram("MC_Pi0_Secondaries_R", particle->R());
818 fHistograms->FillHistogram("MC_Pi0_Secondaries_ZR", particle->Vz(),particle->R());
819 fHistograms->FillHistogram("MC_Pi0_Secondaries_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
820 fHistograms->FillHistogram("MC_Pi0_Secondaries_XY", particle->Vx(),particle->Vy());//only fill from one daughter to avoid multiple filling
821
01b7fdcc 822 if(gammaEtaCut && gammaRCut){
d7d7e825 823 //if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
824 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
825 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
826 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
827 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
828 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
829 }
830 }
831 }
832 else{
833 fHistograms->FillHistogram("MC_Pi0_Eta", particle->Eta());
834 fHistograms->FillHistogram("MC_Pi0_Rapid", rapidity);
835 fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);
836 fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());
837 fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());
838 fHistograms->FillHistogram("MC_Pi0_R", particle->R());
839 fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());
840 fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
841 fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
a0b94e5c 842
01b7fdcc 843 if(gammaEtaCut && gammaRCut){
d7d7e825 844 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
845 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
846 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
847 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
848 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
849 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
850 fHistograms->FillHistogram("MC_Pi0_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
851 }
852 }
853 }
854 }
855
856 if(particle->GetPdgCode()==221){ //Eta
857 fHistograms->FillHistogram("MC_Eta_Eta", particle->Eta());
858 fHistograms->FillHistogram("MC_Eta_Rapid", rapidity);
859 fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);
860 fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());
861 fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());
862 fHistograms->FillHistogram("MC_Eta_R", particle->R());
863 fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());
864 fHistograms->FillHistogram("MC_Eta_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
865 fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
866
01b7fdcc 867 if(gammaEtaCut && gammaRCut){
a0b94e5c 868 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 869 fHistograms->FillHistogram("MC_Eta_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
870 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
871 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
872 fHistograms->FillHistogram("MC_Eta_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
873 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
874 fHistograms->FillHistogram("MC_Eta_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
875 }
876
877 }
878
879 }
880
881 // all motherparticles with 2 gammas as daughters
882 fHistograms->FillHistogram("MC_Mother_R", particle->R());
883 fHistograms->FillHistogram("MC_Mother_ZR", particle->Vz(),particle->R());
884 fHistograms->FillHistogram("MC_Mother_XY", particle->Vx(),particle->Vy());
885 fHistograms->FillHistogram("MC_Mother_Mass", particle->GetCalcMass());
886 fHistograms->FillHistogram("MC_Mother_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
887 fHistograms->FillHistogram("MC_Mother_Energy", particle->Energy());
888 fHistograms->FillHistogram("MC_Mother_Pt", particle->Pt());
889 fHistograms->FillHistogram("MC_Mother_Eta", particle->Eta());
890 fHistograms->FillHistogram("MC_Mother_Rapid", rapidity);
891 fHistograms->FillHistogram("MC_Mother_Phi",tmpPhi);
892 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt",particle->GetMass(),particle->Pt());
a0b94e5c 893
01b7fdcc 894 if(gammaEtaCut && gammaRCut){
a0b94e5c 895 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 896 fHistograms->FillHistogram("MC_Mother_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
897 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
898 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_withinAcceptance",particle->GetMass(),particle->Pt());
899 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
900 fHistograms->FillHistogram("MC_Mother_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
901 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
902 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_ConvGamma_withinAcceptance",particle->GetMass(),particle->Pt());
a0b94e5c 903
d7d7e825 904 }
905
906
907 } // end passed R and eta cut
a0b94e5c 908
d7d7e825 909 } // end if(particle->GetNDaughters() == 2)
910
911 }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++)
a0b94e5c 912
d7d7e825 913} // end ProcessMCData
914
915
916
917void AliAnalysisTaskGammaConversion::FillNtuple(){
918 //Fills the ntuple with the different values
a0b94e5c 919
d7d7e825 920 if(fGammaNtuple == NULL){
921 return;
922 }
923 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
924 for(Int_t i=0;i<numberOfV0s;i++){
925 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};
926 AliESDv0 * cV0 = fV0Reader->GetV0(i);
927 Double_t negPID=0;
928 Double_t posPID=0;
929 fV0Reader->GetPIDProbability(negPID,posPID);
930 values[0]=cV0->GetOnFlyStatus();
931 values[1]=fV0Reader->CheckForPrimaryVertex();
932 values[2]=negPID;
933 values[3]=posPID;
934 values[4]=fV0Reader->GetX();
935 values[5]=fV0Reader->GetY();
936 values[6]=fV0Reader->GetZ();
937 values[7]=fV0Reader->GetXYRadius();
938 values[8]=fV0Reader->GetMotherCandidateNDF();
939 values[9]=fV0Reader->GetMotherCandidateChi2();
940 values[10]=fV0Reader->GetMotherCandidateEnergy();
941 values[11]=fV0Reader->GetMotherCandidateEta();
942 values[12]=fV0Reader->GetMotherCandidatePt();
943 values[13]=fV0Reader->GetMotherCandidateMass();
944 values[14]=fV0Reader->GetMotherCandidateWidth();
945 // values[15]=fV0Reader->GetMotherMCParticle()->Pt(); MOVED TO THE END, HAS TO BE CALLED AFTER HasSameMother NB: still has the same entry in the array
946 values[16]=fV0Reader->GetOpeningAngle();
947 values[17]=fV0Reader->GetNegativeTrackEnergy();
948 values[18]=fV0Reader->GetNegativeTrackPt();
949 values[19]=fV0Reader->GetNegativeTrackEta();
950 values[20]=fV0Reader->GetNegativeTrackPhi();
951 values[21]=fV0Reader->GetPositiveTrackEnergy();
952 values[22]=fV0Reader->GetPositiveTrackPt();
953 values[23]=fV0Reader->GetPositiveTrackEta();
954 values[24]=fV0Reader->GetPositiveTrackPhi();
955 values[25]=fV0Reader->HasSameMCMother();
956 if(values[25] != 0){
957 values[26]=fV0Reader->GetMotherMCParticlePDGCode();
958 values[15]=fV0Reader->GetMotherMCParticle()->Pt();
959 }
960 fTotalNumberOfAddedNtupleEntries++;
961 fGammaNtuple->Fill(values);
962 }
963 fV0Reader->ResetV0IndexNumber();
964
965}
966
967void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){
968 // Process all the V0's without applying any cuts to it
a0b94e5c 969
d7d7e825 970 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
971 for(Int_t i=0;i<numberOfV0s;i++){
972 /*AliESDv0 * cV0 = */fV0Reader->GetV0(i);
a0b94e5c 973
d7d7e825 974 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
cb90a330 975 continue;
d7d7e825 976 }
9640a3d1 977
77880bd8 978 // if( !fV0Reader->GetV0(i)->GetOnFlyStatus()){
979 if( !fV0Reader->CheckV0FinderStatus(i)){
cb90a330 980 continue;
9640a3d1 981 }
982
cb90a330 983
984 if( !((fV0Reader->GetNegativeESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ||
985 !((fV0Reader->GetPositiveESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ){
986 continue;
9640a3d1 987 }
cb90a330 988
ebcfaa7e 989 if( fV0Reader->GetNegativeESDTrack()->GetSign()== fV0Reader->GetPositiveESDTrack()->GetSign()){
990 continue;
991 }
992
993 if( fV0Reader->GetNegativeESDTrack()->GetKinkIndex(0) > 0 ||
994 fV0Reader->GetPositiveESDTrack()->GetKinkIndex(0) > 0) {
995 continue;
996 }
997
d7d7e825 998 if(fDoMCTruth){
a0b94e5c 999
d7d7e825 1000 if(fV0Reader->HasSameMCMother() == kFALSE){
1001 continue;
1002 }
a0b94e5c 1003
d7d7e825 1004 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1005 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 1006
d7d7e825 1007 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1008 continue;
1009 }
26923b22 1010 if(negativeMC->GetPdgCode() == positiveMC->GetPdgCode()){
d7d7e825 1011 continue;
1012 }
a68437fb 1013
26923b22 1014 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){ // id 5 is conversion
a68437fb 1015 continue;
1016 }
a0b94e5c 1017
d7d7e825 1018 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a0b94e5c 1019
d7d7e825 1020 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1021 fHistograms->FillHistogram("ESD_NoCutConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1022 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1023 fHistograms->FillHistogram("ESD_NoCutConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1024 fHistograms->FillHistogram("ESD_NoCutConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1025 fHistograms->FillHistogram("ESD_NoCutConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1026 fHistograms->FillHistogram("ESD_NoCutConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1027 fHistograms->FillHistogram("ESD_NoCutConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1028 fHistograms->FillHistogram("ESD_NoCutConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1029 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
a0b94e5c 1030
d7d7e825 1031 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1032 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
1033
1034 fHistograms->FillHistogram("ESD_NoCutConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1035 fHistograms->FillHistogram("ESD_NoCutConversion_R", fV0Reader->GetXYRadius());
1036 fHistograms->FillHistogram("ESD_NoCutConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1037 fHistograms->FillHistogram("ESD_NoCutConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1038 fHistograms->FillHistogram("ESD_NoCutConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1039 fHistograms->FillHistogram("ESD_NoCutConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1040 fHistograms->FillHistogram("ESD_NoCutConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1041 fHistograms->FillHistogram("ESD_NoCutConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1042
1043 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1044 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1045 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1046 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1047
d7d7e825 1048 //store MCTruth properties
1049 fHistograms->FillHistogram("ESD_NoCutConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1050 fHistograms->FillHistogram("ESD_NoCutConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1051 fHistograms->FillHistogram("ESD_NoCutConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
d7d7e825 1052 }
1053 }
1054 }
1055 fV0Reader->ResetV0IndexNumber();
1056}
1057
1058void AliAnalysisTaskGammaConversion::ProcessV0s(){
1059 // see header file for documentation
037dc2db 1060
1061
d7d7e825 1062 if(fWriteNtuple == kTRUE){
1063 FillNtuple();
1064 }
1065
1066 Int_t nSurvivingV0s=0;
1067 while(fV0Reader->NextV0()){
1068 nSurvivingV0s++;
1069
1070
1071 //-------------------------- filling v0 information -------------------------------------
1072 fHistograms->FillHistogram("ESD_Conversion_R", fV0Reader->GetXYRadius());
1073 fHistograms->FillHistogram("ESD_Conversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1074 fHistograms->FillHistogram("ESD_Conversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1075 fHistograms->FillHistogram("ESD_Conversion_OpeningAngle", fV0Reader->GetOpeningAngle());
1076
1077 fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());
1078 fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());
1079 fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());
1080 fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());
ebcfaa7e 1081 fHistograms->FillHistogram("ESD_E_nTPCClusters", fV0Reader->GetNegativeTracknTPCClusters());
1082 fHistograms->FillHistogram("ESD_E_nITSClusters", fV0Reader->GetNegativeTracknITSClusters());
d7d7e825 1083
1084 fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());
1085 fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());
1086 fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());
1087 fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());
ebcfaa7e 1088 fHistograms->FillHistogram("ESD_P_nTPCClusters", fV0Reader->GetPositiveTracknTPCClusters());
1089 fHistograms->FillHistogram("ESD_P_nITSClusters", fV0Reader->GetPositiveTracknITSClusters());
d7d7e825 1090
1091 fHistograms->FillHistogram("ESD_ConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1092 fHistograms->FillHistogram("ESD_ConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1093 fHistograms->FillHistogram("ESD_ConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1094 fHistograms->FillHistogram("ESD_ConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1095 fHistograms->FillHistogram("ESD_ConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1096 fHistograms->FillHistogram("ESD_ConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1097 fHistograms->FillHistogram("ESD_ConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1098 fHistograms->FillHistogram("ESD_ConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1099 fHistograms->FillHistogram("ESD_ConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1100 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1101
1102 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1103 fHistograms->FillHistogram("ESD_ConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
9640a3d1 1104
1105 fHistograms->FillHistogram("ESD_ConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1106 fHistograms->FillHistogram("ESD_ConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1107 fHistograms->FillHistogram("ESD_ConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1108 fHistograms->FillHistogram("ESD_ConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1109
1110 fHistograms->FillHistogram("ESD_ConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1111 fHistograms->FillHistogram("ESD_ConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1112 fHistograms->FillHistogram("ESD_ConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1113 fHistograms->FillHistogram("ESD_ConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1114
037dc2db 1115
d7d7e825 1116 // begin mapping
1117 Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());
9640a3d1 1118 Int_t zBin = fHistograms->GetZBin(fV0Reader->GetZ());
d7d7e825 1119 Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());
9640a3d1 1120 TVector3 vtxConv(fV0Reader->GetX(),fV0Reader->GetY(), fV0Reader->GetZ());
1121
ebcfaa7e 1122 // Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();
d7d7e825 1123
1124 TString nameESDMappingPhiR="";
e158cbc3 1125 nameESDMappingPhiR.Form("ESD_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 1126 //fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1127
1128 TString nameESDMappingPhi="";
e158cbc3 1129 nameESDMappingPhi.Form("ESD_Conversion_Mapping_Phi%02d",phiBin);
ebcfaa7e 1130 //fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1131
1132 TString nameESDMappingR="";
e158cbc3 1133 nameESDMappingR.Form("ESD_Conversion_Mapping_R%02d",rBin);
ebcfaa7e 1134 //fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1135
1136 TString nameESDMappingPhiInR="";
e158cbc3 1137 nameESDMappingPhiInR.Form("ESD_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 1138 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1139 fHistograms->FillHistogram(nameESDMappingPhiInR, vtxConv.Phi());
1140
1141 TString nameESDMappingZInR="";
1142 nameESDMappingZInR.Form("ESD_Conversion_Mapping_Z_in_R_%02d",rBin);
1143 fHistograms->FillHistogram(nameESDMappingZInR, fV0Reader->GetZ());
1144
1145 TString nameESDMappingPhiInZ="";
1146 nameESDMappingPhiInZ.Form("ESD_Conversion_Mapping_Phi_in_Z_%02d",zBin);
1147 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1148 fHistograms->FillHistogram(nameESDMappingPhiInZ, vtxConv.Phi());
1149
1150 TString nameESDMappingRInZ="";
1151 nameESDMappingRInZ.Form("ESD_Conversion_Mapping_R_in_Z_%02d",zBin);
1152 fHistograms->FillHistogram(nameESDMappingRInZ, fV0Reader->GetXYRadius());
1153
1154 if(fV0Reader->GetMotherCandidatePt() > fLowPtMapping && fV0Reader->GetMotherCandidatePt()< fHighPtMapping){
1155 TString nameESDMappingMidPtPhiInR="";
1156 nameESDMappingMidPtPhiInR.Form("ESD_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
1157 fHistograms->FillHistogram(nameESDMappingMidPtPhiInR, vtxConv.Phi());
1158
1159 TString nameESDMappingMidPtZInR="";
1160 nameESDMappingMidPtZInR.Form("ESD_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
1161 fHistograms->FillHistogram(nameESDMappingMidPtZInR, fV0Reader->GetZ());
1162
1163 TString nameESDMappingMidPtPhiInZ="";
1164 nameESDMappingMidPtPhiInZ.Form("ESD_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
1165 fHistograms->FillHistogram(nameESDMappingMidPtPhiInZ, vtxConv.Phi());
1166
1167 TString nameESDMappingMidPtRInZ="";
1168 nameESDMappingMidPtRInZ.Form("ESD_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
1169 fHistograms->FillHistogram(nameESDMappingMidPtRInZ, fV0Reader->GetXYRadius());
1170
1171 }
1172
1173
d7d7e825 1174 // end mapping
1175
6c84d371 1176 new((*fKFReconstructedGammasTClone)[fKFReconstructedGammasTClone->GetEntriesFast()]) AliKFParticle(*fV0Reader->GetMotherCandidateKFCombination());
037dc2db 1177 fKFReconstructedGammasV0Index.push_back(fV0Reader->GetCurrentV0IndexNumber()-1);
6c84d371 1178 // fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());
d7d7e825 1179 fElectronv1.push_back(fV0Reader->GetCurrentV0()->GetPindex());
1180 fElectronv2.push_back(fV0Reader->GetCurrentV0()->GetNindex());
a0b94e5c 1181
d7d7e825 1182 //----------------------------------- checking for "real" conversions (MC match) --------------------------------------
1183 if(fDoMCTruth){
1184
a68437fb 1185 if(fV0Reader->HasSameMCMother() == kFALSE){
1186 continue;
1187 }
d7d7e825 1188 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
a68437fb 1189 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1190
1191 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1192 continue;
1193 }
1194
1195 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
1196 continue;
1197 }
037dc2db 1198 if(negativeMC->GetUniqueID() == 4 && positiveMC->GetUniqueID() ==4){// fill r distribution for Dalitz decays
1199 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 111){ //pi0
1200 fHistograms->FillHistogram("ESD_TrueDalitzContamination_R", fV0Reader->GetXYRadius());
1201 }
1202 }
a68437fb 1203
1204 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){// check if the daughters come from a conversion
1205 continue;
1206 }
a0b94e5c 1207
d7d7e825 1208 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a68437fb 1209
26923b22 1210 if(fDoCF){
1211 Double_t containerInput[3];
1212 containerInput[0] = fV0Reader->GetMotherCandidatePt();
1213 containerInput[1] = fV0Reader->GetMotherCandidateEta();
1214 containerInput[2] = fV0Reader->GetMotherCandidateMass();
1215 fCFManager->GetParticleContainer()->Fill(containerInput,kStepTrueGamma); // for CF
1216 }
1217
d7d7e825 1218 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1219 fHistograms->FillHistogram("ESD_TrueConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1220 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1221 fHistograms->FillHistogram("ESD_TrueConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1222 fHistograms->FillHistogram("ESD_TrueConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1223 fHistograms->FillHistogram("ESD_TrueConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1224 fHistograms->FillHistogram("ESD_TrueConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1225 fHistograms->FillHistogram("ESD_TrueConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1226 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1227 fHistograms->FillHistogram("ESD_TrueConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1228 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters());
1229 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters());
1230 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters(),fV0Reader->GetMotherCandidateMass());
1231 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters(),fV0Reader->GetMotherCandidateMass());
a0b94e5c 1232
d7d7e825 1233 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1234 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
a0b94e5c 1235
d7d7e825 1236
1237 fHistograms->FillHistogram("ESD_TrueConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1238 fHistograms->FillHistogram("ESD_TrueConversion_R", fV0Reader->GetXYRadius());
1239 fHistograms->FillHistogram("ESD_TrueConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1240 fHistograms->FillHistogram("ESD_TrueConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1241
1242 fHistograms->FillHistogram("ESD_TrueConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1243 fHistograms->FillHistogram("ESD_TrueConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1244 fHistograms->FillHistogram("ESD_TrueConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1245 fHistograms->FillHistogram("ESD_TrueConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1246
1247 fHistograms->FillHistogram("ESD_TrueConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1248 fHistograms->FillHistogram("ESD_TrueConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1249 fHistograms->FillHistogram("ESD_TrueConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1250 fHistograms->FillHistogram("ESD_TrueConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1251
1252
a0b94e5c 1253
d7d7e825 1254 //store MCTruth properties
1255 fHistograms->FillHistogram("ESD_TrueConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1256 fHistograms->FillHistogram("ESD_TrueConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1257 fHistograms->FillHistogram("ESD_TrueConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
a0b94e5c 1258
d7d7e825 1259 //resolution
1260 Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();
1261 Double_t esdpt = fV0Reader->GetMotherCandidatePt();
1262 Double_t resdPt = 0;
4a6157dc 1263 if(mcpt > 0){
d7d7e825 1264 resdPt = ((esdpt - mcpt)/mcpt)*100;
1265 }
4a6157dc 1266 else if(mcpt < 0){
1267 cout<<"Pt of MC particle is negative, this will cause wrong calculation of resPt"<<endl;
1268 }
d7d7e825 1269
1270 fHistograms->FillHistogram("Resolution_dPt", mcpt, resdPt);
1271 fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);
1272 fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);
1273
1274 Double_t resdZ = 0;
1275 if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){
1276 resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100;
1277 }
037dc2db 1278 Double_t resdZAbs = 0;
1279 resdZAbs = (fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz());
1280
1281 fHistograms->FillHistogram("Resolution_dZAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdZAbs);
d7d7e825 1282 fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);
1283 fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());
1284 fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());
1285
1286 Double_t resdR = 0;
1287 if(fV0Reader->GetNegativeMCParticle()->R() != 0){
1288 resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100;
1289 }
037dc2db 1290 Double_t resdRAbs = 0;
1291 resdRAbs = (fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R());
1292
1293 fHistograms->FillHistogram("Resolution_dRAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdRAbs);
d7d7e825 1294 fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);
1295 fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());
1296 fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());
1297 fHistograms->FillHistogram("Resolution_dR_dPt", resdR, resdPt);
037dc2db 1298
1299 Double_t resdPhiAbs=0;
1300 resdPhiAbs=0;
1301 resdPhiAbs= (fV0Reader->GetMotherCandidatePhi()-fV0Reader->GetNegativeMCParticle()->Phi());
1302 fHistograms->FillHistogram("Resolution_dPhiAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdPhiAbs);
1303
d7d7e825 1304 }//if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22)
d7d7e825 1305 }//if(fDoMCTruth)
1306 }//while(fV0Reader->NextV0)
1307 fHistograms->FillHistogram("ESD_NumberOfSurvivingV0s", nSurvivingV0s);
1308 fHistograms->FillHistogram("ESD_NumberOfV0s", fV0Reader->GetNumberOfV0s());
b5832f95 1309 fHistograms->FillHistogram("ESD_NumberOfContributorsVtx", fV0Reader->GetNumberOfContributorsVtx());
cb90a330 1310
1311 fV0Reader->ResetV0IndexNumber();
1312
d7d7e825 1313}
1314
1315void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){
1316 // Fill AOD with reconstructed Gamma
a0b94e5c 1317
6c84d371 1318 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
1319 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
d7d7e825 1320 //Create AOD particle object from AliKFParticle
a0b94e5c 1321
d7d7e825 1322 /* AliKFParticle * gammakf = &fKFReconstructedGammas[gammaIndex];
1323 //You could add directly AliKFParticle objects to the AOD, avoiding dependences with PartCorr
1324 //but this means that I have to work a little bit more in my side.
1325 //AODPWG4Particle objects are simpler and lighter, I think
1326 AliAODPWG4Particle gamma = AliAODPWG4Particle(gammakf->Px(),gammakf->Py(),gammakf->Pz(), gammakf->E());
1327 gamma.SetLabel(-1);//How to get the MC label of the reconstructed gamma?
a0b94e5c 1328 gamma.SetCaloLabel(-1,-1); //How to get the MC label of the 2 electrons that form the gamma?
1329 gamma.SetDetector("CTS"); //tag the gamma as reconstructed in the central barrel
1330 gamma.SetPdg(AliCaloPID::kPhotonConv); //photon id
1331 gamma.SetTag(-1); //Here I usually put a flag saying that montecarlo says it is prompt, decay fragmentation photon, or hadrons or whatever
1332
1333 //Add it to the aod list
1334 Int_t i = fAODBranch->GetEntriesFast();
1335 new((*fAODBranch)[i]) AliAODPWG4Particle(gamma);
d7d7e825 1336 */
6c84d371 1337 // AliKFParticle * gammakf = &fKFReconstructedGammas[gammaIndex];
1338 AliKFParticle * gammakf = (AliKFParticle *)fKFReconstructedGammasTClone->At(gammaIndex);
d7d7e825 1339 AliGammaConversionAODObject aodObject;
1340 aodObject.SetPx(gammakf->GetPx());
1341 aodObject.SetPy(gammakf->GetPy());
1342 aodObject.SetPz(gammakf->GetPz());
1343 aodObject.SetLabel1(fElectronv1[gammaIndex]);
1344 aodObject.SetLabel2(fElectronv2[gammaIndex]);
1345 Int_t i = fAODBranch->GetEntriesFast();
6c84d371 1346 new((*fAODBranch)[i]) AliGammaConversionAODObject(aodObject);
a0b94e5c 1347 }
1348
d7d7e825 1349}
1350
1351
1352void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
1353 // see header file for documentation
1354
6c84d371 1355 // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
1356 // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
037dc2db 1357
1358 if(fKFReconstructedGammasTClone->GetEntriesFast()>fV0Reader->GetNumberOfV0s()){
1359 cout<<"Warning, number of entries in the tclone is bigger than number of v0s"<<endl;
1360 }
1361
6c84d371 1362 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
1363 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
d7d7e825 1364
6c84d371 1365 // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
1366 // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
a0b94e5c 1367
6c84d371 1368 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
1369 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
a0b94e5c 1370
d7d7e825 1371 if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
1372 continue;
1373 }
1374 if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
1375 continue;
1376 }
a0b94e5c 1377
d7d7e825 1378 AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
1379
1380 Double_t massTwoGammaCandidate = 0.;
1381 Double_t widthTwoGammaCandidate = 0.;
1382 Double_t chi2TwoGammaCandidate =10000.;
1383 twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
1384 if(twoGammaCandidate->GetNDF()>0){
1385 chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
1386
6c84d371 1387 if(chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()){
d7d7e825 1388
1389 TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
1390 TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
1391
1392 Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
1393 Double_t rapidity;
1394 if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() == 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() == 0){
1395 rapidity=0;
1396 }
1397 else{
1398 rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
1399 }
1400
1e7846f4 1401 if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
1402 delete twoGammaCandidate;
1403 continue; // minimum opening angle to avoid using ghosttracks
1404 }
1405
d7d7e825 1406 fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
1407 fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
1408 fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
1409 fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
1410 fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
1411 fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
1412 fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
1413 fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
1414 fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
1415 fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
1416 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1417 fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
9640a3d1 1418
ebcfaa7e 1419 // if(fDoNeutralMesonV0MCCheck){
1420 if(fDoMCTruth){
037dc2db 1421 //Kenneth: Checking the eta of the gamma to check the difference between 0.9 and 1.2
1422 Int_t indexKF1 = fKFReconstructedGammasV0Index.at(firstGammaIndex);
1423 if(indexKF1<fV0Reader->GetNumberOfV0s()){
1424 fV0Reader->GetV0(indexKF1);//updates to the correct v0
1425 Double_t eta1 = fV0Reader->GetMotherCandidateEta();
1426 Bool_t isRealPi0=kFALSE;
1427 Int_t gamma1MotherLabel=-1;
1428 if(fV0Reader->HasSameMCMother() == kTRUE){
1429 //cout<<"This v0 is a real v0!!!!"<<endl;
1430 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1431 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1432 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
1433 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
1434 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
1435 gamma1MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
1436 }
1437 }
1438 }
1439 }
1440 Int_t indexKF2 = fKFReconstructedGammasV0Index.at(secondGammaIndex);
1441 if(indexKF1 == indexKF2){
1442 cout<<"index of the two KF particles are the same.... should not happen"<<endl;
1443 }
1444 if(indexKF2<fV0Reader->GetNumberOfV0s()){
1445 fV0Reader->GetV0(indexKF2);
1446 Double_t eta2 = fV0Reader->GetMotherCandidateEta();
1447 Int_t gamma2MotherLabel=-1;
1448 if(fV0Reader->HasSameMCMother() == kTRUE){
1449 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1450 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1451 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
1452 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
1453 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
1454 gamma2MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
1455 }
1456 }
1457 }
1458 }
1459 if(gamma1MotherLabel>=0 && gamma1MotherLabel==gamma2MotherLabel){
1460 if(fV0Reader->CheckIfPi0IsMother(gamma1MotherLabel)){
1461 isRealPi0=kTRUE;
1462 }
1463 }
1464
1465 if(TMath::Abs(eta1)>0.9 && TMath::Abs(eta2)>0.9){
ebcfaa7e 1466 // fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
1467 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
037dc2db 1468 if(isRealPi0){
1469 fHistograms->FillHistogram("ESD_TruePi0_InvMass_1212",massTwoGammaCandidate);
1470 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_1212",openingAngleTwoGammaCandidate);
1471 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 1472 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1473 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
037dc2db 1474 }
1475 }
1476 else if(TMath::Abs(eta1)>0.9 || TMath::Abs(eta2)>0.9){
ebcfaa7e 1477 // fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
1478 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
037dc2db 1479 if(isRealPi0){
1480 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0912",massTwoGammaCandidate);
1481 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0912",openingAngleTwoGammaCandidate);
1482 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 1483 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1484 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
037dc2db 1485 }
1486 }
1487 else{
ebcfaa7e 1488 // fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
1489 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
037dc2db 1490 if(isRealPi0){
1491 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0909",massTwoGammaCandidate);
1492 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0909",openingAngleTwoGammaCandidate);
1493 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 1494 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1495 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
037dc2db 1496 }
1497 }
1498 }
1499 }
1500 }
9640a3d1 1501 if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
1502 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
1503 fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
1504 }
ebcfaa7e 1505
1506 if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
1507 fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
1508 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
1509 }
1510 else if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 || TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
1511 fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
1512 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
1513 }
1514 else{
1515 fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
1516 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
1517 }
d7d7e825 1518 }
1519 }
1520 delete twoGammaCandidate;
d7d7e825 1521 }
1522 }
1523}
1524
1525void AliAnalysisTaskGammaConversion::CalculateBackground(){
1526 // see header file for documentation
5e55d806 1527
1528
1529 TClonesArray * currentEventV0s = fV0Reader->GetCurrentEventGoodV0s();
1530
1531 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
1532 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
1533 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
1534 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
1535 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
1536 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
1537
1538 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
d7d7e825 1539
5e55d806 1540 Double_t massBG =0.;
1541 Double_t widthBG = 0.;
1542 Double_t chi2BG =10000.;
1543 backgroundCandidate->GetMass(massBG,widthBG);
1544 if(backgroundCandidate->GetNDF()>0){
1545 chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
1546 if(chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()){
d7d7e825 1547
5e55d806 1548 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
1549 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
d7d7e825 1550
5e55d806 1551 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
a0b94e5c 1552
5e55d806 1553 Double_t rapidity;
1554 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
1555 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
a0b94e5c 1556
d7d7e825 1557
1558
1559
5e55d806 1560 if(openingAngleBG < fMinOpeningAngleGhostCut ){
1561 delete backgroundCandidate;
1562 continue; // minimum opening angle to avoid using ghosttracks
1563 }
037dc2db 1564
1565 // original
5e55d806 1566 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
1567 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
1568 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
1569 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
1570 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
1571 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
1572 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
1573 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
1574 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
1575 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
1576 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
1577 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
1578
1579 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
1580 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
1581 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
1582 }
037dc2db 1583
1584
1585 // test
1586 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
1587
1588 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
1589 Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
1590
1591 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
1592 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
1593 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
1594 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
1595 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
1596 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
1597 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
1598 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
1599 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
1600 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
1601 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
1602 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
9640a3d1 1603
037dc2db 1604 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
1605 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
1606 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
1607 }
5e55d806 1608 }
d7d7e825 1609 }
5e55d806 1610 delete backgroundCandidate;
d7d7e825 1611 }
d7d7e825 1612 }
1613 }
1614}
1615
1616
d7d7e825 1617void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){
1618 //ProcessGammasForGammaJetAnalysis
a0b94e5c 1619
d7d7e825 1620 Double_t distIsoMin;
a0b94e5c 1621
d7d7e825 1622 CreateListOfChargedParticles();
a0b94e5c 1623
1624
6c84d371 1625 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
1626 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
1627 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
01b7fdcc 1628 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 1629
01b7fdcc 1630 if( momentumVectorCurrentGamma.Pt()> fMinPtForGammaJet){
d7d7e825 1631 distIsoMin=GetMinimumDistanceToCharge(gammaIndex);
a0b94e5c 1632
d7d7e825 1633 if (distIsoMin > fMinIsoConeSize && fLeadingChargedIndex>=0){
01b7fdcc 1634 CalculateJetCone(gammaIndex);
d7d7e825 1635 }
1636 }
1637 }
1638}
1639
1640void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){
1641 // CreateListOfChargedParticles
a0b94e5c 1642
d7d7e825 1643 fESDEvent = fV0Reader->GetESDEvent();
037dc2db 1644 Int_t numberOfESDTracks=0;
d7d7e825 1645 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
1646 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 1647
d7d7e825 1648 if(!curTrack){
1649 continue;
1650 }
a0b94e5c 1651
d7d7e825 1652 if(fEsdTrackCuts->AcceptTrack(curTrack) ){
6c84d371 1653 new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack);
1654 // fChargedParticles.push_back(curTrack);
d7d7e825 1655 fChargedParticlesId.push_back(iTracks);
037dc2db 1656 numberOfESDTracks++;
d7d7e825 1657 }
1658 }
037dc2db 1659 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",numberOfESDTracks);
d7d7e825 1660}
01b7fdcc 1661void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){
6c84d371 1662 // CaculateJetCone
a0b94e5c 1663
d7d7e825 1664 Double_t cone;
1665 Double_t coneSize=0.3;
1666 Double_t ptJet=0;
a0b94e5c 1667
6c84d371 1668 // AliKFParticle * currentGamma = &fKFReconstructedGammas[gammaIndex];
1669 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
98778c17 1670
01b7fdcc 1671 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 1672
6c84d371 1673 AliESDtrack* leadingCharged = (AliESDtrack*)(fChargedParticles->At(fLeadingChargedIndex));
98778c17 1674
d7d7e825 1675 Double_t momLeadingCharged[3];
1676 leadingCharged->GetConstrainedPxPyPz(momLeadingCharged);
a0b94e5c 1677
01b7fdcc 1678 TVector3 momentumVectorLeadingCharged(momLeadingCharged[0],momLeadingCharged[1],momLeadingCharged[2]);
a0b94e5c 1679
01b7fdcc 1680 Double_t phi1=momentumVectorLeadingCharged.Phi();
1681 Double_t eta1=momentumVectorLeadingCharged.Eta();
1682 Double_t phi3=momentumVectorCurrentGamma.Phi();
a0b94e5c 1683
6c84d371 1684 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
1685 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 1686 Int_t chId = fChargedParticlesId[iCh];
1687 if(fLeadingChargedIndex==chId || fLeadingChargedIndex==chId) continue;
1688 Double_t mom[3];
1689 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 1690 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
1691 Double_t phi2=momentumVectorChargedParticle.Phi();
1692 Double_t eta2=momentumVectorChargedParticle.Eta();
a0b94e5c 1693
1694
d7d7e825 1695 cone=100.;
1696 if( TMath::Abs(phi2 - phi1) <= ( TMath::TwoPi()-coneSize) ){
1697 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-phi1),2) );
1698 }else{
1699 if( (phi2 - phi1)> TMath::TwoPi()-coneSize ){
1700 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-TMath::TwoPi()-phi1),2) );
1701 }
1702 if( (phi2 - phi1)< -(TMath::TwoPi()-coneSize) ){
1703 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2+TMath::TwoPi()-phi1),2) );
1704 }
1705 }
a0b94e5c 1706
01b7fdcc 1707 if(cone <coneSize&& momentumVectorChargedParticle.Pt()>fMinPtJetCone ){
1708 ptJet+= momentumVectorChargedParticle.Pt();
1709 Double_t ffzHdrGam = momentumVectorChargedParticle.Pt()/momentumVectorCurrentGamma.Pt();
1710 Double_t imbalanceHdrGam=-momentumVectorChargedParticle.Dot(momentumVectorCurrentGamma)/momentumVectorCurrentGamma.Mag2();
d7d7e825 1711 fHistograms->FillHistogram("ESD_FFzHdrGam",ffzHdrGam);
1712 fHistograms->FillHistogram("ESD_ImbalanceHdrGam",imbalanceHdrGam);
a0b94e5c 1713
d7d7e825 1714 }
a0b94e5c 1715
d7d7e825 1716 Double_t dphiHdrGam=phi3-phi2;
1717 if ( dphiHdrGam < (-TMath::PiOver2())){
1718 dphiHdrGam+=(TMath::TwoPi());
1719 }
a0b94e5c 1720
d7d7e825 1721 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
1722 dphiHdrGam-=(TMath::TwoPi());
1723 }
a0b94e5c 1724
01b7fdcc 1725 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 1726 fHistograms->FillHistogram("ESD_dphiHdrGamIsolated",dphiHdrGam);
1727 }
1728 }//track loop
a0b94e5c 1729
1730
d7d7e825 1731}
1732
1733Double_t AliAnalysisTaskGammaConversion::GetMinimumDistanceToCharge(Int_t indexHighestPtGamma){
1734 // GetMinimumDistanceToCharge
a0b94e5c 1735
d7d7e825 1736 Double_t fIsoMin=100.;
1737 Double_t ptLeadingCharged=-1.;
98778c17 1738
1739 fLeadingChargedIndex=-1;
a0b94e5c 1740
6c84d371 1741 AliKFParticle * gammaHighestPt = (AliKFParticle*)fKFReconstructedGammasTClone->At(indexHighestPtGamma);
01b7fdcc 1742 TVector3 momentumVectorgammaHighestPt(gammaHighestPt->GetPx(),gammaHighestPt->GetPy(),gammaHighestPt->GetPz());
a0b94e5c 1743
01b7fdcc 1744 Double_t phi1=momentumVectorgammaHighestPt.Phi();
1745 Double_t eta1=momentumVectorgammaHighestPt.Eta();
a0b94e5c 1746
6c84d371 1747 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
1748 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 1749 Int_t chId = fChargedParticlesId[iCh];
1750 if(fElectronv1[indexHighestPtGamma]==chId || fElectronv2[indexHighestPtGamma]==chId) continue;
1751 Double_t mom[3];
1752 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 1753 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
1754 Double_t phi2=momentumVectorChargedParticle.Phi();
1755 Double_t eta2=momentumVectorChargedParticle.Eta();
d7d7e825 1756 Double_t iso=pow( (pow( (eta1-eta2),2)+ pow((phi1-phi2),2)),0.5 );
a0b94e5c 1757
01b7fdcc 1758 if(momentumVectorChargedParticle.Pt()>fMinPtIsoCone ){
d7d7e825 1759 if (iso<fIsoMin){
1760 fIsoMin=iso;
1761 }
1762 }
a0b94e5c 1763
d7d7e825 1764 Double_t dphiHdrGam=phi1-phi2;
1765 if ( dphiHdrGam < (-TMath::PiOver2())){
1766 dphiHdrGam+=(TMath::TwoPi());
1767 }
a0b94e5c 1768
d7d7e825 1769 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
1770 dphiHdrGam-=(TMath::TwoPi());
1771 }
01b7fdcc 1772 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 1773 fHistograms->FillHistogram("ESD_dphiHdrGam",dphiHdrGam);
1774 }
a0b94e5c 1775
d7d7e825 1776 if (dphiHdrGam>0.9*TMath::Pi() && dphiHdrGam<1.1*TMath::Pi()){
a0b94e5c 1777 if (momentumVectorChargedParticle.Pt()> ptLeadingCharged && momentumVectorChargedParticle.Pt()>0.1*momentumVectorgammaHighestPt.Pt()){
01b7fdcc 1778 ptLeadingCharged=momentumVectorChargedParticle.Pt();
d7d7e825 1779 fLeadingChargedIndex=iCh;
1780 }
1781 }
a0b94e5c 1782
d7d7e825 1783 }//track loop
1784 fHistograms->FillHistogram("ESD_MinimumIsoDistance",fIsoMin);
1785 return fIsoMin;
a0b94e5c 1786
d7d7e825 1787}
1788
a0b94e5c 1789Int_t AliAnalysisTaskGammaConversion::GetIndexHighestPtGamma(){
1790 //GetIndexHighestPtGamma
1791
d7d7e825 1792 Int_t indexHighestPtGamma=-1;
01b7fdcc 1793 //Double_t
1794 fGammaPtHighest = -100.;
a0b94e5c 1795
6c84d371 1796 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
1797 AliKFParticle * gammaHighestPtCandidate = (AliKFParticle*)fKFReconstructedGammasTClone->At(firstGammaIndex);
a0b94e5c 1798 TVector3 momentumVectorgammaHighestPtCandidate(gammaHighestPtCandidate->GetPx(),gammaHighestPtCandidate->GetPy(),gammaHighestPtCandidate->GetPz());
1799 if (momentumVectorgammaHighestPtCandidate.Pt() > fGammaPtHighest){
1800 fGammaPtHighest=momentumVectorgammaHighestPtCandidate.Pt();
1801 //gammaHighestPt = gammaHighestPtCandidate;
1802 indexHighestPtGamma=firstGammaIndex;
1803 }
d7d7e825 1804 }
a0b94e5c 1805
d7d7e825 1806 return indexHighestPtGamma;
a0b94e5c 1807
d7d7e825 1808}
1809
1810
1811void AliAnalysisTaskGammaConversion::Terminate(Option_t */*option*/)
1812{
1813 // Terminate analysis
1814 //
1815 AliDebug(1,"Do nothing in Terminate");
1816}
1817
1818void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()
1819{
1820 //AOD
1e7846f4 1821 if(fAODBranch==NULL){
1822 fAODBranch = new TClonesArray("AliGammaConversionAODObject", 0);
1823 }
1824 fAODBranch->Delete();
d7d7e825 1825 fAODBranch->SetName(fAODBranchName);
1826 AddAODBranch("TClonesArray", &fAODBranch);
a0b94e5c 1827
d7d7e825 1828 // Create the output container
1829 if(fOutputContainer != NULL){
1830 delete fOutputContainer;
1831 fOutputContainer = NULL;
1832 }
1833 if(fOutputContainer == NULL){
1834 fOutputContainer = new TList();
b58d3c74 1835 fOutputContainer->SetOwner(kTRUE);
d7d7e825 1836 }
1837
1838 //Adding the histograms to the output container
1839 fHistograms->GetOutputContainer(fOutputContainer);
1840
1841
1842 if(fWriteNtuple){
1843 if(fGammaNtuple == NULL){
1844 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);
1845 }
1846 if(fNeutralMesonNtuple == NULL){
1847 fNeutralMesonNtuple = new TNtuple("NeutralMesonNtuple","NeutralMesonNtuple","test");
1848 }
1849 TList * ntupleTList = new TList();
b58d3c74 1850 ntupleTList->SetOwner(kTRUE);
d7d7e825 1851 ntupleTList->SetName("Ntuple");
1852 ntupleTList->Add((TNtuple*)fGammaNtuple);
1853 fOutputContainer->Add(ntupleTList);
1854 }
1855
1856 fOutputContainer->SetName(GetName());
1857}
1858
1859Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(TParticle* const daughter0, TParticle* const daughter1) const{
1860 //helper function
1861 TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());
1862 TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());
1863 return v3D0.Angle(v3D1);
1864}
1865
1866void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){
1867 // see header file for documentation
5e55d806 1868
d7d7e825 1869 vector<Int_t> indexOfGammaParticle;
a0b94e5c 1870
d7d7e825 1871 fStack = fV0Reader->GetMCStack();
a0b94e5c 1872
d7d7e825 1873 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
1874 return; // aborts if the primary vertex does not have contributors.
1875 }
a0b94e5c 1876
d7d7e825 1877 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
1878 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
1879 if(particle->GetPdgCode()==22){ //Gamma
1880 if(particle->GetNDaughters() >= 2){
1881 TParticle* electron=NULL;
1882 TParticle* positron=NULL;
1883 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
1884 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
1885 if(tmpDaughter->GetUniqueID() == 5){
1886 if(tmpDaughter->GetPdgCode() == 11){
1887 electron = tmpDaughter;
1888 }
1889 else if(tmpDaughter->GetPdgCode() == -11){
1890 positron = tmpDaughter;
1891 }
1892 }
1893 }
1894 if(electron!=NULL && positron!=0){
1895 if(electron->R()<160){
1896 indexOfGammaParticle.push_back(iTracks);
1897 }
1898 }
1899 }
1900 }
1901 }
a0b94e5c 1902
d7d7e825 1903 Int_t nFoundGammas=0;
1904 Int_t nNotFoundGammas=0;
a0b94e5c 1905
d7d7e825 1906 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1907 for(Int_t i=0;i<numberOfV0s;i++){
1908 fV0Reader->GetV0(i);
a0b94e5c 1909
d7d7e825 1910 if(fV0Reader->HasSameMCMother() == kFALSE){
1911 continue;
1912 }
a0b94e5c 1913
d7d7e825 1914 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1915 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 1916
d7d7e825 1917 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1918 continue;
1919 }
1920 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
1921 continue;
1922 }
a0b94e5c 1923
d7d7e825 1924 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
1925 //TParticle * v0Gamma = fV0Reader->GetMotherMCParticle();
1926 for(UInt_t mcIndex=0;mcIndex<indexOfGammaParticle.size();mcIndex++){
1927 if(negativeMC->GetFirstMother()==indexOfGammaParticle[mcIndex]){
1928 nFoundGammas++;
1929 }
1930 else{
1931 nNotFoundGammas++;
1932 }
1933 }
1934 }
1935 }
d7d7e825 1936}
1937
1938
1939void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){
1940 // see header file for documantation
a0b94e5c 1941
d7d7e825 1942 fESDEvent = fV0Reader->GetESDEvent();
a0b94e5c 1943
1944
1945 TClonesArray * vESDeNegTemp = new TClonesArray("AliESDtrack",0);
6c84d371 1946 TClonesArray * vESDePosTemp = new TClonesArray("AliESDtrack",0);
1947 TClonesArray * vESDxNegTemp = new TClonesArray("AliESDtrack",0);
1948 TClonesArray * vESDxPosTemp = new TClonesArray("AliESDtrack",0);
1949 TClonesArray * vESDeNegNoJPsi = new TClonesArray("AliESDtrack",0);
1950 TClonesArray * vESDePosNoJPsi = new TClonesArray("AliESDtrack",0);
a0b94e5c 1951
6c84d371 1952 /*
1953 vector <AliESDtrack*> vESDeNegTemp(0);
1954 vector <AliESDtrack*> vESDePosTemp(0);
1955 vector <AliESDtrack*> vESDxNegTemp(0);
1956 vector <AliESDtrack*> vESDxPosTemp(0);
1957 vector <AliESDtrack*> vESDeNegNoJPsi(0);
1958 vector <AliESDtrack*> vESDePosNoJPsi(0);
1959 */
a0b94e5c 1960
1961
d7d7e825 1962 fHistograms->FillTable("Table_Electrons",0);//Count number of Events
a0b94e5c 1963
d7d7e825 1964 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
1965 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 1966
d7d7e825 1967 if(!curTrack){
1968 //print warning here
1969 continue;
1970 }
a0b94e5c 1971
d7d7e825 1972 double p[3];if(!curTrack->GetConstrainedPxPyPz(p))continue;
1973 double r[3];curTrack->GetConstrainedXYZ(r);
a0b94e5c 1974
d7d7e825 1975 TVector3 rXYZ(r);
a0b94e5c 1976
d7d7e825 1977 fHistograms->FillTable("Table_Electrons",4);//Count number of ESD tracks
a0b94e5c 1978
d7d7e825 1979 Bool_t flagKink = kTRUE;
1980 Bool_t flagTPCrefit = kTRUE;
1981 Bool_t flagTRDrefit = kTRUE;
1982 Bool_t flagITSrefit = kTRUE;
1983 Bool_t flagTRDout = kTRUE;
1984 Bool_t flagVertex = kTRUE;
a0b94e5c 1985
1986
d7d7e825 1987 //Cuts ---------------------------------------------------------------
a0b94e5c 1988
d7d7e825 1989 if(curTrack->GetKinkIndex(0) > 0){
1990 fHistograms->FillHistogram("Table_Electrons",5);//Count kink
1991 flagKink = kFALSE;
1992 }
a0b94e5c 1993
d7d7e825 1994 ULong_t trkStatus = curTrack->GetStatus();
a0b94e5c 1995
d7d7e825 1996 ULong_t tpcRefit = (trkStatus & AliESDtrack::kTPCrefit);
a0b94e5c 1997
d7d7e825 1998 if(!tpcRefit){
1999 fHistograms->FillHistogram("Table_Electrons",9);//Count not TPCrefit
2000 flagTPCrefit = kFALSE;
2001 }
a0b94e5c 2002
d7d7e825 2003 ULong_t itsRefit = (trkStatus & AliESDtrack::kITSrefit);
2004 if(!itsRefit){
2005 fHistograms->FillHistogram("Table_Electrons",10);//Count not ITSrefit
2006 flagITSrefit = kFALSE;
2007 }
a0b94e5c 2008
d7d7e825 2009 ULong_t trdRefit = (trkStatus & AliESDtrack::kTRDrefit);
a0b94e5c 2010
d7d7e825 2011 if(!trdRefit){
2012 fHistograms->FillHistogram("Table_Electrons",8); //Count not TRDrefit
2013 flagTRDrefit = kFALSE;
2014 }
a0b94e5c 2015
d7d7e825 2016 ULong_t trdOut = (trkStatus & AliESDtrack::kTRDout);
a0b94e5c 2017
d7d7e825 2018 if(!trdOut) {
2019 fHistograms->FillHistogram("Table_Electrons",7); //Count not TRDout
2020 flagTRDout = kFALSE;
2021 }
a0b94e5c 2022
d7d7e825 2023 double nsigmaToVxt = GetSigmaToVertex(curTrack);
a0b94e5c 2024
d7d7e825 2025 if(nsigmaToVxt > 3){
2026 fHistograms->FillHistogram("Table_Electrons",6); //Count Tracks with number of sigmas > 3
2027 flagVertex = kFALSE;
2028 }
a0b94e5c 2029
d7d7e825 2030 if(! (flagKink && flagTPCrefit && flagITSrefit && flagTRDrefit && flagTRDout && flagVertex ) ) continue;
2031 fHistograms->FillHistogram("Table_Electrons",11);//Count Tracks passed Cuts
a0b94e5c 2032
2033
d7d7e825 2034 Stat_t pid, weight;
2035 GetPID(curTrack, pid, weight);
a0b94e5c 2036
d7d7e825 2037 if(pid!=0){
2038 fHistograms->FillHistogram("Table_Electrons",12); //Count Tracks with pid != 0
2039 }
a0b94e5c 2040
d7d7e825 2041 if(pid == 0){
2042 fHistograms->FillHistogram("Table_Electrons",13); //Count Tracks with pid != 0
2043 }
a0b94e5c 2044
2045
2046
2047
a0b94e5c 2048
2049
d7d7e825 2050 TLorentzVector curElec;
2051 curElec.SetXYZM(p[0],p[1],p[2],fElectronMass);
a0b94e5c 2052
2053
113d8432 2054 if(fDoMCTruth){
2055 Int_t labelMC = TMath::Abs(curTrack->GetLabel());
2056 TParticle* curParticle = fStack->Particle(labelMC);
2057 if(curTrack->GetSign() > 0){
2058 if( pid == 0){
2059 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
2060 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
2061 }
2062 else{
2063 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
2064 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
2065 }
2066 }
2067 }
a0b94e5c 2068
2069
d7d7e825 2070 if(curTrack->GetSign() > 0){
a0b94e5c 2071
6c84d371 2072 // vESDxPosTemp.push_back(curTrack);
2073 new((*vESDxPosTemp)[vESDxPosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 2074
d7d7e825 2075 if( pid == 0){
a0b94e5c 2076
d7d7e825 2077 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
2078 fHistograms->FillHistogram("ESD_ElectronPosPt",curElec.Pt());
113d8432 2079 // fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
d7d7e825 2080 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
113d8432 2081 // fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
6c84d371 2082 // vESDePosTemp.push_back(curTrack);
2083 new((*vESDePosTemp)[vESDePosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 2084
d7d7e825 2085 }
a0b94e5c 2086
d7d7e825 2087 }
2088 else {
5e55d806 2089
6c84d371 2090 new((*vESDxNegTemp)[vESDxNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 2091
d7d7e825 2092 if( pid == 0){
a0b94e5c 2093
d7d7e825 2094 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
2095 fHistograms->FillHistogram("ESD_ElectronNegPt",curElec.Pt());
d7d7e825 2096 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
6c84d371 2097 new((*vESDeNegTemp)[vESDeNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 2098
d7d7e825 2099 }
a0b94e5c 2100
d7d7e825 2101 }
a0b94e5c 2102
d7d7e825 2103 }
a0b94e5c 2104
2105
d7d7e825 2106 Bool_t ePosJPsi = kFALSE;
2107 Bool_t eNegJPsi = kFALSE;
2108 Bool_t ePosPi0 = kFALSE;
2109 Bool_t eNegPi0 = kFALSE;
2110
2111 UInt_t iePosJPsi=0,ieNegJPsi=0,iePosPi0=0,ieNegPi0=0;
a0b94e5c 2112
6c84d371 2113 for(Int_t iNeg=0; iNeg < vESDeNegTemp->GetEntriesFast(); iNeg++){
2114 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetPdgCode() == 11)
2115 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0) > -1){
2116 Int_t labelMother = fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0);
d7d7e825 2117 TParticle* partMother = fStack ->Particle(labelMother);
2118 if (partMother->GetPdgCode() == 111){
2119 ieNegPi0 = iNeg;
2120 eNegPi0 = kTRUE;
2121 }
2122 if(partMother->GetPdgCode() == 443){ //Mother JPsi
2123 fHistograms->FillTable("Table_Electrons",14);
2124 ieNegJPsi = iNeg;
2125 eNegJPsi = kTRUE;
2126 }
2127 else{
6c84d371 2128 // vESDeNegNoJPsi.push_back(vESDeNegTemp[iNeg]);
2129 new((*vESDeNegNoJPsi)[vESDeNegNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDeNegTemp->At(iNeg)));
d7d7e825 2130 // cout<<"ESD No Positivo JPsi "<<endl;
2131 }
a0b94e5c 2132
d7d7e825 2133 }
2134 }
a0b94e5c 2135
6c84d371 2136 for(Int_t iPos=0; iPos < vESDePosTemp->GetEntriesFast(); iPos++){
2137 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetPdgCode() == -11)
2138 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0) > -1){
2139 Int_t labelMother = fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0);
d7d7e825 2140 TParticle* partMother = fStack ->Particle(labelMother);
2141 if (partMother->GetPdgCode() == 111){
2142 iePosPi0 = iPos;
2143 ePosPi0 = kTRUE;
2144 }
2145 if(partMother->GetPdgCode() == 443){ //Mother JPsi
2146 fHistograms->FillTable("Table_Electrons",15);
2147 iePosJPsi = iPos;
2148 ePosJPsi = kTRUE;
2149 }
2150 else{
6c84d371 2151 // vESDePosNoJPsi.push_back(vESDePosTemp[iPos]);
2152 new((*vESDePosNoJPsi)[vESDePosNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDePosTemp->At(iPos)));
d7d7e825 2153 // cout<<"ESD No Negativo JPsi "<<endl;
2154 }
a0b94e5c 2155
d7d7e825 2156 }
2157 }
2158
2159 if( eNegJPsi && ePosJPsi ){
2160 TVector3 tempeNegV,tempePosV;
6c84d371 2161 tempeNegV.SetXYZ(((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Px(),((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Py(),((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Pz());
2162 tempePosV.SetXYZ(((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Px(),((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Py(),((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Pz());
d7d7e825 2163 fHistograms->FillTable("Table_Electrons",16);
2164 fHistograms->FillHistogram("ESD_ElectronPosNegJPsiAngle",tempeNegV.Angle(tempePosV));
6c84d371 2165 fHistograms->FillHistogram("MC_ElectronPosNegJPsiAngle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->GetLabel())),
2166 fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->GetLabel()))));
d7d7e825 2167 }
2168
2169 if( eNegPi0 && ePosPi0 ){
2170 TVector3 tempeNegV,tempePosV;
6c84d371 2171 tempeNegV.SetXYZ(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Px(),((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Py(),((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Pz());
2172 tempePosV.SetXYZ(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Px(),((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Py(),((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Pz());
d7d7e825 2173 fHistograms->FillHistogram("ESD_ElectronPosNegPi0Angle",tempeNegV.Angle(tempePosV));
6c84d371 2174 fHistograms->FillHistogram("MC_ElectronPosNegPi0Angle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->GetLabel())),
2175 fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->GetLabel()))));
d7d7e825 2176 }
a0b94e5c 2177
2178
d7d7e825 2179 FillAngle("ESD_eNegePosAngleBeforeCut",GetTLorentzVector(vESDeNegTemp),GetTLorentzVector(vESDePosTemp));
a0b94e5c 2180
6c84d371 2181 CleanWithAngleCuts(*vESDeNegTemp,*vESDePosTemp,*fKFReconstructedGammasTClone);
d7d7e825 2182
6c84d371 2183 // vector <TLorentzVector> vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectron);
2184 // vector <TLorentzVector> vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectron);
a0b94e5c 2185
6c84d371 2186 TClonesArray vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectronTClone);
2187 TClonesArray vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectronTClone);
a0b94e5c 2188
2189
d7d7e825 2190 FillAngle("ESD_eNegePosAngleAfterCut",vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 2191
2192
2193
2194
d7d7e825 2195 //FillAngle("ESD_eNegePosAngleAfterCut",CurrentTLVeNeg,CurrentTLVePos);
a0b94e5c 2196
2197
d7d7e825 2198 FillElectronInvMass("ESD_InvMass_ePluseMinus",vCurrentTLVeNeg,vCurrentTLVePos);
2199 FillElectronInvMass("ESD_InvMass_xPlusxMinus",GetTLorentzVector(vESDxNegTemp),GetTLorentzVector(vESDxPosTemp));
a0b94e5c 2200
2201
2202
6c84d371 2203 FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusChiC","ESD_InvMass_GammaePluseMinusChiCDiff",*fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 2204
d7d7e825 2205 FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusPi0","ESD_InvMass_GammaePluseMinusPi0Diff",
6c84d371 2206 *fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 2207
d7d7e825 2208 //BackGround
a0b94e5c 2209
d7d7e825 2210 //Like Sign e+e-
2211 ElectronBackground("ESD_ENegBackground",vCurrentTLVeNeg);
2212 ElectronBackground("ESD_EPosBackground",vCurrentTLVePos);
2213 ElectronBackground("ESD_EPosENegBackground",vCurrentTLVeNeg);
2214 ElectronBackground("ESD_EPosENegBackground",vCurrentTLVePos);
a0b94e5c 2215
d7d7e825 2216 // Like Sign e+e- no JPsi
2217 ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDeNegNoJPsi));
2218 ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDePosNoJPsi));
a0b94e5c 2219
d7d7e825 2220 //Mixed Event
a0b94e5c 2221
6c84d371 2222 if( fCurrentEventPosElectronTClone->GetEntriesFast() > 0 && fCurrentEventNegElectronTClone->GetEntriesFast() > 0 && fKFReconstructedGammasCutTClone->GetEntriesFast() > 0 ){
d7d7e825 2223 FillGammaElectronInvMass("ESD_EPosENegGammaBackgroundMX","ESD_EPosENegGammaBackgroundMXDiff",
6c84d371 2224 *fKFReconstructedGammasCutTClone,*fPreviousEventTLVNegElectronTClone,*fPreviousEventTLVPosElectronTClone);
2225 *fPreviousEventTLVNegElectronTClone = vCurrentTLVeNeg;
2226 *fPreviousEventTLVPosElectronTClone = vCurrentTLVePos;
a0b94e5c 2227
d7d7e825 2228 }
a0b94e5c 2229
d7d7e825 2230 /*
2231 //Photons P
2232 Double_t vtx[3];
2233 vtx[0]=0;vtx[1]=0;vtx[2]=0;
2234 for(UInt_t i=0;i<fKFReconstructedGammasChic.size();i++){
a0b94e5c 2235
d7d7e825 2236 // if(fMCGammaChicTempCut[i]->GetMother(0) < 0) continue;
a0b94e5c 2237
2238
2239
d7d7e825 2240 Int_t tempLabel = fStack->Particle(fMCGammaChicTempCut[i]->GetMother(0))->GetPdgCode();
2241 // cout<<" Label Pedro Gonzalez " <<tempLabel <<endl;
a0b94e5c 2242
d7d7e825 2243 // cout<<" Label Distance"<<fKFReconstructedGammasChic[i].GetDistanceFromVertex(vtx)<<endl;
a0b94e5c 2244
d7d7e825 2245 if( tempLabel == 10441 || tempLabel == 20443 || tempLabel == 445 )
a0b94e5c 2246
d7d7e825 2247 fHistograms->FillHistogram("ESD_PhotonsMomentum",fKFReconstructedGammasChic[i].GetMomentum());
a0b94e5c 2248
2249
d7d7e825 2250 }
a0b94e5c 2251
2252
d7d7e825 2253 */
1e7846f4 2254
2255
2256 vESDeNegTemp->Delete();
2257 vESDePosTemp->Delete();
2258 vESDxNegTemp->Delete();
2259 vESDxPosTemp->Delete();
2260 vESDeNegNoJPsi->Delete();
2261 vESDePosNoJPsi->Delete();
2262
2263 delete vESDeNegTemp;
2264 delete vESDePosTemp;
2265 delete vESDxNegTemp;
2266 delete vESDxPosTemp;
2267 delete vESDeNegNoJPsi;
2268 delete vESDePosNoJPsi;
d7d7e825 2269}
2270
6c84d371 2271/*
2272 void AliAnalysisTaskGammaConversion::FillAngle(TString histoName,vector <TLorentzVector> tlVeNeg, vector <TLorentzVector> tlVePos){
d7d7e825 2273 //see header file for documentation
2274 for( UInt_t iNeg=0; iNeg < tlVeNeg.size(); iNeg++){
6c84d371 2275 for (UInt_t iPos=0; iPos < tlVePos.size(); iPos++){
2276 fHistograms->FillHistogram(histoName.Data(),tlVeNeg[iNeg].Vect().Angle(tlVePos[iPos].Vect()));
2277 }
2278 }
2279 }
2280*/
2281void AliAnalysisTaskGammaConversion::FillAngle(TString histoName,TClonesArray const tlVeNeg, TClonesArray const tlVePos){
2282 //see header file for documentation
2283 for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++){
2284 for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){
2285 fHistograms->FillHistogram(histoName.Data(),((TLorentzVector*)(tlVeNeg.At(iNeg)))->Vect().Angle(((TLorentzVector*)(tlVePos.At(iPos)))->Vect()));
d7d7e825 2286 }
2287 }
2288}
6c84d371 2289void AliAnalysisTaskGammaConversion::FillElectronInvMass(TString histoName, TClonesArray const eNeg, TClonesArray const ePos){
d7d7e825 2290 //see header file for documentation
6c84d371 2291 for( Int_t n=0; n < eNeg.GetEntriesFast(); n++){
2292 TLorentzVector en = (*(TLorentzVector*)(eNeg.At(n)));
2293 for (Int_t p=0; p < ePos.GetEntriesFast(); p++){
2294 TLorentzVector ep = (*(TLorentzVector*)(ePos.At(p)));
d7d7e825 2295 TLorentzVector np = ep + en;
2296 fHistograms->FillHistogram(histoName.Data(),np.M());
2297 }
2298 }
d7d7e825 2299}
2300
6c84d371 2301void AliAnalysisTaskGammaConversion::FillGammaElectronInvMass(TString histoMass,TString histoDiff,TClonesArray const fKFGammas,
2302 TClonesArray const tlVeNeg,TClonesArray const tlVePos)
d7d7e825 2303{
2304 //see header file for documentation
a0b94e5c 2305
6c84d371 2306 for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++ ){
a0b94e5c 2307
6c84d371 2308 for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){
a0b94e5c 2309
6c84d371 2310 TLorentzVector xy = *((TLorentzVector *)(tlVePos.At(iPos))) + *((TLorentzVector *)(tlVeNeg.At(iNeg)));
a0b94e5c 2311
6c84d371 2312 for (Int_t iGam=0; iGam < fKFGammas.GetEntriesFast(); iGam++){
a0b94e5c 2313
6c84d371 2314 // AliKFParticle * gammaCandidate = &fKFGammas[iGam];
2315 AliKFParticle * gammaCandidate = (AliKFParticle *)(fKFGammas.At(iGam));
d7d7e825 2316 TLorentzVector g;
a0b94e5c 2317
d7d7e825 2318 g.SetXYZM(gammaCandidate->GetPx(),gammaCandidate->GetPy(),gammaCandidate->GetPz(),fGammaMass);
2319 TLorentzVector xyg = xy + g;
2320 fHistograms->FillHistogram(histoMass.Data(),xyg.M());
2321 fHistograms->FillHistogram(histoDiff.Data(),(xyg.M()-xy.M()));
2322 }
2323 }
2324 }
a0b94e5c 2325
d7d7e825 2326}
6c84d371 2327void AliAnalysisTaskGammaConversion::ElectronBackground(TString hBg, TClonesArray e)
d7d7e825 2328{
2329 // see header file for documentation
6c84d371 2330 for(Int_t i=0; i < e.GetEntriesFast(); i++)
d7d7e825 2331 {
6c84d371 2332 for (Int_t j=i+1; j < e.GetEntriesFast(); j++)
d7d7e825 2333 {
6c84d371 2334 TLorentzVector ee = (*(TLorentzVector*)(e.At(i))) + (*(TLorentzVector*)(e.At(j)));
a0b94e5c 2335
d7d7e825 2336 fHistograms->FillHistogram(hBg.Data(),ee.M());
2337 }
2338 }
2339}
2340
2341
6c84d371 2342void AliAnalysisTaskGammaConversion::CleanWithAngleCuts(TClonesArray const negativeElectrons,
2343 TClonesArray const positiveElectrons,
2344 TClonesArray const gammas){
d7d7e825 2345 // see header file for documentation
a0b94e5c 2346
6c84d371 2347 UInt_t sizeN = negativeElectrons.GetEntriesFast();
2348 UInt_t sizeP = positiveElectrons.GetEntriesFast();
2349 UInt_t sizeG = gammas.GetEntriesFast();
a0b94e5c 2350
2351
2352
d7d7e825 2353 vector <Bool_t> xNegBand(sizeN);
2354 vector <Bool_t> xPosBand(sizeP);
2355 vector <Bool_t> gammaBand(sizeG);
a0b94e5c 2356
2357
d7d7e825 2358 for(UInt_t iNeg=0; iNeg < sizeN; iNeg++) xNegBand[iNeg]=kTRUE;
2359 for(UInt_t iPos=0; iPos < sizeP; iPos++) xPosBand[iPos]=kTRUE;
2360 for(UInt_t iGam=0; iGam < sizeG; iGam++) gammaBand[iGam]=kTRUE;
d7d7e825 2361
a0b94e5c 2362
2363 for(UInt_t iPos=0; iPos < sizeP; iPos++){
2364
6c84d371 2365 Double_t aP[3];
2366 ((AliESDtrack*)(positiveElectrons.At(iPos)))->GetConstrainedPxPyPz(aP);
a0b94e5c 2367
d7d7e825 2368 TVector3 ePosV(aP[0],aP[1],aP[2]);
a0b94e5c 2369
d7d7e825 2370 for(UInt_t iNeg=0; iNeg < sizeN; iNeg++){
a0b94e5c 2371
6c84d371 2372 Double_t aN[3];
2373 ((AliESDtrack*)(negativeElectrons.At(iNeg)))->GetConstrainedPxPyPz(aN);
d7d7e825 2374 TVector3 eNegV(aN[0],aN[1],aN[2]);
a0b94e5c 2375
d7d7e825 2376 if(ePosV.Angle(eNegV) < 0.05){ //e+e- from gamma
2377 xPosBand[iPos]=kFALSE;
2378 xNegBand[iNeg]=kFALSE;
2379 }
a0b94e5c 2380
d7d7e825 2381 for(UInt_t iGam=0; iGam < sizeG; iGam++){
6c84d371 2382 AliKFParticle* gammaCandidate = (AliKFParticle*)gammas.At(iGam);
d7d7e825 2383 TVector3 gammaCandidateVector(gammaCandidate->Px(),gammaCandidate->Py(),gammaCandidate->Pz());
2384 if(ePosV.Angle(gammaCandidateVector) < 0.05 || eNegV.Angle(gammaCandidateVector) < 0.05)
2385 gammaBand[iGam]=kFALSE;
2386 }
2387 }
2388 }
a0b94e5c 2389
2390
2391
2392
d7d7e825 2393 for(UInt_t iPos=0; iPos < sizeP; iPos++){
2394 if(xPosBand[iPos]){
6c84d371 2395 new((*fCurrentEventPosElectronTClone)[fCurrentEventPosElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(positiveElectrons.At(iPos))));
2396 // fCurrentEventPosElectron.push_back(positiveElectrons[iPos]);
d7d7e825 2397 }
2398 }
2399 for(UInt_t iNeg=0;iNeg < sizeN; iNeg++){
2400 if(xNegBand[iNeg]){
6c84d371 2401 new((*fCurrentEventNegElectronTClone)[fCurrentEventNegElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(negativeElectrons.At(iNeg))));
2402 // fCurrentEventNegElectron.push_back(negativeElectrons[iNeg]);
d7d7e825 2403 }
2404 }
2405 for(UInt_t iGam=0; iGam < sizeG; iGam++){
2406 if(gammaBand[iGam]){
6c84d371 2407 new((*fKFReconstructedGammasCutTClone)[fKFReconstructedGammasCutTClone->GetEntriesFast()]) AliKFParticle((*(AliKFParticle*)(gammas.At(iGam))));
2408 //fKFReconstructedGammasCut.push_back(*(AliKFParticle*)gammas->At(iGam));
d7d7e825 2409 }
2410 }
2411}
2412
2413
2414void AliAnalysisTaskGammaConversion::GetPID(AliESDtrack *track, Stat_t &pid, Stat_t &weight)
2415{
2416 // see header file for documentation
2417 pid = -1;
2418 weight = -1;
a0b94e5c 2419
d7d7e825 2420 double wpart[5];
2421 double wpartbayes[5];
a0b94e5c 2422
d7d7e825 2423 //get probability of the diffenrent particle types
2424 track->GetESDpid(wpart);
a0b94e5c 2425
d7d7e825 2426 // Tentative particle type "concentrations"
2427 double c[5]={0.01, 0.01, 0.85, 0.10, 0.05};
a0b94e5c 2428
d7d7e825 2429 //Bayes' formula
2430 double rcc = 0.;
2431 for (int i = 0; i < 5; i++)
2432 {
2433 rcc+=(c[i] * wpart[i]);
2434 }
a0b94e5c 2435
2436
2437
d7d7e825 2438 for (int i=0; i<5; i++) {
4a6157dc 2439 if( rcc>0 || rcc<0){//Kenneth: not sure if the rcc<0 should be there, this is from fixing a coding violation where we are not allowed to say: rcc!=0 (RC19)
d7d7e825 2440 wpartbayes[i] = c[i] * wpart[i] / rcc;
2441 }
2442 }
a0b94e5c 2443
2444
2445
d7d7e825 2446 Float_t max=0.;
2447 int ipid=-1;
2448 //find most probable particle in ESD pid
2449 //0:Electron - 1:Muon - 2:Pion - 3:Kaon - 4:Proton
2450 for (int i = 0; i < 5; i++)
2451 {
2452 if (wpartbayes[i] > max)
2453 {
a0b94e5c 2454 ipid = i;
2455 max = wpartbayes[i];
d7d7e825 2456 }
2457 }
a0b94e5c 2458
d7d7e825 2459 pid = ipid;
2460 weight = max;
2461}
2462double AliAnalysisTaskGammaConversion::GetSigmaToVertex(AliESDtrack* t)
2463{
2464 // Calculates the number of sigma to the vertex.
a0b94e5c 2465
d7d7e825 2466 Float_t b[2];
2467 Float_t bRes[2];
2468 Float_t bCov[3];
2469 t->GetImpactParameters(b,bCov);
2470 if (bCov[0]<=0 || bCov[2]<=0) {
2471 AliDebug(1, "Estimated b resolution lower or equal zero!");
2472 bCov[0]=0; bCov[2]=0;
2473 }
2474 bRes[0] = TMath::Sqrt(bCov[0]);
2475 bRes[1] = TMath::Sqrt(bCov[2]);
a0b94e5c 2476
d7d7e825 2477 // -----------------------------------
2478 // How to get to a n-sigma cut?
2479 //
2480 // The accumulated statistics from 0 to d is
2481 //
2482 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
2483 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
2484 //
2485 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
2486 // Can this be expressed in a different way?
a0b94e5c 2487
d7d7e825 2488 if (bRes[0] == 0 || bRes[1] ==0)
2489 return -1;
a0b94e5c 2490
d7d7e825 2491 double d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
a0b94e5c 2492
d7d7e825 2493 // stupid rounding problem screws up everything:
2494 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
2495 if (TMath::Exp(-d * d / 2) < 1e-10)
2496 return 1000;
a0b94e5c 2497
2498
d7d7e825 2499 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
2500 return d;
2501}
6c84d371 2502
2503//vector <TLorentzVector> AliAnalysisTaskGammaConversion::GetTLorentzVector(vector <AliESDtrack*> esdTrack){
2504TClonesArray AliAnalysisTaskGammaConversion::GetTLorentzVector(TClonesArray *const esdTrack){
01b7fdcc 2505 //Return TLoresntz vector of track?
6c84d371 2506 // vector <TLorentzVector> tlVtrack(0);
2507 TClonesArray array("TLorentzVector",0);
a0b94e5c 2508
6c84d371 2509 for(Int_t itrack=0; itrack < esdTrack->GetEntriesFast(); itrack++){
2510 double p[3];
2511 //esdTrack[itrack]->GetConstrainedPxPyPz(p);
2512 ((AliESDtrack*)(esdTrack->At(itrack)))->GetConstrainedPxPyPz(p);
d7d7e825 2513 TLorentzVector currentTrack;
01b7fdcc 2514 currentTrack.SetXYZM(p[0],p[1],p[2],fElectronMass);
6c84d371 2515 new((array)[array.GetEntriesFast()]) TLorentzVector(currentTrack);
2516 // tlVtrack.push_back(currentTrack);
d7d7e825 2517 }
a0b94e5c 2518
6c84d371 2519 return array;
a0b94e5c 2520
6c84d371 2521 // return tlVtrack;
d7d7e825 2522}
2523