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