]> 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
coverity 15108, the code needs further revision since the loop does not make sense...
[u/mrichter/AliRoot.git] / PWG4 / GammaConv / AliAnalysisTaskGammaConversion.cxx
CommitLineData
d7d7e825 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *
5 * Version 1.1 *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16////////////////////////////////////////////////
17//---------------------------------------------
18// Class used to do analysis on conversion pairs
19//---------------------------------------------
20////////////////////////////////////////////////
21
22// root
23#include <TChain.h>
24
25// analysis
26#include "AliAnalysisTaskGammaConversion.h"
27#include "AliStack.h"
28#include "AliLog.h"
29#include "AliESDtrackCuts.h"
30#include "TNtuple.h"
4a6157dc 31//#include "AliCFManager.h" // for CF
32//#include "AliCFContainer.h" // for CF
9c1cb6f7 33#include "AliESDInputHandler.h"
34#include "AliAnalysisManager.h"
6441e967 35#include "AliAODPWG4Particle.h"
36#include "AliAODPWG4ParticleCorrelation.h"
4a6157dc 37#include "AliGammaConversionAODObject.h"
6441e967 38#include "AliAODConversionParticle.h"
5e55d806 39#include "AliGammaConversionBGHandler.h"
6272370b 40#include "AliESDCaloCluster.h" // for combining PHOS and GammaConv
9c1cb6f7 41#include "AliKFVertex.h"
3c45d101 42#include "AliGenPythiaEventHeader.h"
43#include "AliGenDPMjetEventHeader.h"
44#include "AliGenEventHeader.h"
e60f3265 45#include <AliMCEventHandler.h>
77ac6f3e 46#include "TRandom3.h"
18a1e325 47#include "AliTriggerAnalysis.h"
b6dd6ad2 48#include "AliCentrality.h"
77ac6f3e 49
d707e3cf 50class AliESDTrackCuts;
4a6157dc 51class AliCFContainer;
52class AliCFManager;
d7d7e825 53class AliKFVertex;
54class AliAODHandler;
55class AliAODEvent;
56class ALiESDEvent;
57class AliMCEvent;
58class AliMCEventHandler;
59class AliESDInputHandler;
60class AliAnalysisManager;
61class Riostream;
62class TFile;
63class TInterpreter;
64class TSystem;
65class TROOT;
66
67ClassImp(AliAnalysisTaskGammaConversion)
68
69
70AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion():
71AliAnalysisTaskSE(),
72 fV0Reader(NULL),
73 fStack(NULL),
a0b94e5c 74 fMCTruth(NULL), // for CF
4a6157dc 75 fGCMCEvent(NULL), // for CF
d7d7e825 76 fESDEvent(NULL),
77 fOutputContainer(NULL),
a0b94e5c 78 fCFManager(0x0), // for CF
d7d7e825 79 fHistograms(NULL),
b5832f95 80 fTriggerCINT1B(kFALSE),
d7d7e825 81 fDoMCTruth(kFALSE),
82 fDoNeutralMeson(kFALSE),
6272370b 83 fDoOmegaMeson(kFALSE),
d7d7e825 84 fDoJet(kFALSE),
85 fDoChic(kFALSE),
9c1cb6f7 86 fRecalculateV0ForGamma(kFALSE),
6c84d371 87 fKFReconstructedGammasTClone(NULL),
6272370b 88 fKFReconstructedPi0sTClone(NULL),
9c1cb6f7 89 fKFRecalculatedGammasTClone(NULL),
6c84d371 90 fCurrentEventPosElectronTClone(NULL),
91 fCurrentEventNegElectronTClone(NULL),
92 fKFReconstructedGammasCutTClone(NULL),
93 fPreviousEventTLVNegElectronTClone(NULL),
94 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 95 fElectronv1(),
96 fElectronv2(),
6272370b 97 fGammav1(),
98 fGammav2(),
9c1cb6f7 99 fElectronRecalculatedv1(),
100 fElectronRecalculatedv2(),
d7d7e825 101 fElectronMass(-1),
102 fGammaMass(-1),
103 fPi0Mass(-1),
104 fEtaMass(-1),
105 fGammaWidth(-1),
106 fPi0Width(-1),
107 fEtaWidth(-1),
108 fMinOpeningAngleGhostCut(0.),
9640a3d1 109 fEsdTrackCuts(NULL),
d7d7e825 110 fCalculateBackground(kFALSE),
111 fWriteNtuple(kFALSE),
112 fGammaNtuple(NULL),
113 fNeutralMesonNtuple(NULL),
114 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 115 fChargedParticles(NULL),
d7d7e825 116 fChargedParticlesId(),
117 fGammaPtHighest(0.),
118 fMinPtForGammaJet(1.),
119 fMinIsoConeSize(0.2),
120 fMinPtIsoCone(0.7),
121 fMinPtGamChargedCorr(0.5),
122 fMinPtJetCone(0.5),
123 fLeadingChargedIndex(-1),
9640a3d1 124 fLowPtMapping(1.),
125 fHighPtMapping(3.),
1e7846f4 126 fDoCF(kFALSE),
04bf4381 127 fAODGamma(NULL),
128 fAODPi0(NULL),
129 fAODOmega(NULL),
037dc2db 130 fAODBranchName("GammaConv"),
63477d82 131 fKFCreateAOD(kTRUE),
d765d400 132 fKFForceAOD(kFALSE),
332f1f44 133 fKFDeltaAODFileName(""),
037dc2db 134 fDoNeutralMesonV0MCCheck(kFALSE),
5ce758b0 135 fUseTrackMultiplicityForBG(kTRUE),
136 fMoveParticleAccordingToVertex(kFALSE),
2e2da371 137 fApplyChi2Cut(kFALSE),
77ac6f3e 138 fNRandomEventsForBG(15),
139 fNDegreesPMBackground(15),
111d75df 140 fDoRotation(kTRUE),
141 fCheckBGProbability(kTRUE),
2f426a53 142 fKFReconstructedGammasV0Index(),
18a1e325 143 fRemovePileUp(kFALSE),
144 fSelectV0AND(kFALSE),
6746e1e1 145 fTriggerAnalysis(NULL),
146 fMultiplicity(0),
147 fUseMultiplicity(0),
64f4118c 148 fUseMultiplicityBin(0),
149 fUseCentrality(0),
2ea0dd4a 150 fUseCentralityBin(0),
151 fRandom(0)
d7d7e825 152{
153 // Default constructor
5a34881d 154
155 /* Kenneth: the default constructor should not have any define input/output or the call to SetESDtrackCuts
d7d7e825 156 // Common I/O in slot 0
157 DefineInput (0, TChain::Class());
158 DefineOutput(0, TTree::Class());
159
160 // Your private output
161 DefineOutput(1, TList::Class());
a0b94e5c 162
d7d7e825 163 // Define standard ESD track cuts for Gamma-hadron correlation
164 SetESDtrackCuts();
5a34881d 165 */
d7d7e825 166}
167
168AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name):
169 AliAnalysisTaskSE(name),
170 fV0Reader(NULL),
171 fStack(NULL),
a0b94e5c 172 fMCTruth(NULL), // for CF
4a6157dc 173 fGCMCEvent(NULL), // for CF
d7d7e825 174 fESDEvent(NULL),
175 fOutputContainer(0x0),
a0b94e5c 176 fCFManager(0x0), // for CF
d7d7e825 177 fHistograms(NULL),
b5832f95 178 fTriggerCINT1B(kFALSE),
d7d7e825 179 fDoMCTruth(kFALSE),
180 fDoNeutralMeson(kFALSE),
6272370b 181 fDoOmegaMeson(kFALSE),
d7d7e825 182 fDoJet(kFALSE),
183 fDoChic(kFALSE),
9c1cb6f7 184 fRecalculateV0ForGamma(kFALSE),
6c84d371 185 fKFReconstructedGammasTClone(NULL),
6272370b 186 fKFReconstructedPi0sTClone(NULL),
9c1cb6f7 187 fKFRecalculatedGammasTClone(NULL),
6c84d371 188 fCurrentEventPosElectronTClone(NULL),
189 fCurrentEventNegElectronTClone(NULL),
190 fKFReconstructedGammasCutTClone(NULL),
191 fPreviousEventTLVNegElectronTClone(NULL),
192 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 193 fElectronv1(),
194 fElectronv2(),
6272370b 195 fGammav1(),
196 fGammav2(),
9c1cb6f7 197 fElectronRecalculatedv1(),
198 fElectronRecalculatedv2(),
d7d7e825 199 fElectronMass(-1),
200 fGammaMass(-1),
201 fPi0Mass(-1),
202 fEtaMass(-1),
203 fGammaWidth(-1),
204 fPi0Width(-1),
205 fEtaWidth(-1),
206 fMinOpeningAngleGhostCut(0.),
9640a3d1 207 fEsdTrackCuts(NULL),
d7d7e825 208 fCalculateBackground(kFALSE),
209 fWriteNtuple(kFALSE),
210 fGammaNtuple(NULL),
211 fNeutralMesonNtuple(NULL),
212 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 213 fChargedParticles(NULL),
d7d7e825 214 fChargedParticlesId(),
215 fGammaPtHighest(0.),
216 fMinPtForGammaJet(1.),
217 fMinIsoConeSize(0.2),
218 fMinPtIsoCone(0.7),
219 fMinPtGamChargedCorr(0.5),
220 fMinPtJetCone(0.5),
221 fLeadingChargedIndex(-1),
9640a3d1 222 fLowPtMapping(1.),
223 fHighPtMapping(3.),
c59360eb 224 fDoCF(kFALSE),
225 fAODGamma(NULL),
226 fAODPi0(NULL),
227 fAODOmega(NULL),
228 fAODBranchName("GammaConv"),
63477d82 229 fKFCreateAOD(kTRUE),
c59360eb 230 fKFForceAOD(kFALSE),
231 fKFDeltaAODFileName(""),
232 fDoNeutralMesonV0MCCheck(kFALSE),
5ce758b0 233 fUseTrackMultiplicityForBG(kTRUE),
234 fMoveParticleAccordingToVertex(kFALSE),
2e2da371 235 fApplyChi2Cut(kFALSE),
77ac6f3e 236 fNRandomEventsForBG(15),
237 fNDegreesPMBackground(15),
111d75df 238 fDoRotation(kTRUE),
239 fCheckBGProbability(kTRUE),
2f426a53 240 fKFReconstructedGammasV0Index(),
18a1e325 241 fRemovePileUp(kFALSE),
242 fSelectV0AND(kFALSE),
6746e1e1 243 fTriggerAnalysis(NULL),
244 fMultiplicity(0),
245 fUseMultiplicity(0),
64f4118c 246 fUseMultiplicityBin(0),
247 fUseCentrality(0),
2ea0dd4a 248 fUseCentralityBin(0),
249 fRandom(0)
d7d7e825 250{
143874ce 251 // Common I/O in slot 0, don't define when inheriting from AnalysisTaskSE
252 // DefineInput (0, TChain::Class());
253 // DefineOutput(0, TTree::Class());
d7d7e825 254
255 // Your private output
256 DefineOutput(1, TList::Class());
a0b94e5c 257 DefineOutput(2, AliCFContainer::Class()); // for CF
258
259
d7d7e825 260 // Define standard ESD track cuts for Gamma-hadron correlation
261 SetESDtrackCuts();
9c1cb6f7 262
d7d7e825 263}
264
265AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion()
266{
267 // Remove all pointers
268
269 if(fOutputContainer){
270 fOutputContainer->Clear() ;
271 delete fOutputContainer ;
272 }
273 if(fHistograms){
274 delete fHistograms;
275 }
276 if(fV0Reader){
277 delete fV0Reader;
278 }
a0b94e5c 279
280 // for CF
281 if(fCFManager){
282 delete fCFManager;
283 }
9640a3d1 284
285 if(fEsdTrackCuts){
286 delete fEsdTrackCuts;
287 }
288
04bf4381 289 //Delete AODs
290 if (fAODGamma) {
291 fAODGamma->Clear();
292 delete fAODGamma;
d7d7e825 293 }
04bf4381 294 fAODGamma = NULL;
295
296 if (fAODPi0) {
297 fAODPi0->Clear();
298 delete fAODPi0;
299 }
300 fAODPi0 = NULL;
d7d7e825 301
04bf4381 302 if (fAODOmega) {
303 fAODOmega->Clear();
304 delete fAODOmega;
305 }
306 fAODOmega = NULL;
307
18a1e325 308 if(fTriggerAnalysis) {
309 delete fTriggerAnalysis;
310 }
311
312
04bf4381 313}
d7d7e825 314
12464034 315
d7d7e825 316void AliAnalysisTaskGammaConversion::Init()
317{
318 // Initialization
4a6157dc 319 // AliLog::SetGlobalLogLevel(AliLog::kError);
d7d7e825 320}
321void AliAnalysisTaskGammaConversion::SetESDtrackCuts()
322{
323 // SetESDtrackCuts
9640a3d1 324 if (fEsdTrackCuts!=NULL){
325 delete fEsdTrackCuts;
326 }
d7d7e825 327 fEsdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
a0b94e5c 328 //standard cuts from:
329 //http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup
ebcfaa7e 330
331 // Cuts used up to 3rd of March
332
333 // fEsdTrackCuts->SetMinNClustersTPC(50);
334// fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
335// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
336// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
337// fEsdTrackCuts->SetMaxNsigmaToVertex(3);
338// fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE);
339
340 //------- To be tested-----------
d707e3cf 341 // Cuts used up to 26th of Agost
342// Int_t minNClustersTPC = 70;
343// Double_t maxChi2PerClusterTPC = 4.0;
344// Double_t maxDCAtoVertexXY = 2.4; // cm
345// Double_t maxDCAtoVertexZ = 3.2; // cm
346// fEsdTrackCuts->SetRequireSigmaToVertex(kFALSE);
347// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
348// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
349// // fEsdTrackCuts->SetRequireTPCStandAlone(kTRUE);
350// fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
351// fEsdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
352// fEsdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
353// fEsdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
354// fEsdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
355// fEsdTrackCuts->SetDCAToVertex2D(kTRUE);
356// fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
357// fEsdTrackCuts->SetPtRange(0.15);
358
359// fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
360
361
362// Using standard function for setting Cuts
363 Bool_t selectPrimaries=kTRUE;
406a77e0 364 fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(selectPrimaries);
64f4118c 365 fEsdTrackCuts->SetMaxDCAToVertexZ(2);
d707e3cf 366 fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
367 fEsdTrackCuts->SetPtRange(0.15);
ebcfaa7e 368
369 //----- From Jacek 10.03.03 ------------------/
370// minNClustersTPC = 70;
371// maxChi2PerClusterTPC = 4.0;
372// maxDCAtoVertexXY = 2.4; // cm
373// maxDCAtoVertexZ = 3.2; // cm
374
375// esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
376// esdTrackCuts->SetRequireTPCRefit(kFALSE);
377// esdTrackCuts->SetRequireTPCStandAlone(kTRUE);
378// esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
379// esdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
380// esdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
381// esdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
382// esdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
383// esdTrackCuts->SetDCAToVertex2D(kTRUE);
384
385
386
d7d7e825 387 // fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
037dc2db 388 // fV0Reader->SetESDtrackCuts(fEsdTrackCuts);
d7d7e825 389}
390
c00009fb 391void AliAnalysisTaskGammaConversion::UserExec(Option_t */*option*/)
d7d7e825 392{
393 // Execute analysis for current event
6272370b 394
9c1cb6f7 395 // Load the esdpid from the esdhandler if exists (tender was applied) otherwise set the Bethe Bloch parameters
3c45d101 396 Int_t eventQuality=-1;
9c1cb6f7 397
398 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
399 AliESDInputHandler *esdHandler=0x0;
400 if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){
401 AliV0Reader::SetESDpid(esdHandler->GetESDpid());
402 } else {
403 //load esd pid bethe bloch parameters depending on the existance of the MC handler
404 // yes: MC parameters
405 // no: data parameters
406 if (!AliV0Reader::GetESDpid()){
407 if (fMCEvent ) {
408 AliV0Reader::InitESDpid();
409 } else {
410 AliV0Reader::InitESDpid(1);
411 }
412 }
413 }
414
d765d400 415 //Must set fForceAOD to true for the AOD to get filled. Should only be done when running independent chain / train.
416 if(fKFForceAOD) {
417 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
418 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
419 }
48682642 420
239a56d1 421 // if(fV0Reader == NULL){ // coverty does not permit this test
48682642 422 // Write warning here cuts and so on are default if this ever happens
239a56d1 423 // }
48682642 424
e60f3265 425 if (fMCEvent ) {
426 // To avoid crashes due to unzip errors. Sometimes the trees are not there.
427
428 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
3c45d101 429 if (!mcHandler){
430 AliError("Could not retrive MC event handler!");
431 return;
432
433 eventQuality=0;
434 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
435 }
436 if (!mcHandler->InitOk() ){
437 eventQuality=0;
438 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
439 return;
440 }
441 if (!mcHandler->TreeK() ){
442 eventQuality=0;
443 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
444 return;
445 }
446 if (!mcHandler->TreeTR() ) {
447 eventQuality=0;
448 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
449 return;
450 }
e60f3265 451 }
452
48682642 453 fV0Reader->SetInputAndMCEvent(InputEvent(), MCEvent());
454
455 fV0Reader->Initialize();
456 fDoMCTruth = fV0Reader->GetDoMCTruth();
457
12464034 458 if(fAODGamma) fAODGamma->Delete();
459 if(fAODPi0) fAODPi0->Delete();
460 if(fAODOmega) fAODOmega->Delete();
a0b94e5c 461
6c84d371 462 if(fKFReconstructedGammasTClone == NULL){
463 fKFReconstructedGammasTClone = new TClonesArray("AliKFParticle",0);
464 }
465 if(fCurrentEventPosElectronTClone == NULL){
466 fCurrentEventPosElectronTClone = new TClonesArray("AliESDtrack",0);
467 }
468 if(fCurrentEventNegElectronTClone == NULL){
469 fCurrentEventNegElectronTClone = new TClonesArray("AliESDtrack",0);
470 }
471 if(fKFReconstructedGammasCutTClone == NULL){
472 fKFReconstructedGammasCutTClone = new TClonesArray("AliKFParticle",0);
473 }
474 if(fPreviousEventTLVNegElectronTClone == NULL){
475 fPreviousEventTLVNegElectronTClone = new TClonesArray("TLorentzVector",0);
476 }
477 if(fPreviousEventTLVPosElectronTClone == NULL){
478 fPreviousEventTLVPosElectronTClone = new TClonesArray("TLorentzVector",0);
479 }
480 if(fChargedParticles == NULL){
481 fChargedParticles = new TClonesArray("AliESDtrack",0);
482 }
6272370b 483
484 if(fKFReconstructedPi0sTClone == NULL){
485 fKFReconstructedPi0sTClone = new TClonesArray("AliKFParticle",0);
486 }
9c1cb6f7 487
18a1e325 488 if(fKFRecalculatedGammasTClone == NULL){
9c1cb6f7 489 fKFRecalculatedGammasTClone = new TClonesArray("AliKFParticle",0);
490 }
491
18a1e325 492 if(fTriggerAnalysis== NULL){
493 fTriggerAnalysis = new AliTriggerAnalysis;
494 }
495
6c84d371 496 //clear TClones
1e7846f4 497 fKFReconstructedGammasTClone->Delete();
498 fCurrentEventPosElectronTClone->Delete();
499 fCurrentEventNegElectronTClone->Delete();
500 fKFReconstructedGammasCutTClone->Delete();
501 fPreviousEventTLVNegElectronTClone->Delete();
502 fPreviousEventTLVPosElectronTClone->Delete();
6272370b 503 fKFReconstructedPi0sTClone->Delete();
9c1cb6f7 504 fKFRecalculatedGammasTClone->Delete();
5e55d806 505
d7d7e825 506 //clear vectors
d7d7e825 507 fElectronv1.clear();
508 fElectronv2.clear();
6272370b 509
510 fGammav1.clear();
511 fGammav2.clear();
512
9c1cb6f7 513 fElectronRecalculatedv1.clear();
514 fElectronRecalculatedv2.clear();
515
d7d7e825 516
1e7846f4 517 fChargedParticles->Delete();
5e55d806 518
d7d7e825 519 fChargedParticlesId.clear();
037dc2db 520
521 fKFReconstructedGammasV0Index.clear();
a0b94e5c 522
d7d7e825 523 //Clear the data in the v0Reader
5e55d806 524 // fV0Reader->UpdateEventByEventData();
b5832f95 525
526 //Take Only events with proper trigger
5e55d806 527 /*
b5832f95 528 if(fTriggerCINT1B){
529 if(!fV0Reader->GetESDEvent()->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL")) return;
530 }
5e55d806 531 */
22182c37 532 Bool_t v0A = fTriggerAnalysis->IsOfflineTriggerFired(fV0Reader->GetESDEvent(), AliTriggerAnalysis::kV0A);
533 Bool_t v0C = fTriggerAnalysis->IsOfflineTriggerFired(fV0Reader->GetESDEvent(), AliTriggerAnalysis::kV0C);
534 Bool_t v0AND = v0A && v0C;
535
536 if(fSelectV0AND && !v0AND){
537 eventQuality=5;
538 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
539 if(fDoMCTruth){
540 CheckMesonProcessTypeEventQuality(eventQuality);
541 }
542
543 return;
544 }
10e3319b 545
c8206114 546 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
547 // cout<< "Event not taken"<< endl;
3c45d101 548 eventQuality=1;
549 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
550 if(fDoMCTruth){
551 CheckMesonProcessTypeEventQuality(eventQuality);
552 }
c8206114 553 return; // aborts if the primary vertex does not have contributors.
554 }
555
ada7a878 556
c8206114 557
885114d1 558 if(!fV0Reader->CheckForPrimaryVertexZ() ){
3c45d101 559 eventQuality=2;
560 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
561 if(fDoMCTruth){
562 CheckMesonProcessTypeEventQuality(eventQuality);
563 }
885114d1 564 return;
565 }
2f426a53 566
ada7a878 567 if(fV0Reader->GetESDEvent()->GetPrimaryVertexTracks()->GetNContributors()>0) {
568 fHistograms->FillHistogram("ESD_GlobalPrimaryVtxZ",fV0Reader->GetESDEvent()->GetPrimaryVertex()->GetZ());
569 }else{
570 if(fV0Reader->GetESDEvent()->GetPrimaryVertexSPD()->GetNContributors()>0) {
571 fHistograms->FillHistogram("ESD_SPDPrimaryVtxZ",fV0Reader->GetESDEvent()->GetPrimaryVertex()->GetZ());
572 }
573 }
574
2f426a53 575 if(fRemovePileUp && fV0Reader->GetESDEvent()->IsPileupFromSPD()) {
0a3d9f49 576 eventQuality=4;
45441b78 577 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
2f426a53 578 return;
579 }
18a1e325 580
6746e1e1 581 fMultiplicity = fEsdTrackCuts->CountAcceptedTracks(fV0Reader->GetESDEvent());
582
583 if( CalculateMultiplicityBin() != fUseMultiplicityBin){
584 eventQuality=6;
585 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
586 return;
587 }
18a1e325 588
64f4118c 589 if(fV0Reader->GetIsHeavyIon()){
590 if(fUseCentrality>0){
b6dd6ad2 591 AliCentrality *esdCentrality = fV0Reader->GetESDEvent()->GetCentrality();
64f4118c 592 Int_t centralityC = -1;
593
594 if(fUseCentrality==1){
595 centralityC = esdCentrality->GetCentralityClass10("V0M");
596 if( centralityC != fUseCentralityBin ){
597 eventQuality=7;
598 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
599 return;
600 }
601 }
602
603 if(fUseCentrality==2){
604 centralityC = esdCentrality->GetCentralityClass10("CL1");
605 if( centralityC != fUseCentralityBin ){
606 eventQuality=7;
607 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
608 return;
609 }
610 }
14d7f2f9 611
612////////////////////////////////////// RRnew start /////////////////////////////////////////////////////
613 if(fUseCentrality==3){
614 centralityC = esdCentrality->GetCentralityClass10("V0M");
615 if( (fUseCentralityBin == 0) && (centralityC!=0) ){ // 0-10%
616 eventQuality=7;
617 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
618 return;
619 }
620 if( (fUseCentralityBin == 1) && (centralityC!=1) ){ // 10-20%
621 eventQuality=7;
622 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
623 return;
624 }
625 if( (fUseCentralityBin == 2) && (centralityC!=2) && (centralityC!=3) ){ // 20-40%
626 eventQuality=7;
627 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
628 return;
629 }
630 if( (fUseCentralityBin == 4) && (centralityC!=4) && (centralityC!=5) ){ // 40-60%
631 eventQuality=7;
632 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
633 return;
634 }
635 if( (fUseCentralityBin == 6) && (centralityC!=6) && (centralityC!=7) && (centralityC!=8) ){ // 60-90%
636 eventQuality=7;
637 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
638 return;
639 }
640 }
641
642 if(fUseCentrality==4){
643 centralityC = esdCentrality->GetCentralityClass10("CL1");
644 if( (fUseCentralityBin == 0) && (centralityC!=0) ){ // 0-10%
645 eventQuality=7;
646 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
647 return;
648 }
649 if( (fUseCentralityBin == 1) && (centralityC!=1) ){ // 10-20%
650 eventQuality=7;
651 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
652 return;
653 }
654 if( (fUseCentralityBin == 2) && (centralityC!=2) && (centralityC!=3) ){ // 20-40%
655 eventQuality=7;
656 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
657 return;
658 }
659 if( (fUseCentralityBin == 4) && (centralityC!=4) && (centralityC!=5) ){ // 40-60%
660 eventQuality=7;
661 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
662 return;
663 }
664 if( (fUseCentralityBin == 6) && (centralityC!=6) && (centralityC!=7) && (centralityC!=8) ){ // 60-90%
665 eventQuality=7;
666 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
667 return;
668 }
669 }
670////////////////////////////////////// RRnew end ///////////////////////////////////////////////////////
671
64f4118c 672 }
673 }
3c45d101 674 eventQuality=3;
675 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
885114d1 676
6746e1e1 677
678
679 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",fMultiplicity);
680 if (fV0Reader->GetNumberOfContributorsVtx()>=1){
681 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",fMultiplicity);
682 }
683
684
685
d7d7e825 686 // Process the MC information
687 if(fDoMCTruth){
688 ProcessMCData();
689 }
a0b94e5c 690
d7d7e825 691 //Process the v0 information with no cuts
692 ProcessV0sNoCut();
cb90a330 693
d7d7e825 694 // Process the v0 information
695 ProcessV0s();
cb90a330 696
5e55d806 697
d7d7e825 698 //Fill Gamma AOD
63477d82 699 if(fKFCreateAOD) {
700 FillAODWithConversionGammas() ;
701 }
702
703
d7d7e825 704 // Process reconstructed gammas
705 if(fDoNeutralMeson == kTRUE){
a0b94e5c 706 ProcessGammasForNeutralMesonAnalysis();
77880bd8 707
48682642 708 }
709
77880bd8 710 if(fDoMCTruth == kTRUE){
711 CheckV0Efficiency();
712 }
d7d7e825 713 //Process reconstructed gammas electrons for Chi_c Analysis
6c84d371 714 if(fDoChic == kTRUE){
d7d7e825 715 ProcessGammaElectronsForChicAnalysis();
716 }
717 // Process reconstructed gammas for gamma Jet/hadron correlations
718 if(fDoJet == kTRUE){
719 ProcessGammasForGammaJetAnalysis();
720 }
a0b94e5c 721
5e55d806 722 //calculate background if flag is set
723 if(fCalculateBackground){
724 CalculateBackground();
725 }
48682642 726
727 if(fDoNeutralMeson == kTRUE){
728// ProcessConvPHOSGammasForNeutralMesonAnalysis();
729 if(fDoOmegaMeson == kTRUE){
730 ProcessGammasForOmegaMesonAnalysis();
731 }
732 }
733
5e55d806 734 //Clear the data in the v0Reader
735 fV0Reader->UpdateEventByEventData();
9c1cb6f7 736 if(fRecalculateV0ForGamma==kTRUE){
737 RecalculateV0ForGamma();
738 }
d7d7e825 739 PostData(1, fOutputContainer);
a0b94e5c 740 PostData(2, fCFManager->GetParticleContainer()); // for CF
d7d7e825 741
742}
743
48682642 744// void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t *option){
745// // see header file for documentation
746// // printf(" ConnectInputData %s\n", GetName());
8a685cf3 747
48682642 748// AliAnalysisTaskSE::ConnectInputData(option);
8a685cf3 749
48682642 750// if(fV0Reader == NULL){
751// // Write warning here cuts and so on are default if this ever happens
752// }
753// fV0Reader->Initialize();
754// fDoMCTruth = fV0Reader->GetDoMCTruth();
755// }
d7d7e825 756
3c45d101 757void AliAnalysisTaskGammaConversion::CheckMesonProcessTypeEventQuality(Int_t evtQ){
77ac6f3e 758 // Check meson process type event quality
3c45d101 759 fStack= MCEvent()->Stack();
760 fGCMCEvent=MCEvent();
761
762 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
763 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
764 if (!particle) {
765 //print warning here
766 continue;
767 }
768 if(particle->GetPdgCode()!=111){ //Pi0
769 continue;
770 }
62ae0124 771
772 Double_t rapidity;
773 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
774 rapidity=0;
775 }
776 else{
777 rapidity = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())));
778 }
779
780 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ) continue;
f3026589 781 if(fV0Reader->GetIsHeavyIon()) continue;
62ae0124 782
3c45d101 783 if(evtQ==1){
784 switch(GetProcessType(fGCMCEvent)){
785 case kProcSD:
786 fHistograms->FillHistogram("MC_SD_EvtQ1_Pi0_Pt", particle->Pt());
787 break;
788 case kProcDD:
789 fHistograms->FillHistogram("MC_DD_EvtQ1_Pi0_Pt", particle->Pt());
790 break;
791 case kProcND:
792 fHistograms->FillHistogram("MC_ND_EvtQ1_Pi0_Pt", particle->Pt());
793 break;
794 default:
795 AliError("Unknown Process");
796 }
797 }
798 if(evtQ==2){
799 switch(GetProcessType(fGCMCEvent)){
800 case kProcSD:
801 fHistograms->FillHistogram("MC_SD_EvtQ2_Pi0_Pt", particle->Pt());
802 break;
803 case kProcDD:
804 fHistograms->FillHistogram("MC_DD_EvtQ2_Pi0_Pt", particle->Pt());
805 break;
806 case kProcND:
807 fHistograms->FillHistogram("MC_ND_EvtQ2_Pi0_Pt", particle->Pt());
808 break;
809 default:
810 AliError("Unknown Process");
811 }
812 }
d7d7e825 813
62ae0124 814 if(evtQ==4){
815 switch(GetProcessType(fGCMCEvent)){
816 case kProcSD:
817 fHistograms->FillHistogram("MC_SD_EvtQ4_Pi0_Pt", particle->Pt());
818 break;
819 case kProcDD:
820 fHistograms->FillHistogram("MC_DD_EvtQ4_Pi0_Pt", particle->Pt());
821 break;
822 case kProcND:
823 fHistograms->FillHistogram("MC_ND_EvtQ4_Pi0_Pt", particle->Pt());
824 break;
825 default:
826 AliError("Unknown Process");
827 }
828 }
829
830 if(evtQ==5){
831 switch(GetProcessType(fGCMCEvent)){
832 case kProcSD:
833 fHistograms->FillHistogram("MC_SD_EvtQ5_Pi0_Pt", particle->Pt());
834 break;
835 case kProcDD:
836 fHistograms->FillHistogram("MC_DD_EvtQ5_Pi0_Pt", particle->Pt());
837 break;
838 case kProcND:
839 fHistograms->FillHistogram("MC_ND_EvtQ5_Pi0_Pt", particle->Pt());
840 break;
841 default:
842 AliError("Unknown Process");
843 }
844 }
d7d7e825 845
3c45d101 846 }
847
848}
849
d7d7e825 850void AliAnalysisTaskGammaConversion::ProcessMCData(){
851 // see header file for documentation
48682642 852 //InputEvent(), MCEvent());
853 /* TestAnaMarin
d7d7e825 854 fStack = fV0Reader->GetMCStack();
a0b94e5c 855 fMCTruth = fV0Reader->GetMCTruth(); // for CF
4a6157dc 856 fGCMCEvent = fV0Reader->GetMCEvent(); // for CF
48682642 857 */
858 fStack= MCEvent()->Stack();
859 fGCMCEvent=MCEvent();
a0b94e5c 860
861 // for CF
1e7846f4 862 Double_t containerInput[3];
863 if(fDoCF){
864 if(!fGCMCEvent) cout << "NO MC INFO FOUND" << endl;
865 fCFManager->SetEventInfo(fGCMCEvent);
866 }
a0b94e5c 867 // end for CF
868
869
d7d7e825 870 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
871 return; // aborts if the primary vertex does not have contributors.
872 }
2f426a53 873 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
874 // for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
d7d7e825 875 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
bd6d9fa3 876
111d75df 877
bd6d9fa3 878
d7d7e825 879 if (!particle) {
880 //print warning here
881 continue;
882 }
a0b94e5c 883
d7d7e825 884 ///////////////////////Begin Chic Analysis/////////////////////////////
91b98865 885 if(fDoChic) {
886 if(particle->GetPdgCode() == 443){//Is JPsi
887 if(particle->GetNDaughters()==2){
888 if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 &&
889 TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){
890
891 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
892 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
893 if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9)
894 fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance
a0b94e5c 895
91b98865 896 if( TMath::Abs(daug0->Eta()) < 0.9){
897 if(daug0->GetPdgCode() == -11)
898 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
899 else
900 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
a0b94e5c 901
91b98865 902 }
903 if(TMath::Abs(daug1->Eta()) < 0.9){
904 if(daug1->GetPdgCode() == -11)
905 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
906 else
907 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
908 }
d7d7e825 909 }
910 }
911 }
91b98865 912 // const int CHI_C0 = 10441;
913 // const int CHI_C1 = 20443;
914 // const int CHI_C2 = 445
915 if(particle->GetPdgCode() == 22){//gamma from JPsi
916 if(particle->GetMother(0) > -1){
917 if(fStack->Particle(particle->GetMother(0))->GetPdgCode() == 10441 ||
918 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 20443 ||
919 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 445){
920 if(TMath::Abs(particle->Eta()) < 1.2)
921 fHistograms->FillTable("Table_Electrons",17);// gamma from chic inside accptance
922 }
d7d7e825 923 }
924 }
91b98865 925 if(particle->GetPdgCode() == 10441 || particle->GetPdgCode() == 20443 || particle->GetPdgCode() == 445){
926 if( particle->GetNDaughters() == 2){
927 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
928 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
a0b94e5c 929
91b98865 930 if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){
931 if( daug0->GetPdgCode() == 443){
932 TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter());
933 TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter());
934 if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
935 fHistograms->FillTable("Table_Electrons",18);
a0b94e5c 936
91b98865 937 }//if
938 else if (daug1->GetPdgCode() == 443){
939 TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter());
940 TParticle* daugE1 = fStack->Particle(daug1->GetLastDaughter());
941 if( TMath::Abs(daug0->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
942 fHistograms->FillTable("Table_Electrons",18);
943 }//else if
944 }//gamma o Jpsi
945 }//GetNDaughters
946 }
d7d7e825 947 }
a0b94e5c 948
d7d7e825 949 /////////////////////End Chic Analysis////////////////////////////
a0b94e5c 950
951
dc2883e4 952 // if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
a0b94e5c 953
d7d7e825 954 if(particle->R()>fV0Reader->GetMaxRCut()) continue; // cuts on distance from collision point
955
956 Double_t tmpPhi=particle->Phi();
957
958 if(particle->Phi()> TMath::Pi()){
959 tmpPhi = particle->Phi()-(2*TMath::Pi());
960 }
961
962 Double_t rapidity;
963 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
964 rapidity=0;
965 }
966 else{
967 rapidity = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())));
968 }
969
111d75df 970
971
972 if(iTracks<=fStack->GetNprimary() ){
973 if ( particle->GetPdgCode()== -211 || particle->GetPdgCode()== 211 ||
974 particle->GetPdgCode()== 2212 || particle->GetPdgCode()==-2212 ||
975 particle->GetPdgCode()== 321 || particle->GetPdgCode()==-321 ){
dc2883e4 976 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
3c45d101 977 fHistograms->FillHistogram("MC_PhysicalPrimaryCharged_Pt", particle->Pt());
111d75df 978 }
bd6d9fa3 979 }
111d75df 980
bd6d9fa3 981
982
d7d7e825 983 //process the gammas
984 if (particle->GetPdgCode() == 22){
dc2883e4 985 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
a607c218 986
d7d7e825 987 if(particle->GetMother(0) >-1 && fStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
988 continue; // no photon as mothers!
989 }
a607c218 990
d7d7e825 991 if(particle->GetMother(0) >= fStack->GetNprimary()){
992 continue; // the gamma has a mother, and it is not a primary particle
993 }
f8017b04 994
995 if(particle->GetMother(0) >-1){
996 fHistograms->FillHistogram("MC_DecayAllGamma_Pt", particle->Pt()); // All
997 switch(fStack->Particle(particle->GetMother(0))->GetPdgCode()){
998 case 111: // Pi0
999 fHistograms->FillHistogram("MC_DecayPi0Gamma_Pt", particle->Pt());
1000 break;
1001 case 113: // Rho0
1002 fHistograms->FillHistogram("MC_DecayRho0Gamma_Pt", particle->Pt());
1003 break;
1004 case 221: // Eta
1005 fHistograms->FillHistogram("MC_DecayEtaGamma_Pt", particle->Pt());
1006 break;
1007 case 223: // Omega
1008 fHistograms->FillHistogram("MC_DecayOmegaGamma_Pt", particle->Pt());
1009 break;
1010 case 310: // K_s0
1011 fHistograms->FillHistogram("MC_DecayK0sGamma_Pt", particle->Pt());
1012 break;
1013 case 331: // Eta'
1014 fHistograms->FillHistogram("MC_DecayEtapGamma_Pt", particle->Pt());
1015 break;
86e91247 1016 case 333: // Phi
1017 fHistograms->FillHistogram("MC_DecayPhiGamma_Pt", particle->Pt());
1018 break;
f8017b04 1019 }
1020 }
1021
d7d7e825 1022 fHistograms->FillHistogram("MC_allGamma_Energy", particle->Energy());
1023 fHistograms->FillHistogram("MC_allGamma_Pt", particle->Pt());
1024 fHistograms->FillHistogram("MC_allGamma_Eta", particle->Eta());
1025 fHistograms->FillHistogram("MC_allGamma_Phi", tmpPhi);
1026 fHistograms->FillHistogram("MC_allGamma_Rapid", rapidity);
1027
a0b94e5c 1028 // for CF
1e7846f4 1029 if(fDoCF){
1030 containerInput[0] = particle->Pt();
1031 containerInput[1] = particle->Eta();
1032 if(particle->GetMother(0) >=0){
1033 containerInput[2] = fStack->Particle(particle->GetMother(0))->GetMass();
1034 }
1035 else{
1036 containerInput[2]=-1;
1037 }
a607c218 1038
1e7846f4 1039 fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated); // generated gamma
1040 }
86e91247 1041 if(particle->GetMother(0) < 0 || //Phojet p+p -> Direct Photons have no mother
1042 ((particle->GetMother(0) > -1) &&
14d7f2f9 1043 ((TMath::Abs(fStack->Particle(particle->GetMother(0))->GetPdgCode()) < 10)||
1044 (TMath::Abs(fStack->Particle(particle->GetMother(0))->GetPdgCode()) ==21) )) //Pythia p+p -> Direct Photons have quarksor gluons as mother
86e91247 1045 ){ // direct gamma
d7d7e825 1046 fHistograms->FillHistogram("MC_allDirectGamma_Energy",particle->Energy());
1047 fHistograms->FillHistogram("MC_allDirectGamma_Pt", particle->Pt());
1048 fHistograms->FillHistogram("MC_allDirectGamma_Eta", particle->Eta());
1049 fHistograms->FillHistogram("MC_allDirectGamma_Phi", tmpPhi);
1050 fHistograms->FillHistogram("MC_allDirectGamma_Rapid", rapidity);
1051 }
1052
d7d7e825 1053 // looking for conversion (electron + positron from pairbuilding (= 5) )
1054 TParticle* ePos = NULL;
1055 TParticle* eNeg = NULL;
1056
1057 if(particle->GetNDaughters() >= 2){
1058 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
1059 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
1060 if(tmpDaughter->GetUniqueID() == 5){
1061 if(tmpDaughter->GetPdgCode() == 11){
1062 eNeg = tmpDaughter;
1063 }
1064 else if(tmpDaughter->GetPdgCode() == -11){
1065 ePos = tmpDaughter;
1066 }
1067 }
1068 }
1069 }
1070
1071
1072 if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
1073 continue;
1074 }
1075
1076
1077 Double_t ePosPhi = ePos->Phi();
1078 if(ePos->Phi()> TMath::Pi()) ePosPhi = ePos->Phi()-(2*TMath::Pi());
1079
1080 Double_t eNegPhi = eNeg->Phi();
1081 if(eNeg->Phi()> TMath::Pi()) eNegPhi = eNeg->Phi()-(2*TMath::Pi());
1082
1083
1084 if(ePos->Pt()<fV0Reader->GetPtCut() || eNeg->Pt()<fV0Reader->GetPtCut()){
1085 continue; // no reconstruction below the Pt cut
1086 }
a0b94e5c 1087
d7d7e825 1088 if(TMath::Abs(ePos->Eta())> fV0Reader->GetEtaCut() || TMath::Abs(eNeg->Eta())> fV0Reader->GetEtaCut()){
1089 continue;
1090 }
a0b94e5c 1091
d7d7e825 1092 if(ePos->R()>fV0Reader->GetMaxRCut()){
1093 continue; // cuts on distance from collision point
1094 }
a0b94e5c 1095
1096 if(TMath::Abs(ePos->Vz()) > fV0Reader->GetMaxZCut()){
d7d7e825 1097 continue; // outside material
1098 }
1099
a0b94e5c 1100
d7d7e825 1101 if((TMath::Abs(ePos->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() > ePos->R()){
1102 continue; // line cut to exclude regions where we do not reconstruct
1103 }
1104
a0b94e5c 1105
1106 // for CF
1e7846f4 1107 if(fDoCF){
1108 fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructable); // reconstructable gamma
1109 }
d7d7e825 1110 fHistograms->FillHistogram("MC_ConvGamma_Energy", particle->Energy());
1111 fHistograms->FillHistogram("MC_ConvGamma_Pt", particle->Pt());
1112 fHistograms->FillHistogram("MC_ConvGamma_Eta", particle->Eta());
1113 fHistograms->FillHistogram("MC_ConvGamma_Phi", tmpPhi);
1114 fHistograms->FillHistogram("MC_ConvGamma_Rapid", rapidity);
1115 fHistograms->FillHistogram("MC_ConvGamma_Pt_Eta", particle->Pt(),particle->Eta());
1116
1117 fHistograms->FillHistogram("MC_E_Energy", eNeg->Energy());
1118 fHistograms->FillHistogram("MC_E_Pt", eNeg->Pt());
1119 fHistograms->FillHistogram("MC_E_Eta", eNeg->Eta());
1120 fHistograms->FillHistogram("MC_E_Phi", eNegPhi);
1121
1122 fHistograms->FillHistogram("MC_P_Energy", ePos->Energy());
1123 fHistograms->FillHistogram("MC_P_Pt", ePos->Pt());
1124 fHistograms->FillHistogram("MC_P_Eta", ePos->Eta());
1125 fHistograms->FillHistogram("MC_P_Phi", ePosPhi);
1126
1127
d7d7e825 1128 // begin Mapping
1129 Int_t rBin = fHistograms->GetRBin(ePos->R());
9640a3d1 1130 Int_t zBin = fHistograms->GetZBin(ePos->Vz());
d7d7e825 1131 Int_t phiBin = fHistograms->GetPhiBin(particle->Phi());
9c1cb6f7 1132 Double_t rFMD=30;
64f4118c 1133 Double_t rITSTPCMin=50;
1134 Double_t rITSTPCMax=80;
9640a3d1 1135
1136 TVector3 vtxPos(ePos->Vx(),ePos->Vy(),ePos->Vz());
1137
d7d7e825 1138 TString nameMCMappingPhiR="";
e158cbc3 1139 nameMCMappingPhiR.Form("MC_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 1140 // fHistograms->FillHistogram(nameMCMappingPhiR, ePos->Vz(), particle->Eta());
d7d7e825 1141
1142 TString nameMCMappingPhi="";
e158cbc3 1143 nameMCMappingPhi.Form("MC_Conversion_Mapping_Phi%02d",phiBin);
9640a3d1 1144 // fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());
ebcfaa7e 1145 //fHistograms->FillHistogram(nameMCMappingPhi, ePos->Vz(), particle->Eta());
d7d7e825 1146
1147 TString nameMCMappingR="";
e158cbc3 1148 nameMCMappingR.Form("MC_Conversion_Mapping_R%02d",rBin);
9640a3d1 1149 // fHistograms->FillHistogram(nameMCMappingR, particle->Eta());
ebcfaa7e 1150 //fHistograms->FillHistogram(nameMCMappingR,ePos->Vz(), particle->Eta());
d7d7e825 1151
1152 TString nameMCMappingPhiInR="";
e158cbc3 1153 nameMCMappingPhiInR.Form("MC_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 1154 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
1155 fHistograms->FillHistogram(nameMCMappingPhiInR, vtxPos.Phi());
1156
1157 TString nameMCMappingZInR="";
1158 nameMCMappingZInR.Form("MC_Conversion_Mapping_Z_in_R_%02d",rBin);
1159 fHistograms->FillHistogram(nameMCMappingZInR,ePos->Vz() );
1160
1161
1162 TString nameMCMappingPhiInZ="";
1163 nameMCMappingPhiInZ.Form("MC_Conversion_Mapping_Phi_in_Z_%02d",zBin);
1164 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
1165 fHistograms->FillHistogram(nameMCMappingPhiInZ, vtxPos.Phi());
1166
9c1cb6f7 1167
1168 if(ePos->R()<rFMD){
1169 TString nameMCMappingFMDPhiInZ="";
1170 nameMCMappingFMDPhiInZ.Form("MC_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
1171 fHistograms->FillHistogram(nameMCMappingFMDPhiInZ, vtxPos.Phi());
1172 }
1173
64f4118c 1174 if(ePos->R()>rITSTPCMin && ePos->R()<rITSTPCMax){
1175 TString nameMCMappingITSTPCPhiInZ="";
1176 nameMCMappingITSTPCPhiInZ.Form("MC_Conversion_Mapping_ITSTPC_Phi_in_Z_%02d",zBin);
1177 fHistograms->FillHistogram(nameMCMappingITSTPCPhiInZ, vtxPos.Phi());
1178 }
1179
9640a3d1 1180 TString nameMCMappingRInZ="";
1181 nameMCMappingRInZ.Form("MC_Conversion_Mapping_R_in_Z_%02d",zBin);
1182 fHistograms->FillHistogram(nameMCMappingRInZ,ePos->R() );
1183
1184 if(particle->Pt() > fLowPtMapping && particle->Pt()< fHighPtMapping){
1185 TString nameMCMappingMidPtPhiInR="";
1186 nameMCMappingMidPtPhiInR.Form("MC_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
1187 fHistograms->FillHistogram(nameMCMappingMidPtPhiInR, vtxPos.Phi());
1188
1189 TString nameMCMappingMidPtZInR="";
1190 nameMCMappingMidPtZInR.Form("MC_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
1191 fHistograms->FillHistogram(nameMCMappingMidPtZInR,ePos->Vz() );
1192
1193
1194 TString nameMCMappingMidPtPhiInZ="";
1195 nameMCMappingMidPtPhiInZ.Form("MC_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
1196 fHistograms->FillHistogram(nameMCMappingMidPtPhiInZ, vtxPos.Phi());
9c1cb6f7 1197
1198
1199 if(ePos->R()<rFMD){
1200 TString nameMCMappingMidPtFMDPhiInZ="";
1201 nameMCMappingMidPtFMDPhiInZ.Form("MC_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
1202 fHistograms->FillHistogram(nameMCMappingMidPtFMDPhiInZ, vtxPos.Phi());
1203 }
1204
9640a3d1 1205
1206 TString nameMCMappingMidPtRInZ="";
1207 nameMCMappingMidPtRInZ.Form("MC_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
1208 fHistograms->FillHistogram(nameMCMappingMidPtRInZ,ePos->R() );
1209
1210 }
1211
d7d7e825 1212 //end mapping
1213
1214 fHistograms->FillHistogram("MC_Conversion_R",ePos->R());
1215 fHistograms->FillHistogram("MC_Conversion_ZR",ePos->Vz(),ePos->R());
1216 fHistograms->FillHistogram("MC_Conversion_XY",ePos->Vx(),ePos->Vy());
1217 fHistograms->FillHistogram("MC_Conversion_OpeningAngle",GetMCOpeningAngle(ePos, eNeg));
9640a3d1 1218 fHistograms->FillHistogram("MC_ConvGamma_E_AsymmetryP",particle->P(),eNeg->P()/particle->P());
1219 fHistograms->FillHistogram("MC_ConvGamma_P_AsymmetryP",particle->P(),ePos->P()/particle->P());
1220
1221
d7d7e825 1222 if(particle->GetMother(0) < 0){ // no mother = direct gamma, still inside converted
1223 fHistograms->FillHistogram("MC_ConvDirectGamma_Energy",particle->Energy());
1224 fHistograms->FillHistogram("MC_ConvDirectGamma_Pt", particle->Pt());
1225 fHistograms->FillHistogram("MC_ConvDirectGamma_Eta", particle->Eta());
1226 fHistograms->FillHistogram("MC_ConvDirectGamma_Phi", tmpPhi);
1227 fHistograms->FillHistogram("MC_ConvDirectGamma_Rapid", rapidity);
1228
1229 } // end direct gamma
1230 else{ // mother exits
6c84d371 1231 /* if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0
a0b94e5c 1232 fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S
1233 fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2
1234 ){
1235 fMCGammaChic.push_back(particle);
1236 }
6c84d371 1237 */
d7d7e825 1238 } // end if mother exits
1239 } // end if particle is a photon
1240
1241
1242
1243 // process motherparticles (2 gammas as daughters)
1244 // the motherparticle had already to pass the R and the eta cut, but no line cut.
1245 // the line cut is just valid for the conversions!
b8faf1bc 1246
14d7f2f9 1247 // RR primary Pi0 debug ////////////////////////////////////////////////////////////////////////////////////////////
b8faf1bc 1248 if (particle->GetPdgCode()==111){
1249 if( TMath::Abs(rapidity) < fV0Reader->GetRapidityMesonCut() ){
1250 fHistograms->FillHistogram("MC_Pi0_Pt_vs_Rapid_allDaughters", particle->Pt(),rapidity);
1251 }
1252 }
14d7f2f9 1253 // end RR primary Pi0 debug ////////////////////////////////////////////////////////////////////////////////////////
d7d7e825 1254
1255 if(particle->GetNDaughters() == 2){
1256
1257 TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());
1258 TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());
1259
1260 if(daughter0->GetPdgCode() != 22 || daughter1->GetPdgCode() != 22) continue; //check for gamma gamma daughters
dc2883e4 1261
1262 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ) continue;
1263
d7d7e825 1264 // Check the acceptance for both gammas
01b7fdcc 1265 Bool_t gammaEtaCut = kTRUE;
1266 if(TMath::Abs(daughter0->Eta()) > fV0Reader->GetEtaCut() || TMath::Abs(daughter1->Eta()) > fV0Reader->GetEtaCut() ) gammaEtaCut = kFALSE;
1267 Bool_t gammaRCut = kTRUE;
1268 if(daughter0->R() > fV0Reader->GetMaxRCut() || daughter1->R() > fV0Reader->GetMaxRCut() ) gammaRCut = kFALSE;
d7d7e825 1269
1270
1271
1272 // check for conversions now -> have to pass eta, R and line cut!
1273 Bool_t daughter0Electron = kFALSE;
1274 Bool_t daughter0Positron = kFALSE;
1275 Bool_t daughter1Electron = kFALSE;
1276 Bool_t daughter1Positron = kFALSE;
1277
1278 if(daughter0->GetNDaughters() >= 2){ // first gamma
1279 for(Int_t TrackIndex=daughter0->GetFirstDaughter();TrackIndex<=daughter0->GetLastDaughter();TrackIndex++){
1280 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
1281 if(tmpDaughter->GetUniqueID() == 5){
1282 if(tmpDaughter->GetPdgCode() == 11){
1283 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1284 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1285 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1286 daughter0Electron = kTRUE;
1287 }
1288 }
1289
1290 }
1291 }
1292 else if(tmpDaughter->GetPdgCode() == -11){
1293 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1294 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1295 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1296 daughter0Positron = kTRUE;
1297 }
a68437fb 1298 }
1299 }
d7d7e825 1300 }
1301 }
1302 }
1303 }
1304
1305
1306 if(daughter1->GetNDaughters() >= 2){ // second gamma
1307 for(Int_t TrackIndex=daughter1->GetFirstDaughter();TrackIndex<=daughter1->GetLastDaughter();TrackIndex++){
1308 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
1309 if(tmpDaughter->GetUniqueID() == 5){
1310 if(tmpDaughter->GetPdgCode() == 11){
1311 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1312 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1313 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1314 daughter1Electron = kTRUE;
1315 }
1316 }
1317
1318 }
1319 }
1320 else if(tmpDaughter->GetPdgCode() == -11){
1321 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1322 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1323 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1324 daughter1Positron = kTRUE;
1325 }
1326 }
1327
1328 }
1329
1330 }
1331 }
1332 }
1333 }
1334
a0b94e5c 1335
d7d7e825 1336 if(particle->GetPdgCode()==111){ //Pi0
1337 if( iTracks >= fStack->GetNprimary()){
1338 fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());
1339 fHistograms->FillHistogram("MC_Pi0_Secondaries_Rapid", rapidity);
1340 fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);
1341 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt", particle->Pt());
1342 fHistograms->FillHistogram("MC_Pi0_Secondaries_Energy", particle->Energy());
1343 fHistograms->FillHistogram("MC_Pi0_Secondaries_R", particle->R());
1344 fHistograms->FillHistogram("MC_Pi0_Secondaries_ZR", particle->Vz(),particle->R());
1345 fHistograms->FillHistogram("MC_Pi0_Secondaries_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1346 fHistograms->FillHistogram("MC_Pi0_Secondaries_XY", particle->Vx(),particle->Vy());//only fill from one daughter to avoid multiple filling
1347
01b7fdcc 1348 if(gammaEtaCut && gammaRCut){
d7d7e825 1349 //if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
1350 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1351 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1352 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1353 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1354 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1355 }
1356 }
1357 }
1358 else{
1359 fHistograms->FillHistogram("MC_Pi0_Eta", particle->Eta());
1360 fHistograms->FillHistogram("MC_Pi0_Rapid", rapidity);
1361 fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);
1362 fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());
dc2883e4 1363 fHistograms->FillHistogram("MC_Pi0_Pt_vs_Rapid", particle->Pt(),rapidity);
d7d7e825 1364 fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());
1365 fHistograms->FillHistogram("MC_Pi0_R", particle->R());
1366 fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());
1367 fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1368 fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
63e16c52 1369 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_Fiducial", particle->Pt());
f3026589 1370
1371 if(!fV0Reader->GetIsHeavyIon()) {
1372
1373 switch(GetProcessType(fGCMCEvent)){
1374 case kProcSD:
1375 fHistograms->FillHistogram("MC_SD_EvtQ3_Pi0_Pt", particle->Pt());
1376 break;
1377 case kProcDD:
1378 fHistograms->FillHistogram("MC_DD_EvtQ3_Pi0_Pt", particle->Pt());
1379 break;
1380 case kProcND:
1381 fHistograms->FillHistogram("MC_ND_EvtQ3_Pi0_Pt", particle->Pt());
1382 break;
1383 default:
1384 AliError("Unknown Process");
1385 }
1386 }
01b7fdcc 1387 if(gammaEtaCut && gammaRCut){
d7d7e825 1388 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
1389 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1390 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
63e16c52 1391 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_withinAcceptance_Fiducial", particle->Pt());
1392
d7d7e825 1393 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1394 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1395 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1396 fHistograms->FillHistogram("MC_Pi0_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
10e3319b 1397 fHistograms->FillHistogram("MC_Pi0_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
1398 fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
1399 fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
1400
48682642 1401 Double_t alfa=0.;
77ac6f3e 1402 if((daughter0->Energy()+daughter1->Energy()) > 0.){
48682642 1403 alfa= TMath::Abs((daughter0->Energy()-daughter1->Energy())/(daughter0->Energy()+daughter1->Energy()));
1404 }
1405 fHistograms->FillHistogram("MC_Pi0_alpha",alfa);
1406 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_ConvGamma_withinAcceptance_Fiducial", particle->Pt());
1407
d7d7e825 1408 }
1409 }
1410 }
1411 }
1412
1413 if(particle->GetPdgCode()==221){ //Eta
1414 fHistograms->FillHistogram("MC_Eta_Eta", particle->Eta());
1415 fHistograms->FillHistogram("MC_Eta_Rapid", rapidity);
1416 fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);
1417 fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());
dc2883e4 1418 fHistograms->FillHistogram("MC_Eta_Pt_vs_Rapid", particle->Pt(),rapidity);
d7d7e825 1419 fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());
1420 fHistograms->FillHistogram("MC_Eta_R", particle->R());
1421 fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());
1422 fHistograms->FillHistogram("MC_Eta_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1423 fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
1424
01b7fdcc 1425 if(gammaEtaCut && gammaRCut){
a0b94e5c 1426 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 1427 fHistograms->FillHistogram("MC_Eta_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1428 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1429 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1430 fHistograms->FillHistogram("MC_Eta_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1431 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1432 fHistograms->FillHistogram("MC_Eta_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
10e3319b 1433 fHistograms->FillHistogram("MC_Eta_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
1434 fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
1435 fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
1436
dcdc851f 1437 Double_t alfa=0.;
1438 if((daughter0->Energy()+daughter1->Energy()) > 0.){
1439 alfa= TMath::Abs((daughter0->Energy()-daughter1->Energy())/(daughter0->Energy()+daughter1->Energy()));
1440 }
1441 fHistograms->FillHistogram("MC_Eta_alpha",alfa);
1442
d7d7e825 1443 }
1444
1445 }
1446
1447 }
1448
1449 // all motherparticles with 2 gammas as daughters
1450 fHistograms->FillHistogram("MC_Mother_R", particle->R());
1451 fHistograms->FillHistogram("MC_Mother_ZR", particle->Vz(),particle->R());
1452 fHistograms->FillHistogram("MC_Mother_XY", particle->Vx(),particle->Vy());
1453 fHistograms->FillHistogram("MC_Mother_Mass", particle->GetCalcMass());
1454 fHistograms->FillHistogram("MC_Mother_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1455 fHistograms->FillHistogram("MC_Mother_Energy", particle->Energy());
1456 fHistograms->FillHistogram("MC_Mother_Pt", particle->Pt());
1457 fHistograms->FillHistogram("MC_Mother_Eta", particle->Eta());
1458 fHistograms->FillHistogram("MC_Mother_Rapid", rapidity);
1459 fHistograms->FillHistogram("MC_Mother_Phi",tmpPhi);
1460 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt",particle->GetMass(),particle->Pt());
a0b94e5c 1461
01b7fdcc 1462 if(gammaEtaCut && gammaRCut){
a0b94e5c 1463 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 1464 fHistograms->FillHistogram("MC_Mother_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1465 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1466 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_withinAcceptance",particle->GetMass(),particle->Pt());
1467 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1468 fHistograms->FillHistogram("MC_Mother_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1469 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1470 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_ConvGamma_withinAcceptance",particle->GetMass(),particle->Pt());
a0b94e5c 1471
d7d7e825 1472 }
1473
1474
1475 } // end passed R and eta cut
a0b94e5c 1476
d7d7e825 1477 } // end if(particle->GetNDaughters() == 2)
1478
1479 }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++)
a0b94e5c 1480
d7d7e825 1481} // end ProcessMCData
1482
1483
1484
1485void AliAnalysisTaskGammaConversion::FillNtuple(){
1486 //Fills the ntuple with the different values
a0b94e5c 1487
d7d7e825 1488 if(fGammaNtuple == NULL){
1489 return;
1490 }
1491 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1492 for(Int_t i=0;i<numberOfV0s;i++){
1493 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};
1494 AliESDv0 * cV0 = fV0Reader->GetV0(i);
1495 Double_t negPID=0;
1496 Double_t posPID=0;
1497 fV0Reader->GetPIDProbability(negPID,posPID);
1498 values[0]=cV0->GetOnFlyStatus();
1499 values[1]=fV0Reader->CheckForPrimaryVertex();
1500 values[2]=negPID;
1501 values[3]=posPID;
1502 values[4]=fV0Reader->GetX();
1503 values[5]=fV0Reader->GetY();
1504 values[6]=fV0Reader->GetZ();
1505 values[7]=fV0Reader->GetXYRadius();
1506 values[8]=fV0Reader->GetMotherCandidateNDF();
1507 values[9]=fV0Reader->GetMotherCandidateChi2();
1508 values[10]=fV0Reader->GetMotherCandidateEnergy();
1509 values[11]=fV0Reader->GetMotherCandidateEta();
1510 values[12]=fV0Reader->GetMotherCandidatePt();
1511 values[13]=fV0Reader->GetMotherCandidateMass();
1512 values[14]=fV0Reader->GetMotherCandidateWidth();
1513 // values[15]=fV0Reader->GetMotherMCParticle()->Pt(); MOVED TO THE END, HAS TO BE CALLED AFTER HasSameMother NB: still has the same entry in the array
1514 values[16]=fV0Reader->GetOpeningAngle();
1515 values[17]=fV0Reader->GetNegativeTrackEnergy();
1516 values[18]=fV0Reader->GetNegativeTrackPt();
1517 values[19]=fV0Reader->GetNegativeTrackEta();
1518 values[20]=fV0Reader->GetNegativeTrackPhi();
1519 values[21]=fV0Reader->GetPositiveTrackEnergy();
1520 values[22]=fV0Reader->GetPositiveTrackPt();
1521 values[23]=fV0Reader->GetPositiveTrackEta();
1522 values[24]=fV0Reader->GetPositiveTrackPhi();
1523 values[25]=fV0Reader->HasSameMCMother();
1524 if(values[25] != 0){
1525 values[26]=fV0Reader->GetMotherMCParticlePDGCode();
1526 values[15]=fV0Reader->GetMotherMCParticle()->Pt();
1527 }
1528 fTotalNumberOfAddedNtupleEntries++;
1529 fGammaNtuple->Fill(values);
1530 }
1531 fV0Reader->ResetV0IndexNumber();
1532
1533}
1534
1535void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){
1536 // Process all the V0's without applying any cuts to it
a0b94e5c 1537
d7d7e825 1538 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1539 for(Int_t i=0;i<numberOfV0s;i++){
1540 /*AliESDv0 * cV0 = */fV0Reader->GetV0(i);
a0b94e5c 1541
d7d7e825 1542 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
cb90a330 1543 continue;
d7d7e825 1544 }
9640a3d1 1545
77880bd8 1546 // if( !fV0Reader->GetV0(i)->GetOnFlyStatus()){
1547 if( !fV0Reader->CheckV0FinderStatus(i)){
cb90a330 1548 continue;
9640a3d1 1549 }
1550
cb90a330 1551
1552 if( !((fV0Reader->GetNegativeESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ||
1553 !((fV0Reader->GetPositiveESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ){
1554 continue;
9640a3d1 1555 }
cb90a330 1556
ebcfaa7e 1557 if( fV0Reader->GetNegativeESDTrack()->GetSign()== fV0Reader->GetPositiveESDTrack()->GetSign()){
1558 continue;
1559 }
1560
1561 if( fV0Reader->GetNegativeESDTrack()->GetKinkIndex(0) > 0 ||
1562 fV0Reader->GetPositiveESDTrack()->GetKinkIndex(0) > 0) {
1563 continue;
1564 }
9c1cb6f7 1565 if(TMath::Abs(fV0Reader->GetMotherCandidateEta())> fV0Reader->GetEtaCut()){
1566 continue;
1567 }
1568 if(TMath::Abs(fV0Reader->GetPositiveTrackEta())> fV0Reader->GetEtaCut()){
1569 continue;
1570 }
1571 if(TMath::Abs(fV0Reader->GetNegativeTrackEta())> fV0Reader->GetEtaCut()){
1572 continue;
1573 }
1574 if((TMath::Abs(fV0Reader->GetZ())*fV0Reader->GetLineCutZRSlope())-fV0Reader->GetLineCutZValue() > fV0Reader->GetXYRadius() ){ // cuts out regions where we do not reconstruct
1575 continue;
1576 }
f3026589 1577
1578 fHistograms->FillHistogram("ESD_NoCutAllV0_Pt", fV0Reader->GetMotherCandidatePt());
e0155998 1579
1580// RRnewTOF start ///////////////////////////////////////////////
1581 UInt_t statusPos = fV0Reader->GetPositiveESDTrack()->GetStatus();
1582 UInt_t statusNeg = fV0Reader->GetNegativeESDTrack()->GetStatus();
1583
1584 Double_t t0pos = fV0Reader->GetESDpid()->GetTOFResponse().GetStartTime(fV0Reader->GetPositiveTrackP());
1585 Double_t t0neg = fV0Reader->GetESDpid()->GetTOFResponse().GetStartTime(fV0Reader->GetNegativeTrackP());
1586
1587 Double_t timesPos[5];
1588 fV0Reader->GetPositiveESDTrack()->GetIntegratedTimes(timesPos);
1589 Double_t timesNeg[5];
1590 fV0Reader->GetNegativeESDTrack()->GetIntegratedTimes(timesNeg);
1591
1592 Double_t TOFsignalPos = fV0Reader->GetPositiveTrackTOFsignal();
1593 Double_t TOFsignalNeg = fV0Reader->GetNegativeTrackTOFsignal();
1594
1595 Double_t dTpos = TOFsignalPos - t0pos - timesPos[0];
1596 Double_t dTneg = TOFsignalNeg - t0neg - timesNeg[0];
1597
1598 if(statusPos & AliESDtrack::kTOFpid) fHistograms->FillHistogram("ESD_NoCutConvGamma_EandP_P_dT", fV0Reader->GetPositiveTrackP(), dTpos);
1599 if(statusNeg & AliESDtrack::kTOFpid) fHistograms->FillHistogram("ESD_NoCutConvGamma_EandP_P_dT", fV0Reader->GetNegativeTrackP(), dTneg);
1600// RRnewTOF end /////////////////////////////////////////////////
f3026589 1601
d7d7e825 1602 if(fDoMCTruth){
a0b94e5c 1603
d7d7e825 1604 if(fV0Reader->HasSameMCMother() == kFALSE){
1605 continue;
1606 }
a0b94e5c 1607
d7d7e825 1608 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1609 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 1610
d7d7e825 1611 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1612 continue;
1613 }
26923b22 1614 if(negativeMC->GetPdgCode() == positiveMC->GetPdgCode()){
d7d7e825 1615 continue;
1616 }
a68437fb 1617
26923b22 1618 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){ // id 5 is conversion
a68437fb 1619 continue;
1620 }
a0b94e5c 1621
d7d7e825 1622 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a0b94e5c 1623
d7d7e825 1624 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1625 fHistograms->FillHistogram("ESD_NoCutConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1626 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1627 fHistograms->FillHistogram("ESD_NoCutConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1628 fHistograms->FillHistogram("ESD_NoCutConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1629 fHistograms->FillHistogram("ESD_NoCutConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1630 fHistograms->FillHistogram("ESD_NoCutConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1631 fHistograms->FillHistogram("ESD_NoCutConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1632 fHistograms->FillHistogram("ESD_NoCutConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1633 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
a0b94e5c 1634
d7d7e825 1635 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1636 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
1637
1638 fHistograms->FillHistogram("ESD_NoCutConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1639 fHistograms->FillHistogram("ESD_NoCutConversion_R", fV0Reader->GetXYRadius());
1640 fHistograms->FillHistogram("ESD_NoCutConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1641 fHistograms->FillHistogram("ESD_NoCutConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1642 fHistograms->FillHistogram("ESD_NoCutConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1643 fHistograms->FillHistogram("ESD_NoCutConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1644 fHistograms->FillHistogram("ESD_NoCutConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1645 fHistograms->FillHistogram("ESD_NoCutConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1646
1647 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1648 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1649 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1650 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1651
d7d7e825 1652 //store MCTruth properties
1653 fHistograms->FillHistogram("ESD_NoCutConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1654 fHistograms->FillHistogram("ESD_NoCutConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1655 fHistograms->FillHistogram("ESD_NoCutConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
d7d7e825 1656 }
1657 }
1658 }
1659 fV0Reader->ResetV0IndexNumber();
1660}
1661
1662void AliAnalysisTaskGammaConversion::ProcessV0s(){
1663 // see header file for documentation
037dc2db 1664
1665
d7d7e825 1666 if(fWriteNtuple == kTRUE){
1667 FillNtuple();
1668 }
1669
1670 Int_t nSurvivingV0s=0;
5ce758b0 1671 fV0Reader->ResetNGoodV0s();
d7d7e825 1672 while(fV0Reader->NextV0()){
1673 nSurvivingV0s++;
1674
e60f3265 1675
1676 TVector3 vtxConv(fV0Reader->GetX(),fV0Reader->GetY(), fV0Reader->GetZ());
1677
d7d7e825 1678 //-------------------------- filling v0 information -------------------------------------
1679 fHistograms->FillHistogram("ESD_Conversion_R", fV0Reader->GetXYRadius());
1680 fHistograms->FillHistogram("ESD_Conversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1681 fHistograms->FillHistogram("ESD_Conversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1682 fHistograms->FillHistogram("ESD_Conversion_OpeningAngle", fV0Reader->GetOpeningAngle());
e60f3265 1683
1684 // Specific histograms for beam pipe studies
1685 if( TMath::Abs(fV0Reader->GetZ()) < fV0Reader->GetLineCutZValue() ){
1686 fHistograms->FillHistogram("ESD_Conversion_XY_BeamPipe", fV0Reader->GetX(),fV0Reader->GetY());
1687 fHistograms->FillHistogram("ESD_Conversion_RPhi_BeamPipe", vtxConv.Phi(),fV0Reader->GetXYRadius());
1688 }
1689
1690
d7d7e825 1691 fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());
1692 fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());
1693 fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());
1694 fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());
ebcfaa7e 1695 fHistograms->FillHistogram("ESD_E_nTPCClusters", fV0Reader->GetNegativeTracknTPCClusters());
1696 fHistograms->FillHistogram("ESD_E_nITSClusters", fV0Reader->GetNegativeTracknITSClusters());
40051a3e 1697 Double_t eClsToF= 0;
1698 if(!fV0Reader->GetUseCorrectedTPCClsInfo()){
1699 if(fV0Reader->GetNegativeTracknTPCFClusters()!=0 ){
1700 eClsToF=(Double_t)fV0Reader->GetNegativeTracknTPCClusters()/(Double_t)fV0Reader->GetNegativeTracknTPCFClusters();
1701 }
1702 }else{
1703 eClsToF= fV0Reader->GetNegativeESDTrack()->GetTPCClusterInfo(2,0,fV0Reader->GetFirstTPCRow(fV0Reader->GetXYRadius()));
1704 }
1705 fHistograms->FillHistogram("ESD_E_nTPCClustersToFP", fV0Reader->GetNegativeTrackP(),eClsToF );
1706 fHistograms->FillHistogram("ESD_E_nTPCClustersToFR", fV0Reader->GetXYRadius(),eClsToF );
1707
1708 if(fV0Reader->GetNegativeTracknTPCClusters()!=0 ){
5ce758b0 1709 fHistograms->FillHistogram("ESD_E_TPCchi2", fV0Reader->GetNegativeTrackTPCchi2()/(Double_t)fV0Reader->GetNegativeTracknTPCClusters());
1710 }
1711
1712
1713
d7d7e825 1714 fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());
1715 fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());
1716 fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());
1717 fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());
ebcfaa7e 1718 fHistograms->FillHistogram("ESD_P_nTPCClusters", fV0Reader->GetPositiveTracknTPCClusters());
1719 fHistograms->FillHistogram("ESD_P_nITSClusters", fV0Reader->GetPositiveTracknITSClusters());
40051a3e 1720 Double_t pClsToF= 0;
1721 if(!fV0Reader->GetUseCorrectedTPCClsInfo()){
1722 if(fV0Reader->GetPositiveTracknTPCFClusters()!=0){
1723 pClsToF = (Double_t)fV0Reader->GetPositiveTracknTPCClusters()/(Double_t)fV0Reader->GetPositiveTracknTPCFClusters();
1724 }
1725 }else{
1726 pClsToF= fV0Reader->GetPositiveESDTrack()->GetTPCClusterInfo(2,0,fV0Reader->GetFirstTPCRow(fV0Reader->GetXYRadius()));
1727 }
1728
1729 fHistograms->FillHistogram("ESD_P_nTPCClustersToFP",fV0Reader->GetPositiveTrackP(), pClsToF);
1730 fHistograms->FillHistogram("ESD_P_nTPCClustersToFR",fV0Reader->GetXYRadius(), pClsToF);
1731
1732 if(fV0Reader->GetPositiveTracknTPCClusters()!=0){
5ce758b0 1733 fHistograms->FillHistogram("ESD_P_TPCchi2", fV0Reader->GetPositiveTrackTPCchi2()/(Double_t)fV0Reader->GetPositiveTracknTPCClusters());
1734 }
1735
1736
40051a3e 1737
d7d7e825 1738 fHistograms->FillHistogram("ESD_ConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1739 fHistograms->FillHistogram("ESD_ConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1740 fHistograms->FillHistogram("ESD_ConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1741 fHistograms->FillHistogram("ESD_ConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1742 fHistograms->FillHistogram("ESD_ConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1743 fHistograms->FillHistogram("ESD_ConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1744 fHistograms->FillHistogram("ESD_ConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1745 fHistograms->FillHistogram("ESD_ConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1746 fHistograms->FillHistogram("ESD_ConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1747 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1748
1749 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1750 fHistograms->FillHistogram("ESD_ConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
9640a3d1 1751
1752 fHistograms->FillHistogram("ESD_ConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1753 fHistograms->FillHistogram("ESD_ConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1754 fHistograms->FillHistogram("ESD_ConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1755 fHistograms->FillHistogram("ESD_ConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1756
1757 fHistograms->FillHistogram("ESD_ConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1758 fHistograms->FillHistogram("ESD_ConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1759 fHistograms->FillHistogram("ESD_ConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1760 fHistograms->FillHistogram("ESD_ConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
e0155998 1761
1762 UInt_t statusPos = fV0Reader->GetPositiveESDTrack()->GetStatus(); //moved up here from below RRnewTOF
1763 UInt_t statusNeg = fV0Reader->GetNegativeESDTrack()->GetStatus();
1764// RRnewTOF start ///////////////////////////////////////////////
1765 Double_t t0pos = fV0Reader->GetESDpid()->GetTOFResponse().GetStartTime(fV0Reader->GetPositiveTrackP());
1766 Double_t t0neg = fV0Reader->GetESDpid()->GetTOFResponse().GetStartTime(fV0Reader->GetNegativeTrackP());
1767
1768 Double_t timesPos[5];
1769 fV0Reader->GetPositiveESDTrack()->GetIntegratedTimes(timesPos);
1770 Double_t timesNeg[5];
1771 fV0Reader->GetNegativeESDTrack()->GetIntegratedTimes(timesNeg);
1772
1773 Double_t TOFsignalPos = fV0Reader->GetPositiveTrackTOFsignal();
1774 Double_t TOFsignalNeg = fV0Reader->GetNegativeTrackTOFsignal();
1775
1776 Double_t dTpos = TOFsignalPos - t0pos - timesPos[0];
1777 Double_t dTneg = TOFsignalNeg - t0neg - timesNeg[0];
1778
1779 if(statusPos & AliESDtrack::kTOFpid) fHistograms->FillHistogram("ESD_ConvGamma_EandP_P_dT", fV0Reader->GetPositiveTrackP(), dTpos);
1780 if(statusNeg & AliESDtrack::kTOFpid) fHistograms->FillHistogram("ESD_ConvGamma_EandP_P_dT", fV0Reader->GetNegativeTrackP(), dTneg);
1781// RRnewTOF end /////////////////////////////////////////////////
1782
9c1cb6f7 1783 Double_t negPID=0;
1784 Double_t posPID=0;
1785 fV0Reader->GetPIDProbability(negPID,posPID);
1786 fHistograms->FillHistogram("ESD_ConvGamma_E_EProbP",fV0Reader->GetNegativeTrackP(),negPID);
1787 fHistograms->FillHistogram("ESD_ConvGamma_P_EProbP",fV0Reader->GetPositiveTrackP(),posPID);
1788
1789 Double_t negPIDmupi=0;
1790 Double_t posPIDmupi=0;
1791 fV0Reader->GetPIDProbabilityMuonPion(negPIDmupi,posPIDmupi);
1792 fHistograms->FillHistogram("ESD_ConvGamma_E_mupiProbP",fV0Reader->GetNegativeTrackP(),negPIDmupi);
1793 fHistograms->FillHistogram("ESD_ConvGamma_P_mupiProbP",fV0Reader->GetPositiveTrackP(),posPIDmupi);
1794
48682642 1795 Double_t armenterosQtAlfa[2];
1796 fV0Reader->GetArmenterosQtAlfa(fV0Reader-> GetNegativeKFParticle(),
1797 fV0Reader-> GetPositiveKFParticle(),
1798 fV0Reader->GetMotherCandidateKFCombination(),
1799 armenterosQtAlfa);
1800
1801 fHistograms->FillHistogram("ESD_ConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
1802
037dc2db 1803
d7d7e825 1804 // begin mapping
1805 Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());
9640a3d1 1806 Int_t zBin = fHistograms->GetZBin(fV0Reader->GetZ());
d7d7e825 1807 Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());
9c1cb6f7 1808 Double_t rFMD=30;
64f4118c 1809 Double_t rITSTPCMin=50;
1810 Double_t rITSTPCMax=80;
e60f3265 1811
9640a3d1 1812
ebcfaa7e 1813 // Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();
d7d7e825 1814
1815 TString nameESDMappingPhiR="";
e158cbc3 1816 nameESDMappingPhiR.Form("ESD_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 1817 //fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1818
1819 TString nameESDMappingPhi="";
e158cbc3 1820 nameESDMappingPhi.Form("ESD_Conversion_Mapping_Phi%02d",phiBin);
ebcfaa7e 1821 //fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1822
1823 TString nameESDMappingR="";
e158cbc3 1824 nameESDMappingR.Form("ESD_Conversion_Mapping_R%02d",rBin);
ebcfaa7e 1825 //fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1826
1827 TString nameESDMappingPhiInR="";
e158cbc3 1828 nameESDMappingPhiInR.Form("ESD_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 1829 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1830 fHistograms->FillHistogram(nameESDMappingPhiInR, vtxConv.Phi());
1831
1832 TString nameESDMappingZInR="";
1833 nameESDMappingZInR.Form("ESD_Conversion_Mapping_Z_in_R_%02d",rBin);
1834 fHistograms->FillHistogram(nameESDMappingZInR, fV0Reader->GetZ());
1835
1836 TString nameESDMappingPhiInZ="";
1837 nameESDMappingPhiInZ.Form("ESD_Conversion_Mapping_Phi_in_Z_%02d",zBin);
1838 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1839 fHistograms->FillHistogram(nameESDMappingPhiInZ, vtxConv.Phi());
1840
9c1cb6f7 1841 if(fV0Reader->GetXYRadius()<rFMD){
1842 TString nameESDMappingFMDPhiInZ="";
1843 nameESDMappingFMDPhiInZ.Form("ESD_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
1844 fHistograms->FillHistogram(nameESDMappingFMDPhiInZ, vtxConv.Phi());
1845 }
1846
64f4118c 1847 if(fV0Reader->GetXYRadius()>rITSTPCMin && fV0Reader->GetXYRadius()<rITSTPCMax){
1848 TString nameESDMappingITSTPCPhiInZ="";
1849 nameESDMappingITSTPCPhiInZ.Form("ESD_Conversion_Mapping_ITSTPC_Phi_in_Z_%02d",zBin);
1850 fHistograms->FillHistogram(nameESDMappingITSTPCPhiInZ, vtxConv.Phi());
1851 }
9c1cb6f7 1852
9640a3d1 1853 TString nameESDMappingRInZ="";
1854 nameESDMappingRInZ.Form("ESD_Conversion_Mapping_R_in_Z_%02d",zBin);
1855 fHistograms->FillHistogram(nameESDMappingRInZ, fV0Reader->GetXYRadius());
1856
1857 if(fV0Reader->GetMotherCandidatePt() > fLowPtMapping && fV0Reader->GetMotherCandidatePt()< fHighPtMapping){
1858 TString nameESDMappingMidPtPhiInR="";
1859 nameESDMappingMidPtPhiInR.Form("ESD_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
1860 fHistograms->FillHistogram(nameESDMappingMidPtPhiInR, vtxConv.Phi());
1861
1862 TString nameESDMappingMidPtZInR="";
1863 nameESDMappingMidPtZInR.Form("ESD_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
1864 fHistograms->FillHistogram(nameESDMappingMidPtZInR, fV0Reader->GetZ());
1865
1866 TString nameESDMappingMidPtPhiInZ="";
1867 nameESDMappingMidPtPhiInZ.Form("ESD_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
1868 fHistograms->FillHistogram(nameESDMappingMidPtPhiInZ, vtxConv.Phi());
9c1cb6f7 1869 if(fV0Reader->GetXYRadius()<rFMD){
1870 TString nameESDMappingMidPtFMDPhiInZ="";
1871 nameESDMappingMidPtFMDPhiInZ.Form("ESD_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
1872 fHistograms->FillHistogram(nameESDMappingMidPtFMDPhiInZ, vtxConv.Phi());
1873 }
1874
9640a3d1 1875
1876 TString nameESDMappingMidPtRInZ="";
1877 nameESDMappingMidPtRInZ.Form("ESD_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
1878 fHistograms->FillHistogram(nameESDMappingMidPtRInZ, fV0Reader->GetXYRadius());
1879
1880 }
1881
1882
d7d7e825 1883 // end mapping
1884
6c84d371 1885 new((*fKFReconstructedGammasTClone)[fKFReconstructedGammasTClone->GetEntriesFast()]) AliKFParticle(*fV0Reader->GetMotherCandidateKFCombination());
037dc2db 1886 fKFReconstructedGammasV0Index.push_back(fV0Reader->GetCurrentV0IndexNumber()-1);
6c84d371 1887 // fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());
d7d7e825 1888 fElectronv1.push_back(fV0Reader->GetCurrentV0()->GetPindex());
1889 fElectronv2.push_back(fV0Reader->GetCurrentV0()->GetNindex());
a0b94e5c 1890
d7d7e825 1891 //----------------------------------- checking for "real" conversions (MC match) --------------------------------------
1892 if(fDoMCTruth){
f3026589 1893 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1894 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1895
1896 if(fV0Reader->HasSameMCMother() == kFALSE){
1897 fHistograms->FillHistogram("ESD_TrueConvCombinatorial_R", fV0Reader->GetXYRadius());
1898 fHistograms->FillHistogram("ESD_TrueConvCombinatorial_Pt", fV0Reader->GetMotherCandidatePt());
1899 fHistograms->FillHistogram("ESD_TrueConvCombinatorialDaughter_Pt", negativeMC->Pt(),positiveMC->Pt());
1900 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
1901 fHistograms->FillHistogram("ESD_TrueConvCombinatorialElec_R", fV0Reader->GetXYRadius());
1902 fHistograms->FillHistogram("ESD_TrueConvCombinatorialElec_Pt", fV0Reader->GetMotherCandidatePt());
1903 }
1904 if(TMath::Abs(negativeMC->GetPdgCode())==211 && TMath::Abs(positiveMC->GetPdgCode())==211){
1905 fHistograms->FillHistogram("ESD_TrueConvCombinatorialPi_R", fV0Reader->GetXYRadius());
1906 fHistograms->FillHistogram("ESD_TrueConvCombinatorialPi_Pt", fV0Reader->GetMotherCandidatePt());
1907 fHistograms->FillHistogram("ESD_TrueConvCombinatorialPiDaughter_Pt", negativeMC->Pt(),positiveMC->Pt());
1908 }
1909 if((TMath::Abs(negativeMC->GetPdgCode())==211 && TMath::Abs(positiveMC->GetPdgCode())==2211) ||
1910 (TMath::Abs(negativeMC->GetPdgCode())==2212 && TMath::Abs(positiveMC->GetPdgCode())==211)){
1911 fHistograms->FillHistogram("ESD_TrueConvCombinatorialPiP_R", fV0Reader->GetXYRadius());
1912 fHistograms->FillHistogram("ESD_TrueConvCombinatorialPiP_Pt", fV0Reader->GetMotherCandidatePt());
1913 fHistograms->FillHistogram("ESD_TrueConvCombinatorialPiPDaughter_Pt", negativeMC->Pt(),positiveMC->Pt());
1914 }
1915 if((TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==211) ||
1916 (TMath::Abs(negativeMC->GetPdgCode())==211 && TMath::Abs(positiveMC->GetPdgCode())==11)){
1917 fHistograms->FillHistogram("ESD_TrueConvCombinatorialElecPi_R", fV0Reader->GetXYRadius());
1918 fHistograms->FillHistogram("ESD_TrueConvCombinatorialElecPi_Pt", fV0Reader->GetMotherCandidatePt());
1919 }
e0155998 1920 if( (statusPos & AliESDtrack::kTOFpid) && ( TMath::Abs(positiveMC->GetPdgCode()) != 11 ) )
1921 fHistograms->FillHistogram("ESD_TrueConvCombinatorial_DaughtersNotElec_P_dT", fV0Reader->GetPositiveTrackP(), dTpos);//RRnewTOF
1922 if( (statusNeg & AliESDtrack::kTOFpid) && ( TMath::Abs(negativeMC->GetPdgCode()) != 11 ) )
1923 fHistograms->FillHistogram("ESD_TrueConvCombinatorial_DaughtersNotElec_P_dT", fV0Reader->GetNegativeTrackP(), dTneg);//RRnewTOF
f3026589 1924 continue;
1925 }
1926 // Moved up to check true electron background
1927 // TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1928 // TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1929
1930 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1931 fHistograms->FillHistogram("ESD_TrueConvHadronicBck_R", fV0Reader->GetXYRadius());
1932 fHistograms->FillHistogram("ESD_TrueConvHadronicBck_Pt", fV0Reader->GetMotherCandidatePt());
1933 fHistograms->FillHistogram("ESD_TrueConvHadronicBckDaughter_Pt", negativeMC->Pt(),positiveMC->Pt());
e0155998 1934 if(statusPos & AliESDtrack::kTOFpid)
1935 fHistograms->FillHistogram("ESD_TrueConvHadronicBck_Daughters_P_dT", fV0Reader->GetPositiveTrackP(), dTpos);//RRnewTOF
1936 if(statusNeg & AliESDtrack::kTOFpid)
1937 fHistograms->FillHistogram("ESD_TrueConvHadronicBck_Daughters_P_dT", fV0Reader->GetNegativeTrackP(), dTneg);//RRnewTOF
f3026589 1938 if((TMath::Abs(negativeMC->GetPdgCode())==211 && TMath::Abs(positiveMC->GetPdgCode())==2211) ||
1939 (TMath::Abs(negativeMC->GetPdgCode())==2212 && TMath::Abs(positiveMC->GetPdgCode())==211)){
1940 fHistograms->FillHistogram("ESD_TrueConvLambda_R", fV0Reader->GetXYRadius());
1941 fHistograms->FillHistogram("ESD_TrueConvLambda_Pt", fV0Reader->GetMotherCandidatePt());
1942 }
1943 if(TMath::Abs(negativeMC->GetPdgCode())==211 && TMath::Abs(positiveMC->GetPdgCode())==211){
1944 fHistograms->FillHistogram("ESD_TrueConvMeson_R", fV0Reader->GetXYRadius());
1945 fHistograms->FillHistogram("ESD_TrueConvMeson_Pt", fV0Reader->GetMotherCandidatePt());
1946 }
1947 continue;
1948 }
1949
a68437fb 1950
f3026589 1951 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
1952 continue;
1953 }
8de8695b 1954
e0155998 1955 //UInt_t statusPos = fV0Reader->GetPositiveESDTrack()->GetStatus(); moved higher
1956 //UInt_t statusNeg = fV0Reader->GetNegativeESDTrack()->GetStatus();
8de8695b 1957 UChar_t itsPixelPos = fV0Reader->GetPositiveESDTrack()->GetITSClusterMap();
1958 UChar_t itsPixelNeg = fV0Reader->GetNegativeESDTrack()->GetITSClusterMap();
1959
1960 // Using the UniqueID Phojet does not get the Dalitz right
1961 // if( (negativeMC->GetUniqueID() == 4 && positiveMC->GetUniqueID() ==4) ||
1962 // (negativeMC->GetUniqueID() == 0 && positiveMC->GetUniqueID() ==0) ){// fill r distribution for Dalitz decays
1963 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 111){ //pi0
1964 fHistograms->FillHistogram("ESD_TrueDalitzContamination_R", fV0Reader->GetXYRadius());
1965 fHistograms->FillHistogram("ESD_TrueConvDalitzPi0_R", fV0Reader->GetXYRadius());
1966 //--------Histos for HFE
1967
1968 if(statusPos & AliESDtrack::kTOFpid){
1969 fHistograms->FillHistogram("ESD_TrueConvDalitzPi0_SinglePos_R", fV0Reader->GetXYRadius());
1970 if( TESTBIT(itsPixelPos, 0) ){
1971 fHistograms->FillHistogram("ESD_TrueConvDalitzPi0_SinglePos_kFirst_R", fV0Reader->GetXYRadius());
1972 }
037dc2db 1973 }
8de8695b 1974 if(statusNeg & AliESDtrack::kTOFpid){
1975 fHistograms->FillHistogram("ESD_TrueConvDalitzPi0_SingleNeg_R", fV0Reader->GetXYRadius());
1976 if( TESTBIT(itsPixelNeg, 0) ){
1977 fHistograms->FillHistogram("ESD_TrueConvDalitzPi0_SingleNeg_kFirst_R", fV0Reader->GetXYRadius());
1978 }
64f4118c 1979 }
8de8695b 1980 //--------------------------------------------------------
64f4118c 1981
037dc2db 1982 }
8de8695b 1983 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 221){ //eta
1984 fHistograms->FillHistogram("ESD_TrueConvDalitzEta_R", fV0Reader->GetXYRadius());
1985 }
1986
1987 //}
a68437fb 1988
1989 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){// check if the daughters come from a conversion
1990 continue;
1991 }
a0b94e5c 1992
d7d7e825 1993 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a68437fb 1994
26923b22 1995 if(fDoCF){
1996 Double_t containerInput[3];
1997 containerInput[0] = fV0Reader->GetMotherCandidatePt();
1998 containerInput[1] = fV0Reader->GetMotherCandidateEta();
1999 containerInput[2] = fV0Reader->GetMotherCandidateMass();
2000 fCFManager->GetParticleContainer()->Fill(containerInput,kStepTrueGamma); // for CF
2001 }
2002
e0155998 2003// RRnewTOF start ///////////////////////////////////////////////
2004 if(statusPos & AliESDtrack::kTOFpid) fHistograms->FillHistogram("ESD_TrueConvGamma_EandP_P_dT", fV0Reader->GetPositiveTrackP(), dTpos);
2005 if(statusNeg & AliESDtrack::kTOFpid) fHistograms->FillHistogram("ESD_TrueConvGamma_EandP_P_dT", fV0Reader->GetNegativeTrackP(), dTneg);
2006// RRnewTOF end /////////////////////////////////////////////////
2007
d7d7e825 2008 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
2009 fHistograms->FillHistogram("ESD_TrueConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
2010 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
2011 fHistograms->FillHistogram("ESD_TrueConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
2012 fHistograms->FillHistogram("ESD_TrueConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
2013 fHistograms->FillHistogram("ESD_TrueConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
2014 fHistograms->FillHistogram("ESD_TrueConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
2015 fHistograms->FillHistogram("ESD_TrueConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
2016 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
2017 fHistograms->FillHistogram("ESD_TrueConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
2018 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters());
2019 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters());
2020 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters(),fV0Reader->GetMotherCandidateMass());
2021 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters(),fV0Reader->GetMotherCandidateMass());
a0b94e5c 2022
d7d7e825 2023 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
2024 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
a0b94e5c 2025
40051a3e 2026 fHistograms->FillHistogram("ESD_TrueConversion_E_nTPCClustersToFR", fV0Reader->GetXYRadius(),eClsToF );
2027 fHistograms->FillHistogram("ESD_TrueConversion_P_nTPCClustersToFR",fV0Reader->GetXYRadius(), pClsToF);
86e91247 2028
d7d7e825 2029 fHistograms->FillHistogram("ESD_TrueConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
2030 fHistograms->FillHistogram("ESD_TrueConversion_R", fV0Reader->GetXYRadius());
2031 fHistograms->FillHistogram("ESD_TrueConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
2032 fHistograms->FillHistogram("ESD_TrueConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 2033
8de8695b 2034 //----Histos for HFE--------------------------------------
2035
2036 if(statusPos & AliESDtrack::kTOFpid){
2037 fHistograms->FillHistogram("ESD_TrueConversion_SinglePos_R", positiveMC->R(),fV0Reader->GetPositiveMCParticle()->Pt());
2038 if( TESTBIT(itsPixelPos, 0) ){
2039 fHistograms->FillHistogram("ESD_TrueConversion_SinglePos_kFirst_R", positiveMC->R(),fV0Reader->GetPositiveMCParticle()->Pt());
2040 }
2041 }
2042 if(statusNeg & AliESDtrack::kTOFpid){
2043 fHistograms->FillHistogram("ESD_TrueConversion_SingleNeg_R", negativeMC->R(),fV0Reader->GetNegativeMCParticle()->Pt());
2044 if( TESTBIT(itsPixelNeg, 0) ){
2045 fHistograms->FillHistogram("ESD_TrueConversion_SingleNeg_kFirst_R", negativeMC->R(),fV0Reader->GetNegativeMCParticle()->Pt());
2046 }
2047 }
2048 //--------------------------------------------------------
2049
9640a3d1 2050 fHistograms->FillHistogram("ESD_TrueConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
2051 fHistograms->FillHistogram("ESD_TrueConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
2052 fHistograms->FillHistogram("ESD_TrueConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
2053 fHistograms->FillHistogram("ESD_TrueConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
96ade8ef 2054 if (fV0Reader->GetMotherCandidateP() != 0) {
2055 fHistograms->FillHistogram("ESD_TrueConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
2056 fHistograms->FillHistogram("ESD_TrueConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
2057 } else { cout << "Error::fV0Reader->GetNegativeTrackP() == 0 !!!" << endl; }
2058
9640a3d1 2059 fHistograms->FillHistogram("ESD_TrueConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
2060 fHistograms->FillHistogram("ESD_TrueConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
70ef88b5 2061 fHistograms->FillHistogram("ESD_TrueConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
2062
9640a3d1 2063
a0b94e5c 2064
d7d7e825 2065 //store MCTruth properties
2066 fHistograms->FillHistogram("ESD_TrueConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
2067 fHistograms->FillHistogram("ESD_TrueConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
2068 fHistograms->FillHistogram("ESD_TrueConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
48682642 2069
d7d7e825 2070 //resolution
2071 Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();
2072 Double_t esdpt = fV0Reader->GetMotherCandidatePt();
9c1cb6f7 2073 Double_t resdPt = 0.;
4a6157dc 2074 if(mcpt > 0){
9c1cb6f7 2075 resdPt = ((esdpt - mcpt)/mcpt)*100.;
d7d7e825 2076 }
4a6157dc 2077 else if(mcpt < 0){
2078 cout<<"Pt of MC particle is negative, this will cause wrong calculation of resPt"<<endl;
2079 }
d7d7e825 2080
ca6d4600 2081 fHistograms->FillHistogram("Resolution_Gamma_dPt_Pt", mcpt, resdPt);
d7d7e825 2082 fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);
2083 fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);
ca6d4600 2084 fHistograms->FillHistogram("Resolution_Gamma_dPt_Phi", fV0Reader->GetMotherCandidatePhi(), resdPt);
d7d7e825 2085
9c1cb6f7 2086 Double_t resdZ = 0.;
d7d7e825 2087 if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){
9c1cb6f7 2088 resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100.;
d7d7e825 2089 }
9c1cb6f7 2090 Double_t resdZAbs = 0.;
037dc2db 2091 resdZAbs = (fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz());
2092
2093 fHistograms->FillHistogram("Resolution_dZAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdZAbs);
d7d7e825 2094 fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);
2095 fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());
2096 fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());
ca6d4600 2097
2098 // new for dPt_Pt-histograms for Electron and Positron
2099 Double_t mcEpt = fV0Reader->GetNegativeMCParticle()->Pt();
9c1cb6f7 2100 Double_t resEdPt = 0.;
ca6d4600 2101 if (mcEpt > 0){
9c1cb6f7 2102 resEdPt = ((fV0Reader->GetNegativeTrackPt()-mcEpt)/mcEpt)*100.;
ca6d4600 2103 }
2104 UInt_t statusN = fV0Reader->GetNegativeESDTrack()->GetStatus();
2105 // AliESDtrack * negTrk = fV0Reader->GetNegativeESDTrack();
2106 UInt_t kTRDoutN = (statusN & AliESDtrack::kTRDout);
2107
77ac6f3e 2108 Int_t nITSclsE= fV0Reader->GetNegativeTracknITSClusters();
ca6d4600 2109 // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
77ac6f3e 2110 switch(nITSclsE){
ca6d4600 2111 case 0: // 0 ITS clusters
2112 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS0", mcEpt, resEdPt);
2113 break;
2114 case 1: // 1 ITS cluster
2115 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS1", mcEpt, resEdPt);
2116 break;
2117 case 2: // 2 ITS clusters
2118 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS2", mcEpt, resEdPt);
2119 break;
2120 case 3: // 3 ITS clusters
2121 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS3", mcEpt, resEdPt);
2122 break;
2123 case 4: // 4 ITS clusters
2124 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS4", mcEpt, resEdPt);
2125 break;
2126 case 5: // 5 ITS clusters
2127 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS5", mcEpt, resEdPt);
2128 break;
2129 case 6: // 6 ITS clusters
2130 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS6", mcEpt, resEdPt);
2131 break;
2132 }
2133 //Filling histograms with respect to Electron resolution
2134 fHistograms->FillHistogram("Resolution_E_dPt_Pt", mcEpt, resEdPt);
2135 fHistograms->FillHistogram("Resolution_E_dPt_Phi", fV0Reader->GetNegativeTrackPhi(), resEdPt);
2136 if(kTRDoutN){
2137 fHistograms->FillHistogram("Resolution_E_nTRDtracklets_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
2138 fHistograms->FillHistogram("Resolution_E_nTRDtracklets_MCPt", mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
2139 fHistograms->FillHistogram("Resolution_E_nTRDclusters_ESDPt",fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDncls());
2140 fHistograms->FillHistogram("Resolution_E_nTRDclusters_MCPt",mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDncls());
2141 fHistograms->FillHistogram("Resolution_E_TRDsignal_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDsignal());
2142 }
2143
2144 Double_t mcPpt = fV0Reader->GetPositiveMCParticle()->Pt();
2145 Double_t resPdPt = 0;
2146 if (mcPpt > 0){
9c1cb6f7 2147 resPdPt = ((fV0Reader->GetPositiveTrackPt()-mcPpt)/mcPpt)*100.;
ca6d4600 2148 }
2149
2150 UInt_t statusP = fV0Reader->GetPositiveESDTrack()->GetStatus();
2151// AliESDtrack * posTr= fV0Reader->GetPositiveESDTrack();
2152 UInt_t kTRDoutP = (statusP & AliESDtrack::kTRDout);
2153
77ac6f3e 2154 Int_t nITSclsP = fV0Reader->GetPositiveTracknITSClusters();
ca6d4600 2155 // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
77ac6f3e 2156 switch(nITSclsP){
ca6d4600 2157 case 0: // 0 ITS clusters
2158 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS0", mcPpt, resPdPt);
2159 break;
2160 case 1: // 1 ITS cluster
2161 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS1", mcPpt, resPdPt);
2162 break;
2163 case 2: // 2 ITS clusters
2164 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS2", mcPpt, resPdPt);
2165 break;
2166 case 3: // 3 ITS clusters
2167 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS3", mcPpt, resPdPt);
2168 break;
2169 case 4: // 4 ITS clusters
2170 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS4", mcPpt, resPdPt);
2171 break;
2172 case 5: // 5 ITS clusters
2173 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS5", mcPpt, resPdPt);
2174 break;
2175 case 6: // 6 ITS clusters
2176 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS6", mcPpt, resPdPt);
2177 break;
2178 }
2179 //Filling histograms with respect to Positron resolution
2180 fHistograms->FillHistogram("Resolution_P_dPt_Pt", mcPpt, resPdPt);
2181 fHistograms->FillHistogram("Resolution_P_dPt_Phi", fV0Reader->GetPositiveTrackPhi(), resPdPt);
2182 if(kTRDoutP){
2183 fHistograms->FillHistogram("Resolution_P_nTRDtracklets_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
2184 fHistograms->FillHistogram("Resolution_P_nTRDtracklets_MCPt", mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
2185 fHistograms->FillHistogram("Resolution_P_nTRDclusters_ESDPt",fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDncls());
2186 fHistograms->FillHistogram("Resolution_P_nTRDclusters_MCPt",mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDncls());
2187 fHistograms->FillHistogram("Resolution_P_TRDsignal_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDsignal());
2188 }
2189
2190
9c1cb6f7 2191 Double_t resdR = 0.;
d7d7e825 2192 if(fV0Reader->GetNegativeMCParticle()->R() != 0){
9c1cb6f7 2193 resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100.;
d7d7e825 2194 }
9c1cb6f7 2195 Double_t resdRAbs = 0.;
037dc2db 2196 resdRAbs = (fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R());
2197
2198 fHistograms->FillHistogram("Resolution_dRAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdRAbs);
d7d7e825 2199 fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);
2200 fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());
2201 fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());
ca6d4600 2202 fHistograms->FillHistogram("Resolution_R_dPt", fV0Reader->GetNegativeMCParticle()->R(), resdPt);
2203
9c1cb6f7 2204 Double_t resdPhiAbs=0.;
2205 resdPhiAbs=0.;
037dc2db 2206 resdPhiAbs= (fV0Reader->GetMotherCandidatePhi()-fV0Reader->GetNegativeMCParticle()->Phi());
2207 fHistograms->FillHistogram("Resolution_dPhiAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdPhiAbs);
2208
d7d7e825 2209 }//if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22)
d7d7e825 2210 }//if(fDoMCTruth)
2211 }//while(fV0Reader->NextV0)
2212 fHistograms->FillHistogram("ESD_NumberOfSurvivingV0s", nSurvivingV0s);
2213 fHistograms->FillHistogram("ESD_NumberOfV0s", fV0Reader->GetNumberOfV0s());
b5832f95 2214 fHistograms->FillHistogram("ESD_NumberOfContributorsVtx", fV0Reader->GetNumberOfContributorsVtx());
cb90a330 2215
2216 fV0Reader->ResetV0IndexNumber();
d7d7e825 2217}
2218
81900183 2219
6441e967 2220
81900183 2221///_____________________________________________________________________________________
6272370b 2222void AliAnalysisTaskGammaConversion::ProcessGammasForOmegaMesonAnalysis(){
2223 // omega meson analysis pi0+gamma decay
2224 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
48682642 2225 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
6272370b 2226 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
48682642 2227
6272370b 2228 AliKFParticle * omegaCandidateGammaDaughter = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2229 if(fGammav1[firstPi0Index]==firstGammaIndex || fGammav2[firstPi0Index]==firstGammaIndex){
2230 continue;
2231 }
2232
48682642 2233 AliKFParticle omegaCandidate(*omegaCandidatePi0Daughter,*omegaCandidateGammaDaughter);
6272370b 2234 Double_t massOmegaCandidate = 0.;
2235 Double_t widthOmegaCandidate = 0.;
2236
48682642 2237 omegaCandidate.GetMass(massOmegaCandidate,widthOmegaCandidate);
6272370b 2238
e9aea39f 2239 if ( massOmegaCandidate > 733 && massOmegaCandidate < 833 ) {
6441e967 2240 //AddOmegaToAOD(&omegaCandidate, massOmegaCandidate, firstPi0Index, firstGammaIndex);
c59360eb 2241 }
2242
48682642 2243 fHistograms->FillHistogram("ESD_Omega_InvMass_vs_Pt",massOmegaCandidate ,omegaCandidate.GetPt());
6272370b 2244 fHistograms->FillHistogram("ESD_Omega_InvMass",massOmegaCandidate);
48682642 2245
2246 //delete omegaCandidate;
2247
2248 }// end of omega reconstruction in pi0+gamma channel
2249
2250 if(fDoJet == kTRUE){
2251 AliKFParticle* negPiKF=NULL;
2252 AliKFParticle* posPiKF=NULL;
2253
2254 // look at the pi+pi+pi0 channel
2255 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
2256 AliESDtrack* posTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
2257 if (posTrack->GetSign()<0) continue;
c0d9051e 2258 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion))>2.) continue;
48682642 2259 if (posPiKF) delete posPiKF; posPiKF=NULL;
2260 posPiKF = new AliKFParticle( *(posTrack) ,211);
2261
2262 for(Int_t jCh=0;jCh<fChargedParticles->GetEntriesFast();jCh++){
2263 AliESDtrack* negTrack = (AliESDtrack*)(fChargedParticles->At(jCh));
2264 if( negTrack->GetSign()>0) continue;
c0d9051e 2265 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion))>2.) continue;
48682642 2266 if (negPiKF) delete negPiKF; negPiKF=NULL;
2267 negPiKF = new AliKFParticle( *(negTrack) ,-211);
2268 AliKFParticle omegaCandidatePipPinPi0(*omegaCandidatePi0Daughter,*posPiKF,*negPiKF);
2269 Double_t massOmegaCandidatePipPinPi0 = 0.;
2270 Double_t widthOmegaCandidatePipPinPi0 = 0.;
2271
2272 omegaCandidatePipPinPi0.GetMass(massOmegaCandidatePipPinPi0,widthOmegaCandidatePipPinPi0);
ef2e2748 2273
2274 if ( massOmegaCandidatePipPinPi0 > 733 && massOmegaCandidatePipPinPi0 < 833 ) {
6441e967 2275 // AddOmegaToAOD(&omegaCandidatePipPinPi0, massOmegaCandidatePipPinPi0, -1, -1);
ef2e2748 2276 }
48682642 2277
2278 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass_vs_Pt",massOmegaCandidatePipPinPi0 ,omegaCandidatePipPinPi0.GetPt());
2279 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass",massOmegaCandidatePipPinPi0);
2280
2281 // delete omegaCandidatePipPinPi0;
2282 }
2283 }
48682642 2284
45936300 2285 if (posPiKF) delete posPiKF; posPiKF=NULL; if (negPiKF) delete negPiKF; negPiKF=NULL;
48682642 2286
45936300 2287 } // checking ig gammajet because in that case the chargedparticle list is created
6272370b 2288
6272370b 2289 }
48682642 2290
2291 if(fCalculateBackground){
5ce758b0 2292
2293 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2294
2295 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
2296 Int_t mbin = 0;
2297 if(fUseTrackMultiplicityForBG == kTRUE){
2298 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2299 }
2300 else{
2301 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2302 }
2303
2304 AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
2305
48682642 2306 // Background calculation for the omega
2307 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
5ce758b0 2308 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2309
2310 if(fMoveParticleAccordingToVertex == kTRUE){
2311 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
2312 }
48682642 2313 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2314 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
5ce758b0 2315
2316 if(fMoveParticleAccordingToVertex == kTRUE){
2317 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2318 }
2319
48682642 2320 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
2321 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
2322 AliKFParticle * omegaBckCandidate = new AliKFParticle(*omegaCandidatePi0Daughter,previousGoodV0);
2323 Double_t massOmegaBckCandidate = 0.;
2324 Double_t widthOmegaBckCandidate = 0.;
2325
2326 omegaBckCandidate->GetMass(massOmegaBckCandidate,widthOmegaBckCandidate);
ef2e2748 2327
2328
48682642 2329 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass_vs_Pt",massOmegaBckCandidate ,omegaBckCandidate->GetPt());
2330 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass",massOmegaBckCandidate);
2331
2332 delete omegaBckCandidate;
2333 }
2334 }
2335 }
2336 } // end of checking if background calculation is available
6272370b 2337}
d7d7e825 2338
04bf4381 2339
04bf4381 2340
2341
2342
d7d7e825 2343void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
2344 // see header file for documentation
2345
6c84d371 2346 // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
2347 // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
037dc2db 2348
5ce758b0 2349 fESDEvent = fV0Reader->GetESDEvent();
2350
037dc2db 2351 if(fKFReconstructedGammasTClone->GetEntriesFast()>fV0Reader->GetNumberOfV0s()){
2352 cout<<"Warning, number of entries in the tclone is bigger than number of v0s"<<endl;
2353 }
2354
6c84d371 2355 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2356 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
d7d7e825 2357
6c84d371 2358 // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
2359 // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
a0b94e5c 2360
6c84d371 2361 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2362 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
a0b94e5c 2363
d7d7e825 2364 if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
2365 continue;
2366 }
2367 if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
2368 continue;
2369 }
a0b94e5c 2370
d7d7e825 2371 AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
2372
2373 Double_t massTwoGammaCandidate = 0.;
2374 Double_t widthTwoGammaCandidate = 0.;
2375 Double_t chi2TwoGammaCandidate =10000.;
2376 twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
d707e3cf 2377 // if(twoGammaCandidate->GetNDF()>0){
2378 // chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
2379 chi2TwoGammaCandidate = twoGammaCandidate->GetChi2();
d7d7e825 2380
70ef88b5 2381 fHistograms->FillHistogram("ESD_Mother_Chi2",chi2TwoGammaCandidate);
2e2da371 2382 if((chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
d7d7e825 2383
2384 TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
2385 TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
2386
2387 Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
2388 Double_t rapidity;
96ade8ef 2389 if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() <= 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() <= 0){
2390 cout << "Error: |Pz| > E !!!! " << endl;
d7d7e825 2391 rapidity=0;
2392 }
2393 else{
2394 rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
2395 }
dc2883e4 2396
2397 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut()){
2398 delete twoGammaCandidate;
2399 continue; // rapidity cut
2400 }
2401
2402
48682642 2403 Double_t alfa=0.0;
2404 if( (twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()) != 0){
2405 alfa=TMath::Abs((twoGammaDecayCandidateDaughter0->GetE()-twoGammaDecayCandidateDaughter1->GetE())
2406 /(twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()));
2407 }
d7d7e825 2408
1e7846f4 2409 if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
2410 delete twoGammaCandidate;
2411 continue; // minimum opening angle to avoid using ghosttracks
2412 }
2413
67381a40 2414 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 2415 fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
2416 fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
2417 fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
2418 fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
2419 fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
2420 fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
2421 fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
2422 fHistograms->FillHistogram("ESD_Mother_alfa", alfa);
9da64c46 2423 if(massTwoGammaCandidate>0.1 && massTwoGammaCandidate<0.15){
2424 fHistograms->FillHistogram("ESD_Mother_alfa_Pi0", alfa);
2425 }
dcdc851f 2426 if(massTwoGammaCandidate>0.5 && massTwoGammaCandidate<0.57){
2427 fHistograms->FillHistogram("ESD_Mother_alfa_Eta", alfa);
2428 }
2429
6de3471d 2430 fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
2431 fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
2432 fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
2433 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2434 fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
9c1cb6f7 2435 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2436 }
3c45d101 2437 if(alfa<0.1){
2438 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_E_alpha",massTwoGammaCandidate ,twoGammaCandidate->GetE());
2439 }
6de3471d 2440
9640a3d1 2441
9c1cb6f7 2442 if(fCalculateBackground){
5ce758b0 2443 /* Kenneth, just for testing*/
2444 AliGammaConversionBGHandler * bgHandlerTest = fV0Reader->GetBGHandler();
2445
2446 Int_t zbin= bgHandlerTest->GetZBinIndex(fV0Reader->GetVertexZ());
2447 Int_t mbin=0;
2448 Int_t multKAA=0;
2449 if(fUseTrackMultiplicityForBG == kTRUE){
2450 multKAA=fV0Reader->CountESDTracks();
2451 mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2452 }
2453 else{// means we use #v0s for multiplicity
2454 multKAA=fV0Reader->GetNGoodV0s();
2455 mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2456 }
2457 // cout<<"Filling bin number "<<zbin<<" and "<<mbin<<endl;
2458 // cout<<"Corresponding to z = "<<fV0Reader->GetVertexZ()<<" and m = "<<multKAA<<endl;
6de3471d 2459 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2460 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass",zbin,mbin),massTwoGammaCandidate);
2461 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass_vs_Pt",zbin,mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2462 /* end Kenneth, just for testing*/
2463 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2464 }
5ce758b0 2465 }
2466 /* if(fCalculateBackground){
9c1cb6f7 2467 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2468 Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2469 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
5ce758b0 2470 }*/
ebcfaa7e 2471 // if(fDoNeutralMesonV0MCCheck){
2472 if(fDoMCTruth){
037dc2db 2473 //Kenneth: Checking the eta of the gamma to check the difference between 0.9 and 1.2
2474 Int_t indexKF1 = fKFReconstructedGammasV0Index.at(firstGammaIndex);
2475 if(indexKF1<fV0Reader->GetNumberOfV0s()){
2476 fV0Reader->GetV0(indexKF1);//updates to the correct v0
2477 Double_t eta1 = fV0Reader->GetMotherCandidateEta();
2478 Bool_t isRealPi0=kFALSE;
10e3319b 2479 Bool_t isRealEta=kFALSE;
037dc2db 2480 Int_t gamma1MotherLabel=-1;
2481 if(fV0Reader->HasSameMCMother() == kTRUE){
2482 //cout<<"This v0 is a real v0!!!!"<<endl;
2483 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2484 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
2485 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
2486 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
2487 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2488 gamma1MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
2489 }
2490 }
23340c07 2491 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() ==111){
2492 gamma1MotherLabel=-111;
2493 }
d440c756 2494 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() ==221){
2495 gamma1MotherLabel=-221;
2496 }
037dc2db 2497 }
2498 }
2499 Int_t indexKF2 = fKFReconstructedGammasV0Index.at(secondGammaIndex);
2500 if(indexKF1 == indexKF2){
2501 cout<<"index of the two KF particles are the same.... should not happen"<<endl;
2502 }
2503 if(indexKF2<fV0Reader->GetNumberOfV0s()){
2504 fV0Reader->GetV0(indexKF2);
2505 Double_t eta2 = fV0Reader->GetMotherCandidateEta();
2506 Int_t gamma2MotherLabel=-1;
2507 if(fV0Reader->HasSameMCMother() == kTRUE){
2508 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2509 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
2510 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
2511 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
2512 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2513 gamma2MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
2514 }
2515 }
23340c07 2516 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() ==111){
2517 gamma2MotherLabel=-111;
2518 }
d440c756 2519 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() ==221){
2520 gamma2MotherLabel=-221;
2521 }
2522
23340c07 2523 }
037dc2db 2524 }
2525 if(gamma1MotherLabel>=0 && gamma1MotherLabel==gamma2MotherLabel){
2526 if(fV0Reader->CheckIfPi0IsMother(gamma1MotherLabel)){
2527 isRealPi0=kTRUE;
2528 }
10e3319b 2529 if(fV0Reader->CheckIfEtaIsMother(gamma1MotherLabel)){
2530 isRealEta=kTRUE;
2531 }
2532
037dc2db 2533 }
6de3471d 2534 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2535 if(TMath::Abs(eta1)>0.9 && TMath::Abs(eta2)>0.9){
2536 // fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
2537 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2538 if(isRealPi0 || isRealEta){
2539 fHistograms->FillHistogram("ESD_TruePi0_InvMass_1212",massTwoGammaCandidate);
2540 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_1212",openingAngleTwoGammaCandidate);
2541 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2542 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2543 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2544 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2545 }
6de3471d 2546
2547 if(!isRealPi0 && !isRealEta){
2548 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
86e91247 2549 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2550 }else{
86e91247 2551 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2552 }
d440c756 2553 if(gamma1MotherLabel==-111 || gamma2MotherLabel==-111 || gamma1MotherLabel==-221 || gamma2MotherLabel==-221){
86e91247 2554 fHistograms->FillHistogram("ESD_TruePi0DalitzCont_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2555 }
bd6d9fa3 2556 }
bd6d9fa3 2557
6de3471d 2558 }
2559 else if(TMath::Abs(eta1)>0.9 || TMath::Abs(eta2)>0.9){
ebcfaa7e 2560 // fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
2561 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2562
2563 if(isRealPi0 || isRealEta){
2564 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0912",massTwoGammaCandidate);
2565 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0912",openingAngleTwoGammaCandidate);
2566 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2567 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2568 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2569 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2570 }
6de3471d 2571 if(!isRealPi0 && !isRealEta){
2572 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
86e91247 2573 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2574 }else{
86e91247 2575 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2576 }
d440c756 2577 if(gamma1MotherLabel==-111 || gamma2MotherLabel==-111 || gamma1MotherLabel==-221 || gamma2MotherLabel==-221){
86e91247 2578 fHistograms->FillHistogram("ESD_TruePi0DalitzCont_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2579 }
bd6d9fa3 2580 }
2581 }
6de3471d 2582 else{
ebcfaa7e 2583 // fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
2584 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2585 if(isRealPi0 || isRealEta){
2586 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0909",massTwoGammaCandidate);
2587 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0909",openingAngleTwoGammaCandidate);
2588 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2589 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2590 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2591 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2592 if(gamma1MotherLabel > fV0Reader->GetMCStack()->GetNprimary()){
86e91247 2593 fHistograms->FillHistogram("ESD_TruePi0Sec_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2594 }
d707e3cf 2595 }
6de3471d 2596 if(!isRealPi0 && !isRealEta){
2597 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
86e91247 2598 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2599 }else{
86e91247 2600 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2601 }
d440c756 2602 if(gamma1MotherLabel==-111 || gamma2MotherLabel==-111 || gamma1MotherLabel==-221 || gamma2MotherLabel==-221 ){
86e91247 2603 fHistograms->FillHistogram("ESD_TruePi0DalitzCont_InvMass_vs_Pt",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2604 }
bd6d9fa3 2605 }
037dc2db 2606 }
2607 }
2608 }
2609 }
2610 }
6de3471d 2611 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2612 if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
2613 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2614 fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
2615 }
2616
2617 if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
2618 fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
2619 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2620 }
2621 else if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 || TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
2622 fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
2623 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2624 }
2625 else{
2626 fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
2627 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2628 }
ebcfaa7e 2629
6de3471d 2630 Double_t lowMassPi0=0.1;
2631 Double_t highMassPi0=0.15;
95203575 2632 if ( ( massTwoGammaCandidate > lowMassPi0) && (massTwoGammaCandidate < highMassPi0) ){
6de3471d 2633 new((*fKFReconstructedPi0sTClone)[fKFReconstructedPi0sTClone->GetEntriesFast()]) AliKFParticle(*twoGammaCandidate);
2634 fGammav1.push_back(firstGammaIndex);
2635 fGammav2.push_back(secondGammaIndex);
81900183 2636 }
2637
2638 if( fKFCreateAOD ) {
2639 lowMassPi0=0.1;
2640 highMassPi0=0.16;
2641 Double_t lowMassEta=0.5;
2642 Double_t highMassEta=0.6;
2643
2644 if ( ( massTwoGammaCandidate > lowMassPi0) && (massTwoGammaCandidate < highMassPi0) ){
2645 AddPionToAOD(twoGammaCandidate, massTwoGammaCandidate, firstGammaIndex, secondGammaIndex);
2646 } else if ( ( massTwoGammaCandidate > lowMassEta) && (massTwoGammaCandidate < highMassEta) ){
95203575 2647 AddPionToAOD(twoGammaCandidate, massTwoGammaCandidate, firstGammaIndex, secondGammaIndex);
2648 }
81900183 2649 } // if create aod
2650
ebcfaa7e 2651 }
2e2da371 2652 }
6de3471d 2653 delete twoGammaCandidate;
d7d7e825 2654 }
2655 }
2656}
2657
81900183 2658///__________________________________________________________________________________
2659void AliAnalysisTaskGammaConversion::AddToAODBranch(TClonesArray * branch, AliKFParticle * kfParticle, Double_t mass, Int_t daughter1, Int_t daughter2) {
18922343 2660 //Fill AOD with particles
2661 if(branch) {
2662 new((*branch)[branch->GetEntriesFast()]) AliAODConversionParticle(kfParticle, daughter1, daughter2);
2663 AliAODConversionParticle * aodParticle = dynamic_cast<AliAODConversionParticle*>(branch->Last());
2664 if(aodParticle) {
2665 aodParticle->SetChi2(kfParticle->Chi2());
2666 aodParticle->SetIMass(mass);
2667 }
81900183 2668 }
2669}
2670
2671///__________________________________________________________________________________
2672void AliAnalysisTaskGammaConversion::AddPionToAOD(AliKFParticle * kfParticle, Double_t mass, Int_t daughter1, Int_t daughter2) {
18922343 2673 //Add pions to AOD
81900183 2674 AddToAODBranch(fAODPi0, kfParticle, mass, daughter1, daughter2);
2675 TagDaughter(daughter1);
2676 TagDaughter(daughter2);
2677
2678}
2679
2680
2681///__________________________________________________________________________________
2682void AliAnalysisTaskGammaConversion::TagDaughter(Int_t gammaIndex) {
18922343 2683 //Set conversion tag on pion daughters
81900183 2684 AliAODConversionParticle * daughter = dynamic_cast<AliAODConversionParticle*>(fAODGamma->At(gammaIndex));
2685 if(daughter) {
2686 daughter->SetTag(kTRUE);
22eae5f9 2687 } else {
81900183 2688 AliError("Daughter not in gamma tree!!");
2689 }
2690}
2691
2692///___________________________________________________________________________________
2693void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){
2694 // Fill AOD with reconstructed Gamma
81900183 2695 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
2696 AliKFParticle * gammakf = dynamic_cast<AliKFParticle*>(fKFReconstructedGammasTClone->At(gammaIndex));
18922343 2697 if(gammakf) {
2698 AddToAODBranch(fAODGamma, gammakf, gammakf->GetMass(), fElectronv1[gammaIndex], fElectronv2[gammaIndex]);
2699 }
22eae5f9 2700 }
04bf4381 2701}
2702
48682642 2703/*
77ac6f3e 2704 void AliAnalysisTaskGammaConversion::ProcessConvPHOSGammasForNeutralMesonAnalysis(){
2705
6272370b 2706 // see header file for documentation
2707 // Analyse Pi0 with one photon from Phos and 1 photon from conversions
2708
2709
2710
2711 Double_t vtx[3];
2712 vtx[0] = fV0Reader->GetPrimaryVertex()->GetX();
2713 vtx[1] = fV0Reader->GetPrimaryVertex()->GetY();
2714 vtx[2] = fV0Reader->GetPrimaryVertex()->GetZ();
2715
2716
2717 // Loop over all CaloClusters and consider only the PHOS ones:
2718 AliESDCaloCluster *clu;
2719 TLorentzVector pPHOS;
2720 TLorentzVector gammaPHOS;
2721 TLorentzVector gammaGammaConv;
2722 TLorentzVector pi0GammaConvPHOS;
2723 TLorentzVector gammaGammaConvBck;
2724 TLorentzVector pi0GammaConvPHOSBck;
2725
2726
2727 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
77ac6f3e 2728 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2729 if ( !clu->IsPHOS() || clu->E()<0.1 ) continue;
2730 clu ->GetMomentum(pPHOS ,vtx);
2731 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2732 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2733 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz(),0.);
2734 gammaPHOS.SetXYZM(pPHOS.Px(),pPHOS.Py(),pPHOS.Pz(),0.);
2735 pi0GammaConvPHOS=gammaGammaConv+gammaPHOS;
2736 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS",pi0GammaConvPHOS.M());
2737 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvPHOS",pi0GammaConvPHOS.M(),pi0GammaConvPHOS.Pt());
6272370b 2738
77ac6f3e 2739 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz());
2740 TVector3 v3D1(gammaPHOS.Px(),gammaPHOS.Py(),gammaPHOS.Pz());
2741 Double_t opanConvPHOS= v3D0.Angle(v3D1);
2742 if ( opanConvPHOS < 0.35){
2743 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanLow",pi0GammaConvPHOS.M());
2744 }else{
2745 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanHigh",pi0GammaConvPHOS.M());
2746 }
6272370b 2747
77ac6f3e 2748 }
2749
2750 // Now the LorentVector pPHOS is obtained and can be paired with the converted proton
6272370b 2751 }
2752 //==== End of the PHOS cluster selection ============
2753 TLorentzVector pEMCAL;
2754 TLorentzVector gammaEMCAL;
2755 TLorentzVector pi0GammaConvEMCAL;
2756 TLorentzVector pi0GammaConvEMCALBck;
2757
77ac6f3e 2758 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
2759 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2760 if ( !clu->IsEMCAL() || clu->E()<0.1 ) continue;
2761 if (clu->GetNCells() <= 1) continue;
2762 if ( clu->GetTOF()*1e9 < 550 || clu->GetTOF()*1e9 > 750) continue;
6272370b 2763
77ac6f3e 2764 clu ->GetMomentum(pEMCAL ,vtx);
2765 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2766 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2767 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),
2768 twoGammaDecayCandidateDaughter0->Py(),
2769 twoGammaDecayCandidateDaughter0->Pz(),0.);
2770 gammaEMCAL.SetXYZM(pEMCAL.Px(),pEMCAL.Py(),pEMCAL.Pz(),0.);
2771 pi0GammaConvEMCAL=gammaGammaConv+gammaEMCAL;
2772 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL",pi0GammaConvEMCAL.M());
2773 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL",pi0GammaConvEMCAL.M(),pi0GammaConvEMCAL.Pt());
2774 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),
2775 twoGammaDecayCandidateDaughter0->Py(),
2776 twoGammaDecayCandidateDaughter0->Pz());
2777 TVector3 v3D1(gammaEMCAL.Px(),gammaEMCAL.Py(),gammaEMCAL.Pz());
6272370b 2778
77ac6f3e 2779
2780 Double_t opanConvEMCAL= v3D0.Angle(v3D1);
2781 if ( opanConvEMCAL < 0.35){
2782 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanLow",pi0GammaConvEMCAL.M());
2783 }else{
2784 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanHigh",pi0GammaConvEMCAL.M());
2785 }
2786
2787 }
2788 if(fCalculateBackground){
2789 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2790 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
2791 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2792 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2793 gammaGammaConvBck.SetXYZM(previousGoodV0.Px(),
2794 previousGoodV0.Py(),
2795 previousGoodV0.Pz(),0.);
2796 pi0GammaConvEMCALBck=gammaGammaConvBck+gammaEMCAL;
2797 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M());
2798 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M(),
2799 pi0GammaConvEMCALBck.Pt());
2800 }
2801 }
48682642 2802
77ac6f3e 2803 // Now the LorentVector pEMCAL is obtained and can be paired with the converted proton
2804 } // end of checking if background photons are available
2805 }
6272370b 2806 //==== End of the PHOS cluster selection ============
77ac6f3e 2807
2808 }
48682642 2809*/
6272370b 2810
77ac6f3e 2811void AliAnalysisTaskGammaConversion::MoveParticleAccordingToVertex(AliKFParticle * particle,const AliGammaConversionBGHandler::GammaConversionVertex *vertex){
5ce758b0 2812 //see header file for documentation
2813
77ac6f3e 2814 Double_t dx = vertex->fX - fESDEvent->GetPrimaryVertex()->GetX();
2815 Double_t dy = vertex->fY - fESDEvent->GetPrimaryVertex()->GetY();
2816 Double_t dz = vertex->fZ - fESDEvent->GetPrimaryVertex()->GetZ();
5ce758b0 2817
2818 // cout<<"dx, dy, dz: ["<<dx<<","<<dy<<","<<dz<<"]"<<endl;
77ac6f3e 2819 particle->X() = particle->GetX() - dx;
2820 particle->Y() = particle->GetY() - dy;
2821 particle->Z() = particle->GetZ() - dz;
5ce758b0 2822}
2823
111d75df 2824void AliAnalysisTaskGammaConversion::RotateKFParticle(AliKFParticle * kfParticle,Double_t angle){
5510ebe1 2825 // Before rotate needs to be moved to position 0,0,0, ; move back after rotation
2826 Double_t dx = fESDEvent->GetPrimaryVertex()->GetX()-0.;
2827 Double_t dy = fESDEvent->GetPrimaryVertex()->GetY()-0.;
2828 Double_t dz = fESDEvent->GetPrimaryVertex()->GetZ()-0.;
2829
2830 kfParticle->X() = kfParticle->GetX() - dx;
2831 kfParticle->Y() = kfParticle->GetY() - dy;
2832 kfParticle->Z() = kfParticle->GetZ() - dz;
2833
2834
77ac6f3e 2835 // Rotate the kf particle
111d75df 2836 Double_t c = cos(angle);
2837 Double_t s = sin(angle);
2838
5510ebe1 2839 Double_t mA[8][ 8];
2840 for( Int_t i=0; i<8; i++ ){
2841 for( Int_t j=0; j<8; j++){
77ac6f3e 2842 mA[i][j] = 0;
111d75df 2843 }
2844 }
5510ebe1 2845 for( int i=0; i<8; i++ ){
77ac6f3e 2846 mA[i][i] = 1;
111d75df 2847 }
77ac6f3e 2848 mA[0][0] = c; mA[0][1] = s;
2849 mA[1][0] = -s; mA[1][1] = c;
2850 mA[3][3] = c; mA[3][4] = s;
2851 mA[4][3] = -s; mA[4][4] = c;
111d75df 2852
5510ebe1 2853 Double_t mAC[8][8];
2854 Double_t mAp[8];
111d75df 2855
5510ebe1 2856 for( Int_t i=0; i<8; i++ ){
77ac6f3e 2857 mAp[i] = 0;
5510ebe1 2858 for( Int_t k=0; k<8; k++){
77ac6f3e 2859 mAp[i]+=mA[i][k] * kfParticle->GetParameter(k);
111d75df 2860 }
2861 }
2862
5510ebe1 2863 for( Int_t i=0; i<8; i++){
77ac6f3e 2864 kfParticle->Parameter(i) = mAp[i];
111d75df 2865 }
2866
5510ebe1 2867 for( Int_t i=0; i<8; i++ ){
2868 for( Int_t j=0; j<8; j++ ){
77ac6f3e 2869 mAC[i][j] = 0;
5510ebe1 2870 for( Int_t k=0; k<8; k++ ){
77ac6f3e 2871 mAC[i][j]+= mA[i][k] * kfParticle->GetCovariance(k,j);
111d75df 2872 }
2873 }
2874 }
2875
5510ebe1 2876 for( Int_t i=0; i<8; i++ ){
111d75df 2877 for( Int_t j=0; j<=i; j++ ){
2878 Double_t xx = 0;
5510ebe1 2879 for( Int_t k=0; k<8; k++){
77ac6f3e 2880 xx+= mAC[i][k]*mA[j][k];
111d75df 2881 }
2882 kfParticle->Covariance(i,j) = xx;
2883 }
2884 }
5510ebe1 2885
2886 Double_t dx1 = 0.-fESDEvent->GetPrimaryVertex()->GetX();
2887 Double_t dy1 = 0.-fESDEvent->GetPrimaryVertex()->GetY();
2888 Double_t dz1 = 0.-fESDEvent->GetPrimaryVertex()->GetZ();
2889
2890 kfParticle->X() = kfParticle->GetX() - dx1;
2891 kfParticle->Y() = kfParticle->GetY() - dy1;
2892 kfParticle->Z() = kfParticle->GetZ() - dz1;
2893
111d75df 2894}
2895
2896
d7d7e825 2897void AliAnalysisTaskGammaConversion::CalculateBackground(){
2898 // see header file for documentation
5e55d806 2899
2900
2901 TClonesArray * currentEventV0s = fV0Reader->GetCurrentEventGoodV0s();
2902
5ce758b0 2903 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2904
2905 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
2906 Int_t mbin = 0;
2907 if(fUseTrackMultiplicityForBG == kTRUE){
2908 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2909 }
2910 else{
2911 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2912 }
2913
111d75df 2914 if(fDoRotation == kTRUE){
5510ebe1 2915
111d75df 2916 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2917 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2918 for(Int_t iCurrent2=iCurrent+1;iCurrent2<currentEventV0s->GetEntriesFast();iCurrent2++){
77ac6f3e 2919 for(Int_t nRandom=0;nRandom<fNRandomEventsForBG;nRandom++){
111d75df 2920
2921 AliKFParticle currentEventGoodV02 = *(AliKFParticle *)(currentEventV0s->At(iCurrent2));
3bfbe89a 2922
2923 if(fCheckBGProbability == kTRUE){
2924 Double_t massBGprob =0.;
2925 Double_t widthBGprob = 0.;
2926 AliKFParticle *backgroundCandidateProb = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
2927 backgroundCandidateProb->GetMass(massBGprob,widthBGprob);
2928 if(massBGprob>0.1 && massBGprob<0.14){
2ea0dd4a 2929 if(fRandom.Rndm()>bgHandler->GetBGProb(zbin,mbin)){
3bfbe89a 2930 delete backgroundCandidateProb;
2931 continue;
2932 }
2933 }
2934 delete backgroundCandidateProb;
2935 }
111d75df 2936
77ac6f3e 2937 Double_t nRadiansPM = fNDegreesPMBackground*TMath::Pi()/180;
2ea0dd4a 2938
2939 Double_t rotationValue = fRandom.Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
5510ebe1 2940
111d75df 2941 RotateKFParticle(&currentEventGoodV02,rotationValue);
5ce758b0 2942
111d75df 2943 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
5e55d806 2944
5ce758b0 2945 Double_t massBG =0.;
2946 Double_t widthBG = 0.;
2947 Double_t chi2BG =10000.;
2948 backgroundCandidate->GetMass(massBG,widthBG);
111d75df 2949 // if(backgroundCandidate->GetNDF()>0){
2950 chi2BG = backgroundCandidate->GetChi2();
2951 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
2952
2953 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2954 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2955
2956 Double_t openingAngleBG = currentEventGoodV0.GetAngle(currentEventGoodV02);
2957
2958 Double_t rapidity;
2959 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
2960 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2961
dc2883e4 2962 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2963 delete backgroundCandidate;
2964 continue; // rapidity cut
2965 }
2966
111d75df 2967
2968 Double_t alfa=0.0;
2969 if( (currentEventGoodV0.GetE()+currentEventGoodV02.GetE()) != 0){
2970 alfa=TMath::Abs((currentEventGoodV0.GetE()-currentEventGoodV02.GetE())
2971 /(currentEventGoodV0.GetE()+currentEventGoodV02.GetE()));
2972 }
2973
2974
2975 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2976 delete backgroundCandidate;
2977 continue; // minimum opening angle to avoid using ghosttracks
2978 }
2979
2980 // original
111d75df 2981 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 2982 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2983 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2984 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2985 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2986 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2987 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2988 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2989 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2990 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2991 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2992 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2993 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
111d75df 2994 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
6de3471d 2995
dcdc851f 2996 if(massBG>0.1 && massBG<0.15){
2997 fHistograms->FillHistogram("ESD_Background_alfa_Pi0", alfa);
2998 }
2999 if(massBG>0.5 && massBG<0.57){
3000 fHistograms->FillHistogram("ESD_Background_alfa_Eta", alfa);
3001 }
6de3471d 3002
3003 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
3004 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
3005 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
3006 }
3007
3008 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
3009 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
3010 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
3011 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
3012 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
3013 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
3014 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
3015 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3016 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3017 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
3018 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3019 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
3020
3021 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
3022 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3023 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
3024 }
111d75df 3025 }
3c45d101 3026 if(alfa<0.1){
3027 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
3028 }
3029
2e2da371 3030 }
111d75df 3031 //}
5ce758b0 3032 delete backgroundCandidate;
3033 }
3034 }
3035 }
3036 }
111d75df 3037 else{ // means no rotation
3038 AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
3039
3040 if(fUseTrackMultiplicityForBG){
3041 // cout<<"Using charged track multiplicity for background calculation"<<endl;
3042 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
5ce758b0 3043
111d75df 3044 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);//fV0Reader->GetBGGoodV0s(nEventsInBG);
3045
5ce758b0 3046 if(fMoveParticleAccordingToVertex == kTRUE){
3047 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
3048 }
3049
3050 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
3051 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
3052 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
3053 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
111d75df 3054 AliKFParticle previousGoodV0test = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
5ce758b0 3055
111d75df 3056 //cout<<"Primary Vertex event: ["<<fESDEvent->GetPrimaryVertex()->GetX()<<","<<fESDEvent->GetPrimaryVertex()->GetY()<<","<<fESDEvent->GetPrimaryVertex()->GetZ()<<"]"<<endl;
3057 //cout<<"BG prim Vertex event: ["<<bgEventVertex->fX<<","<<bgEventVertex->fY<<","<<bgEventVertex->fZ<<"]"<<endl;
3058
3059 //cout<<"XYZ of particle before transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
5ce758b0 3060 if(fMoveParticleAccordingToVertex == kTRUE){
3061 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
3062 }
111d75df 3063 //cout<<"XYZ of particle after transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
5ce758b0 3064
3065 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
111d75df 3066
5ce758b0 3067 Double_t massBG =0.;
3068 Double_t widthBG = 0.;
3069 Double_t chi2BG =10000.;
3070 backgroundCandidate->GetMass(massBG,widthBG);
5510ebe1 3071
111d75df 3072 // if(backgroundCandidate->GetNDF()>0){
3073 // chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
3074 chi2BG = backgroundCandidate->GetChi2();
3075 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
3076
3077 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
3078 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
3079
3080 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
3081
3082 Double_t rapidity;
3083
3084 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() <= 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() <= 0){
3085 cout << "Error: |Pz| > E !!!! " << endl;
3086 rapidity=0;
3087 } else {
3088 rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
3089 }
dc2883e4 3090 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
3091 delete backgroundCandidate;
3092 continue; // rapidity cut
3093 }
3094
3095
111d75df 3096 Double_t alfa=0.0;
3097 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
3098 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
3099 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
3100 }
3101
3102
3103 if(openingAngleBG < fMinOpeningAngleGhostCut ){
3104 delete backgroundCandidate;
3105 continue; // minimum opening angle to avoid using ghosttracks
3106 }
3107
3108 // original
111d75df 3109 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 3110 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
3111 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
3112 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
3113 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
3114 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
3115 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
3116 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
3117 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3118 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3119 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
3120 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
3121 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
111d75df 3122 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
6de3471d 3123
dcdc851f 3124 if(massBG>0.1 && massBG<0.15){
3125 fHistograms->FillHistogram("ESD_Background_alfa_Pi0", alfa);
3126 }
3127 if(massBG>0.5 && massBG<0.57){
3128 fHistograms->FillHistogram("ESD_Background_alfa_Eta", alfa);
3129 }
6de3471d 3130
3131 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
3132 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
3133 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
3134 }
3135
3136 // test
3137 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
3138 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
3139 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
3140 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
3141 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
3142 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
3143 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
3144 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3145 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3146 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
3147 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3148 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
3149
3150 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
3151 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3152 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
3153 }
3154 // }
111d75df 3155 }
3c45d101 3156 if(alfa<0.1){
3157 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
3158 }
3159
111d75df 3160 }
3161 delete backgroundCandidate;
3162 }
3163 }
3164 }
3165 }
3166 else{ // means using #V0s for multiplicity
3167
3168 // cout<<"Using the v0 multiplicity to calculate background"<<endl;
3169
3170 fHistograms->FillHistogram("ESD_Background_z_m",zbin,mbin);
3171 fHistograms->FillHistogram("ESD_Mother_multpilicityVSv0s",fV0Reader->CountESDTracks(),fV0Reader->GetNumberOfV0s());
3172
3173 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
3174 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);// fV0Reader->GetBGGoodV0s(nEventsInBG);
3175 if(previousEventV0s){
3176
3177 if(fMoveParticleAccordingToVertex == kTRUE){
3178 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
3179 }
3180
3181 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
3182 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
3183 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
3184 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
3185
3186 if(fMoveParticleAccordingToVertex == kTRUE){
3187 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
5ce758b0 3188 }
111d75df 3189
3190 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
3191 Double_t massBG =0.;
3192 Double_t widthBG = 0.;
3193 Double_t chi2BG =10000.;
3194 backgroundCandidate->GetMass(massBG,widthBG);
5510ebe1 3195
111d75df 3196 /* if(backgroundCandidate->GetNDF()>0){
3197 chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
3198 {//remember to remove
3199 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
3200 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
3201
3202 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
3203 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle_nochi2", openingAngleBG);
3204 }
3205 */
2e2da371 3206 chi2BG = backgroundCandidate->GetChi2();
3207 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
5ce758b0 3208 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
3209 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
3210
3211 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
3212
3213 Double_t rapidity;
3214 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
3215 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
3216
dc2883e4 3217 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
3218 delete backgroundCandidate;
3219 continue; // rapidity cut
3220 }
3221
3222
5ce758b0 3223 Double_t alfa=0.0;
3224 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
3225 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
3226 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
3227 }
3228
3229
3230 if(openingAngleBG < fMinOpeningAngleGhostCut ){
3231 delete backgroundCandidate;
3232 continue; // minimum opening angle to avoid using ghosttracks
3233 }
3234
5ce758b0 3235 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 3236 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
3237 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
3238 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
3239 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
3240 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
3241 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
3242 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
3243 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3244 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3245 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
3246 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
3247 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
3248
3249
5ce758b0 3250 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
3c45d101 3251
dcdc851f 3252 if(massBG>0.1 && massBG<0.15){
3253 fHistograms->FillHistogram("ESD_Background_alfa_Pi0", alfa);
3254 }
3255 if(massBG>0.5 && massBG<0.57){
3256 fHistograms->FillHistogram("ESD_Background_alfa_Eta", alfa);
3257 }
3258
6de3471d 3259 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
3260 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
3261 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
3262 }
3263
3264 if(massBG>0.5 && massBG<0.6){
3265 fHistograms->FillHistogram("ESD_Background_alfa_pt0506",momentumVectorbackgroundCandidate.Pt(),alfa);
3266 }
3267 if(massBG>0.3 && massBG<0.4){
3268 fHistograms->FillHistogram("ESD_Background_alfa_pt0304",momentumVectorbackgroundCandidate.Pt(),alfa);
3269 }
3270
3271 // test
3272 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
3273 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
3274 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
3275 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
3276 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
3277 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
3278 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
3279 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3280 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3281 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
3282 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3283 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
3284
3285 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
3286 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3287 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
3288 }
5ce758b0 3289 }
5ce758b0 3290
6de3471d 3291 if(alfa<0.1){
3292 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
5ce758b0 3293 }
2e2da371 3294 // }
5ce758b0 3295 }
111d75df 3296 delete backgroundCandidate;
3297 }
5ce758b0 3298 }
3299 }
d7d7e825 3300 }
111d75df 3301 } // end else (means use #v0s as multiplicity)
3302 } // end no rotation
d7d7e825 3303}
3304
3305
d7d7e825 3306void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){
3307 //ProcessGammasForGammaJetAnalysis
a0b94e5c 3308
d7d7e825 3309 Double_t distIsoMin;
a0b94e5c 3310
d7d7e825 3311 CreateListOfChargedParticles();
a0b94e5c 3312
3313
6c84d371 3314 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
3315 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
3316 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
01b7fdcc 3317 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 3318
01b7fdcc 3319 if( momentumVectorCurrentGamma.Pt()> fMinPtForGammaJet){
d7d7e825 3320 distIsoMin=GetMinimumDistanceToCharge(gammaIndex);
a0b94e5c 3321
d7d7e825 3322 if (distIsoMin > fMinIsoConeSize && fLeadingChargedIndex>=0){
01b7fdcc 3323 CalculateJetCone(gammaIndex);
d7d7e825 3324 }
3325 }
3326 }
3327}
3328
6272370b 3329//____________________________________________________________________
77ac6f3e 3330Bool_t AliAnalysisTaskGammaConversion::IsGoodImpPar(const AliESDtrack *const track)
6272370b 3331{
3332//
3333// check whether particle has good DCAr(Pt) impact
3334// parameter. Only for TPC+ITS tracks (7*sigma cut)
3335// Origin: Andrea Dainese
3336//
3337
3338Float_t d0z0[2],covd0z0[3];
3339track->GetImpactParameters(d0z0,covd0z0);
3340Float_t sigma= 0.0050+0.0060/TMath::Power(track->Pt(),0.9);
3341Float_t d0max = 7.*sigma;
3342if(TMath::Abs(d0z0[0]) < d0max) return kTRUE;
3343
3344return kFALSE;
3345}
3346
3347
d7d7e825 3348void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){
3349 // CreateListOfChargedParticles
a0b94e5c 3350
d7d7e825 3351 fESDEvent = fV0Reader->GetESDEvent();
037dc2db 3352 Int_t numberOfESDTracks=0;
d7d7e825 3353 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
3354 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 3355
d7d7e825 3356 if(!curTrack){
3357 continue;
3358 }
d707e3cf 3359 // Not needed if Standard function used.
3360// if(!IsGoodImpPar(curTrack)){
3361// continue;
3362// }
a0b94e5c 3363
d7d7e825 3364 if(fEsdTrackCuts->AcceptTrack(curTrack) ){
6c84d371 3365 new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack);
3366 // fChargedParticles.push_back(curTrack);
d7d7e825 3367 fChargedParticlesId.push_back(iTracks);
037dc2db 3368 numberOfESDTracks++;
d7d7e825 3369 }
3370 }
6746e1e1 3371// Moved to UserExec using CountAcceptedTracks function. runjet is not needed by default
3372// fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",numberOfESDTracks);
3373// cout<<"esdtracks::"<< numberOfESDTracks<<endl;
3374// if (fV0Reader->GetNumberOfContributorsVtx()>=1){
3375// fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",numberOfESDTracks);
3376// }
d7d7e825 3377}
9c1cb6f7 3378void AliAnalysisTaskGammaConversion::RecalculateV0ForGamma(){
77ac6f3e 3379 //recalculates v0 for gamma
3380
9c1cb6f7 3381 Double_t massE=0.00051099892;
3382 TLorentzVector curElecPos;
3383 TLorentzVector curElecNeg;
3384 TLorentzVector curGamma;
3385
3386 TLorentzVector curGammaAt;
3387 TLorentzVector curElecPosAt;
3388 TLorentzVector curElecNegAt;
3389 AliKFVertex primVtxGamma(*(fESDEvent->GetPrimaryVertex()));
3390 AliKFVertex primVtxImprovedGamma = primVtxGamma;
3391
3392 const AliESDVertex *vtxT3D=fESDEvent->GetPrimaryVertex();
3393
3394 Double_t xPrimaryVertex=vtxT3D->GetXv();
3395 Double_t yPrimaryVertex=vtxT3D->GetYv();
3396 Double_t zPrimaryVertex=vtxT3D->GetZv();
70ef88b5 3397 // Float_t primvertex[3]={xPrimaryVertex,yPrimaryVertex,zPrimaryVertex};
9c1cb6f7 3398
3399 Float_t nsigmaTPCtrackPos;
3400 Float_t nsigmaTPCtrackNeg;
3401 Float_t nsigmaTPCtrackPosToPion;
3402 Float_t nsigmaTPCtrackNegToPion;
3403 AliKFParticle* negKF=NULL;
3404 AliKFParticle* posKF=NULL;
3405
3406 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
3407 AliESDtrack* posTrack = fESDEvent->GetTrack(iTracks);
3408 if(!posTrack){
3409 continue;
3410 }
3411 if (posKF) delete posKF; posKF=NULL;
3412 if(posTrack->GetSign()<0) continue;
3413 if(!(posTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
3414 if(posTrack->GetKinkIndex(0)>0 ) continue;
3415 if(posTrack->GetNcls(1)<50)continue;
3416 Double_t momPos[3];
3417 // posTrack->GetConstrainedPxPyPz(momPos);
3418 posTrack->GetPxPyPz(momPos);
3419 AliESDtrack *ptrk=fESDEvent->GetTrack(iTracks);
3420 curElecPos.SetXYZM(momPos[0],momPos[1],momPos[2],massE);
3421 if(TMath::Abs(curElecPos.Eta())<0.9) continue;
3422 posKF = new AliKFParticle( *(posTrack),-11);
3423
3424 nsigmaTPCtrackPos = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kElectron);
3425 nsigmaTPCtrackPosToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion);
3426
3427 if ( nsigmaTPCtrackPos>5.|| nsigmaTPCtrackPos<-2.){
3428 continue;
3429 }
3430
3431 if(pow((momPos[0]*momPos[0]+momPos[1]*momPos[1]+momPos[2]*momPos[2]),0.5)>0.5 && nsigmaTPCtrackPosToPion<1){
3432 continue;
3433 }
3434
3435
3436
3437 for(Int_t jTracks = 0; jTracks < fESDEvent->GetNumberOfTracks(); jTracks++){
3438 AliESDtrack* negTrack = fESDEvent->GetTrack(jTracks);
3439 if(!negTrack){
3440 continue;
3441 }
3442 if (negKF) delete negKF; negKF=NULL;
3443 if(negTrack->GetSign()>0) continue;
3444 if(!(negTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
3445 if(negTrack->GetKinkIndex(0)>0 ) continue;
3446 if(negTrack->GetNcls(1)<50)continue;
3447 Double_t momNeg[3];
3448 // negTrack->GetConstrainedPxPyPz(momNeg);
3449 negTrack->GetPxPyPz(momNeg);
3450
3451 nsigmaTPCtrackNeg = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kElectron);
3452 nsigmaTPCtrackNegToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion);
3453 if ( nsigmaTPCtrackNeg>5. || nsigmaTPCtrackNeg<-2.){
3454 continue;
3455 }
3456 if(pow((momNeg[0]*momNeg[0]+momNeg[1]*momNeg[1]+momNeg[2]*momNeg[2]),0.5)>0.5 && nsigmaTPCtrackNegToPion<1){
3457 continue;
3458 }
3459 AliESDtrack *ntrk=fESDEvent->GetTrack(jTracks);
3460 curElecNeg.SetXYZM(momNeg[0],momNeg[1],momNeg[2],massE);
3461 if(TMath::Abs(curElecNeg.Eta())<0.9) continue;
3462 negKF = new AliKFParticle( *(negTrack) ,11);
3463
3464 Double_t b=fESDEvent->GetMagneticField();
3465 Double_t xn, xp, dca=ntrk->GetDCA(ptrk,b,xn,xp);
3466 AliExternalTrackParam nt(*ntrk), pt(*ptrk);
3467 nt.PropagateTo(xn,b); pt.PropagateTo(xp,b);
3468
3469
3470 //--- Like in ITSV0Finder
3471 AliExternalTrackParam ntAt0(*ntrk), ptAt0(*ptrk);
3472 Double_t xxP,yyP,alphaP;
3473 Double_t rP[3];
3474
3475 // if (!ptAt0.GetGlobalXYZat(ptAt0->GetX(),xxP,yyP,zzP)) continue;
3476 if (!ptAt0.GetXYZAt(ptAt0.GetX(),b,rP)) continue;
3477 xxP=rP[0];
3478 yyP=rP[1];
3479 alphaP = TMath::ATan2(yyP,xxP);
3480
3481
3482 ptAt0.Propagate(alphaP,0,b);
3483 Float_t ptfacP = (1.+100.*TMath::Abs(ptAt0.GetC(b)));
3484
3485 // Double_t distP = ptAt0.GetY();
3486 Double_t normP = ptfacP*TMath::Sqrt(ptAt0.GetSigmaY2());
3487 Double_t normdist0P = TMath::Abs(ptAt0.GetY()/normP);
3488 Double_t normdist1P = TMath::Abs((ptAt0.GetZ()-zPrimaryVertex)/(ptfacP*TMath::Sqrt(ptAt0.GetSigmaZ2())));
3489 Double_t normdistP = TMath::Sqrt(normdist0P*normdist0P+normdist1P*normdist1P);
3490
3491
3492 Double_t xxN,yyN,alphaN;
3493 Double_t rN[3];
3494 // if (!ntAt0.GetGlobalXYZat(ntAt0->GetX(),xxN,yyN,zzN)) continue;
3495 if (!ntAt0.GetXYZAt(ntAt0.GetX(),b,rN)) continue;
3496 xxN=rN[0];
3497 yyN=rN[1];
3498
3499 alphaN = TMath::ATan2(yyN,xxN);
3500
3501 ntAt0.Propagate(alphaN,0,b);
3502
3503 Float_t ptfacN = (1.+100.*TMath::Abs(ntAt0.GetC(b)));
3504 // Double_t distN = ntAt0.GetY();
3505 Double_t normN = ptfacN*TMath::Sqrt(ntAt0.GetSigmaY2());
3506 Double_t normdist0N = TMath::Abs(ntAt0.GetY()/normN);
3507 Double_t normdist1N = TMath::Abs((ntAt0.GetZ()-zPrimaryVertex)/(ptfacN*TMath::Sqrt(ntAt0.GetSigmaZ2())));
3508 Double_t normdistN = TMath::Sqrt(normdist0N*normdist0N+normdist1N*normdist1N);
3509
3510 //-----------------------------
3511
3512 Double_t momNegAt[3];
3513 nt.GetPxPyPz(momNegAt);
3514 curElecNegAt.SetXYZM(momNegAt[0],momNegAt[1],momNegAt[2],massE);
3515
3516 Double_t momPosAt[3];
3517 pt.GetPxPyPz(momPosAt);
3518 curElecPosAt.SetXYZM(momPosAt[0],momPosAt[1],momPosAt[2],massE);
3519 if(dca>1){
3520 continue;
3521 }
3522
3523 // Double_t dneg= negTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
3524 // Double_t dpos= posTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
3525 AliESDv0 vertex(nt,jTracks,pt,iTracks);
3526
3527
3528 Float_t cpa=vertex.GetV0CosineOfPointingAngle(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex);
3529
9c1cb6f7 3530
70ef88b5 3531
9c1cb6f7 3532 // cout<< "v0Rr::"<< v0Rr<<endl;
3533 // if (pvertex.GetRr()<0.5){
3534 // continue;
3535 //}
3536 // cout<<"vertex.GetChi2V0()"<<vertex.GetChi2V0()<<endl;
3537 if(cpa<0.9)continue;
3538 // if (vertex.GetChi2V0() > 30) continue;
3539 // cout<<"xp+xn::"<<xp<<" "<<xn<<endl;
3540 if ((xn+xp) < 0.4) continue;
3541 if (TMath::Abs(ntrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05)
3542 if (TMath::Abs(ptrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05) continue;
3543
3544 //cout<<"pass"<<endl;
3545
3546 AliKFParticle v0GammaC;
3547 v0GammaC+=(*negKF);
3548 v0GammaC+=(*posKF);
3549 v0GammaC.SetMassConstraint(0,0.001);
3550 primVtxImprovedGamma+=v0GammaC;
3551 v0GammaC.SetProductionVertex(primVtxImprovedGamma);
3552
3553
3554 curGamma=curElecNeg+curElecPos;
3555 curGammaAt=curElecNegAt+curElecPosAt;
3556
3557 // invariant mass versus pt of K0short
3558
3559 Double_t chi2V0GammaC=100000.;
3560 if( v0GammaC.GetNDF() != 0) {
3561 chi2V0GammaC = v0GammaC.GetChi2()/v0GammaC.GetNDF();
3562 }else{
3563 cout<< "ERROR::v0K0C.GetNDF()" << endl;
3564 }
3565
3566 if(chi2V0GammaC<200 &&chi2V0GammaC>0 ){
3567 if(fHistograms != NULL){
3568 fHistograms->FillHistogram("ESD_RecalculateV0_InvMass",v0GammaC.GetMass());
3569 fHistograms->FillHistogram("ESD_RecalculateV0_Pt",v0GammaC.GetPt());
3570 fHistograms->FillHistogram("ESD_RecalculateV0_E_dEdxP",curElecNegAt.P(),negTrack->GetTPCsignal());
3571 fHistograms->FillHistogram("ESD_RecalculateV0_P_dEdxP",curElecPosAt.P(),posTrack->GetTPCsignal());
3572 fHistograms->FillHistogram("ESD_RecalculateV0_cpa",cpa);
3573 fHistograms->FillHistogram("ESD_RecalculateV0_dca",dca);
9c1cb6f7 3574 fHistograms->FillHistogram("ESD_RecalculateV0_normdistP",normdistP);
3575 fHistograms->FillHistogram("ESD_RecalculateV0_normdistN",normdistN);
3576
3577 new((*fKFRecalculatedGammasTClone)[fKFRecalculatedGammasTClone->GetEntriesFast()]) AliKFParticle(v0GammaC);
3578 fElectronRecalculatedv1.push_back(iTracks);
3579 fElectronRecalculatedv2.push_back(jTracks);
3580 }
3581 }
3582 }
3583 }
3584
3585 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();firstGammaIndex++){
3586 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();secondGammaIndex++){
3587 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(firstGammaIndex);
3588 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(secondGammaIndex);
3589
3590 if(fElectronRecalculatedv1[firstGammaIndex]==fElectronRecalculatedv1[secondGammaIndex]){
3591 continue;
3592 }
3593 if( fElectronRecalculatedv2[firstGammaIndex]==fElectronRecalculatedv2[secondGammaIndex]){
3594 continue;
3595 }
3596
3597 AliKFParticle twoGammaCandidate(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
239a56d1 3598 if(fHistograms != NULL){
3599 fHistograms->FillHistogram("ESD_RecalculateGG_InvMass",twoGammaCandidate.GetMass());
3600 fHistograms->FillHistogram("ESD_RecalculateGG_InvMass_vs_Pt",twoGammaCandidate.GetMass(),twoGammaCandidate.GetPt());
3601 }
9c1cb6f7 3602 }
3603 }
3604}
01b7fdcc 3605void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){
6c84d371 3606 // CaculateJetCone
a0b94e5c 3607
d7d7e825 3608 Double_t cone;
3609 Double_t coneSize=0.3;
3610 Double_t ptJet=0;
a0b94e5c 3611
6c84d371 3612 // AliKFParticle * currentGamma = &fKFReconstructedGammas[gammaIndex];
3613 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
98778c17 3614
01b7fdcc 3615 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 3616
6c84d371 3617 AliESDtrack* leadingCharged = (AliESDtrack*)(fChargedParticles->At(fLeadingChargedIndex));
98778c17 3618
d7d7e825 3619 Double_t momLeadingCharged[3];
3620 leadingCharged->GetConstrainedPxPyPz(momLeadingCharged);
a0b94e5c 3621
01b7fdcc 3622 TVector3 momentumVectorLeadingCharged(momLeadingCharged[0],momLeadingCharged[1],momLeadingCharged[2]);
a0b94e5c 3623
01b7fdcc 3624 Double_t phi1=momentumVectorLeadingCharged.Phi();
3625 Double_t eta1=momentumVectorLeadingCharged.Eta();
3626 Double_t phi3=momentumVectorCurrentGamma.Phi();
a0b94e5c 3627
6c84d371 3628 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
3629 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 3630 Int_t chId = fChargedParticlesId[iCh];
3631 if(fLeadingChargedIndex==chId || fLeadingChargedIndex==chId) continue;
3632 Double_t mom[3];
3633 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 3634 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
3635 Double_t phi2=momentumVectorChargedParticle.Phi();
3636 Double_t eta2=momentumVectorChargedParticle.Eta();
a0b94e5c 3637
3638
d7d7e825 3639 cone=100.;
3640 if( TMath::Abs(phi2 - phi1) <= ( TMath::TwoPi()-coneSize) ){
3641 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-phi1),2) );
3642 }else{
3643 if( (phi2 - phi1)> TMath::TwoPi()-coneSize ){
3644 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-TMath::TwoPi()-phi1),2) );
3645 }
3646 if( (phi2 - phi1)< -(TMath::TwoPi()-coneSize) ){
3647 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2+TMath::TwoPi()-phi1),2) );
3648 }
3649 }
a0b94e5c 3650
01b7fdcc 3651 if(cone <coneSize&& momentumVectorChargedParticle.Pt()>fMinPtJetCone ){
3652 ptJet+= momentumVectorChargedParticle.Pt();
3653 Double_t ffzHdrGam = momentumVectorChargedParticle.Pt()/momentumVectorCurrentGamma.Pt();
3654 Double_t imbalanceHdrGam=-momentumVectorChargedParticle.Dot(momentumVectorCurrentGamma)/momentumVectorCurrentGamma.Mag2();
d7d7e825 3655 fHistograms->FillHistogram("ESD_FFzHdrGam",ffzHdrGam);
3656 fHistograms->FillHistogram("ESD_ImbalanceHdrGam",imbalanceHdrGam);
a0b94e5c 3657
d7d7e825 3658 }
a0b94e5c 3659
d7d7e825 3660 Double_t dphiHdrGam=phi3-phi2;
3661 if ( dphiHdrGam < (-TMath::PiOver2())){
3662 dphiHdrGam+=(TMath::TwoPi());
3663 }
a0b94e5c 3664
d7d7e825 3665 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
3666 dphiHdrGam-=(TMath::TwoPi());
3667 }
a0b94e5c 3668
01b7fdcc 3669 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 3670 fHistograms->FillHistogram("ESD_dphiHdrGamIsolated",dphiHdrGam);
3671 }
3672 }//track loop
a0b94e5c 3673
3674
d7d7e825 3675}
3676
3677Double_t AliAnalysisTaskGammaConversion::GetMinimumDistanceToCharge(Int_t indexHighestPtGamma){
3678 // GetMinimumDistanceToCharge
a0b94e5c 3679
d7d7e825 3680 Double_t fIsoMin=100.;
3681 Double_t ptLeadingCharged=-1.;
98778c17 3682
3683 fLeadingChargedIndex=-1;
a0b94e5c 3684
6c84d371 3685 AliKFParticle * gammaHighestPt = (AliKFParticle*)fKFReconstructedGammasTClone->At(indexHighestPtGamma);
01b7fdcc 3686 TVector3 momentumVectorgammaHighestPt(gammaHighestPt->GetPx(),gammaHighestPt->GetPy(),gammaHighestPt->GetPz());
a0b94e5c 3687
01b7fdcc 3688 Double_t phi1=momentumVectorgammaHighestPt.Phi();
3689 Double_t eta1=momentumVectorgammaHighestPt.Eta();
a0b94e5c 3690
6c84d371 3691 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
3692 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 3693 Int_t chId = fChargedParticlesId[iCh];
3694 if(fElectronv1[indexHighestPtGamma]==chId || fElectronv2[indexHighestPtGamma]==chId) continue;
3695 Double_t mom[3];
3696 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 3697 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
3698 Double_t phi2=momentumVectorChargedParticle.Phi();
3699 Double_t eta2=momentumVectorChargedParticle.Eta();
d7d7e825 3700 Double_t iso=pow( (pow( (eta1-eta2),2)+ pow((phi1-phi2),2)),0.5 );
a0b94e5c 3701
01b7fdcc 3702 if(momentumVectorChargedParticle.Pt()>fMinPtIsoCone ){
d7d7e825 3703 if (iso<fIsoMin){
3704 fIsoMin=iso;
3705 }
3706 }
a0b94e5c 3707
d7d7e825 3708 Double_t dphiHdrGam=phi1-phi2;
3709 if ( dphiHdrGam < (-TMath::PiOver2())){
3710 dphiHdrGam+=(TMath::TwoPi());
3711 }
a0b94e5c 3712
d7d7e825 3713 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
3714 dphiHdrGam-=(TMath::TwoPi());
3715 }
01b7fdcc 3716 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 3717 fHistograms->FillHistogram("ESD_dphiHdrGam",dphiHdrGam);
3718 }
a0b94e5c 3719
d7d7e825 3720 if (dphiHdrGam>0.9*TMath::Pi() && dphiHdrGam<1.1*TMath::Pi()){
a0b94e5c 3721 if (momentumVectorChargedParticle.Pt()> ptLeadingCharged && momentumVectorChargedParticle.Pt()>0.1*momentumVectorgammaHighestPt.Pt()){
01b7fdcc 3722 ptLeadingCharged=momentumVectorChargedParticle.Pt();
d7d7e825 3723 fLeadingChargedIndex=iCh;
3724 }
3725 }
a0b94e5c 3726
d7d7e825 3727 }//track loop
3728 fHistograms->FillHistogram("ESD_MinimumIsoDistance",fIsoMin);
3729 return fIsoMin;
a0b94e5c 3730
d7d7e825 3731}
3732
a0b94e5c 3733Int_t AliAnalysisTaskGammaConversion::GetIndexHighestPtGamma(){
3734 //GetIndexHighestPtGamma
3735
d7d7e825 3736 Int_t indexHighestPtGamma=-1;
01b7fdcc 3737 //Double_t
3738 fGammaPtHighest = -100.;
a0b94e5c 3739
6c84d371 3740 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
3741 AliKFParticle * gammaHighestPtCandidate = (AliKFParticle*)fKFReconstructedGammasTClone->At(firstGammaIndex);
a0b94e5c 3742 TVector3 momentumVectorgammaHighestPtCandidate(gammaHighestPtCandidate->GetPx(),gammaHighestPtCandidate->GetPy(),gammaHighestPtCandidate->GetPz());
3743 if (momentumVectorgammaHighestPtCandidate.Pt() > fGammaPtHighest){
3744 fGammaPtHighest=momentumVectorgammaHighestPtCandidate.Pt();
3745 //gammaHighestPt = gammaHighestPtCandidate;
3746 indexHighestPtGamma=firstGammaIndex;
3747 }
d7d7e825 3748 }
a0b94e5c 3749
d7d7e825 3750 return indexHighestPtGamma;
a0b94e5c 3751
d7d7e825 3752}
3753
3754
3755void AliAnalysisTaskGammaConversion::Terminate(Option_t */*option*/)
3756{
3757 // Terminate analysis
3758 //
3759 AliDebug(1,"Do nothing in Terminate");
3760}
3761
3762void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()
3763{
332f1f44 3764
63477d82 3765 if(fKFCreateAOD) {
3766
3767 //AOD
81900183 3768 if(!fAODGamma) fAODGamma = new TClonesArray("AliAODConversionParticle", 0);
63477d82 3769 else fAODGamma->Delete();
3770 fAODGamma->SetName(Form("%s_gamma", fAODBranchName.Data()));
3771
81900183 3772 if(!fAODPi0) fAODPi0 = new TClonesArray("AliAODConversionParticle", 0);
63477d82 3773 else fAODPi0->Delete();
3774 fAODPi0->SetName(Form("%s_Pi0", fAODBranchName.Data()));
3775
81900183 3776 if(!fAODOmega) fAODOmega = new TClonesArray("AliAODConversionParticle", 0);
63477d82 3777 else fAODOmega->Delete();
3778 fAODOmega->SetName(Form("%s_Omega", fAODBranchName.Data()));
3779
3780 //If delta AOD file name set, add in separate file. Else add in standard aod file.
3781 if(GetDeltaAODFileName().Length() > 0) {
3782 AddAODBranch("TClonesArray", &fAODGamma, GetDeltaAODFileName().Data());
3783 AddAODBranch("TClonesArray", &fAODPi0, GetDeltaAODFileName().Data());
3784 AddAODBranch("TClonesArray", &fAODOmega, GetDeltaAODFileName().Data());
3785 AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(GetDeltaAODFileName().Data());
3786 } else {
3787 AddAODBranch("TClonesArray", &fAODGamma);
3788 AddAODBranch("TClonesArray", &fAODPi0);
3789 AddAODBranch("TClonesArray", &fAODOmega);
3790 }
332f1f44 3791 }
04bf4381 3792
d7d7e825 3793 // Create the output container
3794 if(fOutputContainer != NULL){
3795 delete fOutputContainer;
3796 fOutputContainer = NULL;
3797 }
3798 if(fOutputContainer == NULL){
3799 fOutputContainer = new TList();
b58d3c74 3800 fOutputContainer->SetOwner(kTRUE);
d7d7e825 3801 }
3802
3803 //Adding the histograms to the output container
3804 fHistograms->GetOutputContainer(fOutputContainer);
3805
3806
3807 if(fWriteNtuple){
3808 if(fGammaNtuple == NULL){
3809 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);
3810 }
3811 if(fNeutralMesonNtuple == NULL){
3812 fNeutralMesonNtuple = new TNtuple("NeutralMesonNtuple","NeutralMesonNtuple","test");
3813 }
3814 TList * ntupleTList = new TList();
b58d3c74 3815 ntupleTList->SetOwner(kTRUE);
d7d7e825 3816 ntupleTList->SetName("Ntuple");
3817 ntupleTList->Add((TNtuple*)fGammaNtuple);
3818 fOutputContainer->Add(ntupleTList);
3819 }
3820
3821 fOutputContainer->SetName(GetName());
143874ce 3822
3823 PostData(1, fOutputContainer);
3824 PostData(2, fCFManager->GetParticleContainer()); // for CF
3825
d7d7e825 3826}
3827
77ac6f3e 3828Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(const TParticle* const daughter0, const TParticle* const daughter1) const{
d7d7e825 3829 //helper function
3830 TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());
3831 TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());
3832 return v3D0.Angle(v3D1);
3833}
3834
3835void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){
3836 // see header file for documentation
5e55d806 3837
d7d7e825 3838 vector<Int_t> indexOfGammaParticle;
a0b94e5c 3839
d7d7e825 3840 fStack = fV0Reader->GetMCStack();
a0b94e5c 3841
d7d7e825 3842 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
3843 return; // aborts if the primary vertex does not have contributors.
3844 }
a0b94e5c 3845
d7d7e825 3846 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
3847 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
3848 if(particle->GetPdgCode()==22){ //Gamma
3849 if(particle->GetNDaughters() >= 2){
3850 TParticle* electron=NULL;
3851 TParticle* positron=NULL;
3852 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
3853 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
3854 if(tmpDaughter->GetUniqueID() == 5){
3855 if(tmpDaughter->GetPdgCode() == 11){
3856 electron = tmpDaughter;
3857 }
3858 else if(tmpDaughter->GetPdgCode() == -11){
3859 positron = tmpDaughter;
3860 }
3861 }
3862 }
3863 if(electron!=NULL && positron!=0){
3864 if(electron->R()<160){
3865 indexOfGammaParticle.push_back(iTracks);
3866 }
3867 }
3868 }
3869 }
3870 }
a0b94e5c 3871
d7d7e825 3872 Int_t nFoundGammas=0;
3873 Int_t nNotFoundGammas=0;
a0b94e5c 3874
d7d7e825 3875 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
3876 for(Int_t i=0;i<numberOfV0s;i++){
3877 fV0Reader->GetV0(i);
a0b94e5c 3878
d7d7e825 3879 if(fV0Reader->HasSameMCMother() == kFALSE){
3880 continue;
3881 }
a0b94e5c 3882
d7d7e825 3883 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
3884 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 3885
d7d7e825 3886 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
3887 continue;
3888 }
3889 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
3890 continue;
3891 }
a0b94e5c 3892
d7d7e825 3893 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
3894 //TParticle * v0Gamma = fV0Reader->GetMotherMCParticle();
3895 for(UInt_t mcIndex=0;mcIndex<indexOfGammaParticle.size();mcIndex++){
3896 if(negativeMC->GetFirstMother()==indexOfGammaParticle[mcIndex]){
3897 nFoundGammas++;
3898 }
3899 else{
3900 nNotFoundGammas++;
3901 }
3902 }
3903 }
3904 }
d7d7e825 3905}
3906
3907
3908void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){
3909 // see header file for documantation
a0b94e5c 3910
d7d7e825 3911 fESDEvent = fV0Reader->GetESDEvent();
a0b94e5c 3912
3913
3914 TClonesArray * vESDeNegTemp = new TClonesArray("AliESDtrack",0);
6c84d371 3915 TClonesArray * vESDePosTemp = new TClonesArray("AliESDtrack",0);
3916 TClonesArray * vESDxNegTemp = new TClonesArray("AliESDtrack",0);
3917 TClonesArray * vESDxPosTemp = new TClonesArray("AliESDtrack",0);
3918 TClonesArray * vESDeNegNoJPsi = new TClonesArray("AliESDtrack",0);
3919 TClonesArray * vESDePosNoJPsi = new TClonesArray("AliESDtrack",0);
a0b94e5c 3920
6c84d371 3921 /*
3922 vector <AliESDtrack*> vESDeNegTemp(0);
3923 vector <AliESDtrack*> vESDePosTemp(0);
3924 vector <AliESDtrack*> vESDxNegTemp(0);
3925 vector <AliESDtrack*> vESDxPosTemp(0);
3926 vector <AliESDtrack*> vESDeNegNoJPsi(0);
3927 vector <AliESDtrack*> vESDePosNoJPsi(0);
3928 */
a0b94e5c 3929
3930
d7d7e825 3931 fHistograms->FillTable("Table_Electrons",0);//Count number of Events
a0b94e5c 3932
d7d7e825 3933 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
3934 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 3935
d7d7e825 3936 if(!curTrack){
3937 //print warning here
3938 continue;
3939 }
a0b94e5c 3940
d7d7e825 3941 double p[3];if(!curTrack->GetConstrainedPxPyPz(p))continue;
3942 double r[3];curTrack->GetConstrainedXYZ(r);
a0b94e5c 3943
d7d7e825 3944 TVector3 rXYZ(r);
a0b94e5c 3945
d7d7e825 3946 fHistograms->FillTable("Table_Electrons",4);//Count number of ESD tracks
a0b94e5c 3947
d7d7e825 3948 Bool_t flagKink = kTRUE;
3949 Bool_t flagTPCrefit = kTRUE;
3950 Bool_t flagTRDrefit = kTRUE;
3951 Bool_t flagITSrefit = kTRUE;
3952 Bool_t flagTRDout = kTRUE;
3953 Bool_t flagVertex = kTRUE;
a0b94e5c 3954
3955
d7d7e825 3956 //Cuts ---------------------------------------------------------------
a0b94e5c 3957
d7d7e825 3958 if(curTrack->GetKinkIndex(0) > 0){
3959 fHistograms->FillHistogram("Table_Electrons",5);//Count kink
3960 flagKink = kFALSE;
3961 }
a0b94e5c 3962
d7d7e825 3963 ULong_t trkStatus = curTrack->GetStatus();
a0b94e5c 3964
d7d7e825 3965 ULong_t tpcRefit = (trkStatus & AliESDtrack::kTPCrefit);
a0b94e5c 3966
d7d7e825 3967 if(!tpcRefit){
3968 fHistograms->FillHistogram("Table_Electrons",9);//Count not TPCrefit
3969 flagTPCrefit = kFALSE;
3970 }
a0b94e5c 3971
d7d7e825 3972 ULong_t itsRefit = (trkStatus & AliESDtrack::kITSrefit);
3973 if(!itsRefit){
3974 fHistograms->FillHistogram("Table_Electrons",10);//Count not ITSrefit
3975 flagITSrefit = kFALSE;
3976 }
a0b94e5c 3977
d7d7e825 3978 ULong_t trdRefit = (trkStatus & AliESDtrack::kTRDrefit);
a0b94e5c 3979
d7d7e825 3980 if(!trdRefit){
3981 fHistograms->FillHistogram("Table_Electrons",8); //Count not TRDrefit
3982 flagTRDrefit = kFALSE;
3983 }
a0b94e5c 3984
d7d7e825 3985 ULong_t trdOut = (trkStatus & AliESDtrack::kTRDout);
a0b94e5c 3986
d7d7e825 3987 if(!trdOut) {
3988 fHistograms->FillHistogram("Table_Electrons",7); //Count not TRDout
3989 flagTRDout = kFALSE;
3990 }
a0b94e5c 3991
d7d7e825 3992 double nsigmaToVxt = GetSigmaToVertex(curTrack);
a0b94e5c 3993
d7d7e825 3994 if(nsigmaToVxt > 3){
3995 fHistograms->FillHistogram("Table_Electrons",6); //Count Tracks with number of sigmas > 3
3996 flagVertex = kFALSE;
3997 }
a0b94e5c 3998
d7d7e825 3999 if(! (flagKink && flagTPCrefit && flagITSrefit && flagTRDrefit && flagTRDout && flagVertex ) ) continue;
4000 fHistograms->FillHistogram("Table_Electrons",11);//Count Tracks passed Cuts
a0b94e5c 4001
4002
d7d7e825 4003 Stat_t pid, weight;
4004 GetPID(curTrack, pid, weight);
a0b94e5c 4005
d7d7e825 4006 if(pid!=0){
4007 fHistograms->FillHistogram("Table_Electrons",12); //Count Tracks with pid != 0
4008 }
a0b94e5c 4009
d7d7e825 4010 if(pid == 0){
4011 fHistograms->FillHistogram("Table_Electrons",13); //Count Tracks with pid != 0
4012 }
a0b94e5c 4013
4014
4015
4016
a0b94e5c 4017
4018
d7d7e825 4019 TLorentzVector curElec;
4020 curElec.SetXYZM(p[0],p[1],p[2],fElectronMass);
a0b94e5c 4021
4022
113d8432 4023 if(fDoMCTruth){
4024 Int_t labelMC = TMath::Abs(curTrack->GetLabel());
4025 TParticle* curParticle = fStack->Particle(labelMC);
4026 if(curTrack->GetSign() > 0){
4027 if( pid == 0){
4028 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
4029 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
4030 }
4031 else{
4032 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
4033 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
4034 }
4035 }
4036 }
a0b94e5c 4037
4038
d7d7e825 4039 if(curTrack->GetSign() > 0){
a0b94e5c 4040
6c84d371 4041 // vESDxPosTemp.push_back(curTrack);
4042 new((*vESDxPosTemp)[vESDxPosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 4043
d7d7e825 4044 if( pid == 0){
a0b94e5c 4045
d7d7e825 4046 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
4047 fHistograms->FillHistogram("ESD_ElectronPosPt",curElec.Pt());
113d8432 4048 // fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
d7d7e825 4049 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
113d8432 4050 // fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
6c84d371 4051 // vESDePosTemp.push_back(curTrack);
4052 new((*vESDePosTemp)[vESDePosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 4053
d7d7e825 4054 }
a0b94e5c 4055
d7d7e825 4056 }
4057 else {
5e55d806 4058
6c84d371 4059 new((*vESDxNegTemp)[vESDxNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 4060
d7d7e825 4061 if( pid == 0){
a0b94e5c 4062
d7d7e825 4063 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
4064 fHistograms->FillHistogram("ESD_ElectronNegPt",curElec.Pt());
d7d7e825 4065 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
6c84d371 4066 new((*vESDeNegTemp)[vESDeNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 4067
d7d7e825 4068 }
a0b94e5c 4069
d7d7e825 4070 }
a0b94e5c 4071
d7d7e825 4072 }
a0b94e5c 4073
4074
d7d7e825 4075 Bool_t ePosJPsi = kFALSE;
4076 Bool_t eNegJPsi = kFALSE;
4077 Bool_t ePosPi0 = kFALSE;
4078 Bool_t eNegPi0 = kFALSE;
4079
4080 UInt_t iePosJPsi=0,ieNegJPsi=0,iePosPi0=0,ieNegPi0=0;
a0b94e5c 4081
6c84d371 4082 for(Int_t iNeg=0; iNeg < vESDeNegTemp->GetEntriesFast(); iNeg++){
4083 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetPdgCode() == 11)
4084 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0) > -1){
4085 Int_t labelMother = fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0);
d7d7e825 4086 TParticle* partMother = fStack ->Particle(labelMother);
4087 if (partMother->GetPdgCode() == 111){
4088 ieNegPi0 = iNeg;
4089 eNegPi0 = kTRUE;
4090 }
4091 if(partMother->GetPdgCode() == 443){ //Mother JPsi
4092 fHistograms->FillTable("Table_Electrons",14);
4093 ieNegJPsi = iNeg;
4094 eNegJPsi = kTRUE;
4095 }
4096 else{
6c84d371 4097 // vESDeNegNoJPsi.push_back(vESDeNegTemp[iNeg]);
4098 new((*vESDeNegNoJPsi)[vESDeNegNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDeNegTemp->At(iNeg)));
d7d7e825 4099 // cout<<"ESD No Positivo JPsi "<<endl;
4100 }
a0b94e5c 4101
d7d7e825 4102 }
4103 }
a0b94e5c 4104
6c84d371 4105 for(Int_t iPos=0; iPos < vESDePosTemp->GetEntriesFast(); iPos++){
4106 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetPdgCode() == -11)
4107 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0) > -1){
4108 Int_t labelMother = fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0);
d7d7e825 4109 TParticle* partMother = fStack ->Particle(labelMother);
4110 if (partMother->GetPdgCode() == 111){
4111 iePosPi0 = iPos;
4112 ePosPi0 = kTRUE;
4113 }
4114 if(partMother->GetPdgCode() == 443){ //Mother JPsi
4115 fHistograms->FillTable("Table_Electrons",15);
4116 iePosJPsi = iPos;
4117 ePosJPsi = kTRUE;
4118 }
4119 else{
6c84d371 4120 // vESDePosNoJPsi.push_back(vESDePosTemp[iPos]);
4121 new((*vESDePosNoJPsi)[vESDePosNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDePosTemp->At(iPos)));
d7d7e825 4122 // cout<<"ESD No Negativo JPsi "<<endl;
4123 }
a0b94e5c 4124
d7d7e825 4125 }
4126 }
4127
4128 if( eNegJPsi && ePosJPsi ){
4129 TVector3 tempeNegV,tempePosV;
6c84d371 4130 tempeNegV.SetXYZ(((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Px(),((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Py(),((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Pz());
4131 tempePosV.SetXYZ(((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Px(),((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Py(),((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Pz());
d7d7e825 4132 fHistograms->FillTable("Table_Electrons",16);
4133 fHistograms->FillHistogram("ESD_ElectronPosNegJPsiAngle",tempeNegV.Angle(tempePosV));
6c84d371 4134 fHistograms->FillHistogram("MC_ElectronPosNegJPsiAngle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->GetLabel())),
4135 fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->GetLabel()))));
d7d7e825 4136 }
4137
4138 if( eNegPi0 && ePosPi0 ){
4139 TVector3 tempeNegV,tempePosV;
6c84d371 4140 tempeNegV.SetXYZ(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Px(),((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Py(),((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Pz());
4141 tempePosV.SetXYZ(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Px(),((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Py(),((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Pz());
d7d7e825 4142 fHistograms->FillHistogram("ESD_ElectronPosNegPi0Angle",tempeNegV.Angle(tempePosV));
6c84d371 4143 fHistograms->FillHistogram("MC_ElectronPosNegPi0Angle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->GetLabel())),
4144 fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->GetLabel()))));
d7d7e825 4145 }
a0b94e5c 4146
4147
d7d7e825 4148 FillAngle("ESD_eNegePosAngleBeforeCut",GetTLorentzVector(vESDeNegTemp),GetTLorentzVector(vESDePosTemp));
a0b94e5c 4149
6c84d371 4150 CleanWithAngleCuts(*vESDeNegTemp,*vESDePosTemp,*fKFReconstructedGammasTClone);
d7d7e825 4151
6c84d371 4152 // vector <TLorentzVector> vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectron);
4153 // vector <TLorentzVector> vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectron);
a0b94e5c 4154
6c84d371 4155 TClonesArray vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectronTClone);
4156 TClonesArray vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectronTClone);
a0b94e5c 4157
4158
d7d7e825 4159 FillAngle("ESD_eNegePosAngleAfterCut",vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 4160
4161
4162
4163
d7d7e825 4164 //FillAngle("ESD_eNegePosAngleAfterCut",CurrentTLVeNeg,CurrentTLVePos);
a0b94e5c 4165
4166
d7d7e825 4167 FillElectronInvMass("ESD_InvMass_ePluseMinus",vCurrentTLVeNeg,vCurrentTLVePos);
4168 FillElectronInvMass("ESD_InvMass_xPlusxMinus",GetTLorentzVector(vESDxNegTemp),GetTLorentzVector(vESDxPosTemp));
a0b94e5c 4169
4170
4171
6c84d371 4172 FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusChiC","ESD_InvMass_GammaePluseMinusChiCDiff",*fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 4173
d7d7e825 4174 FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusPi0","ESD_InvMass_GammaePluseMinusPi0Diff",
6c84d371 4175 *fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 4176
d7d7e825 4177 //BackGround
a0b94e5c 4178
d7d7e825 4179 //Like Sign e+e-
4180 ElectronBackground("ESD_ENegBackground",vCurrentTLVeNeg);
4181 ElectronBackground("ESD_EPosBackground",vCurrentTLVePos);
4182 ElectronBackground("ESD_EPosENegBackground",vCurrentTLVeNeg);
4183 ElectronBackground("ESD_EPosENegBackground",vCurrentTLVePos);
a0b94e5c 4184
d7d7e825 4185 // Like Sign e+e- no JPsi
4186 ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDeNegNoJPsi));
4187 ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDePosNoJPsi));
a0b94e5c 4188
d7d7e825 4189 //Mixed Event
a0b94e5c 4190
6c84d371 4191 if( fCurrentEventPosElectronTClone->GetEntriesFast() > 0 && fCurrentEventNegElectronTClone->GetEntriesFast() > 0 && fKFReconstructedGammasCutTClone->GetEntriesFast() > 0 ){
d7d7e825 4192 FillGammaElectronInvMass("ESD_EPosENegGammaBackgroundMX","ESD_EPosENegGammaBackgroundMXDiff",
6c84d371 4193 *fKFReconstructedGammasCutTClone,*fPreviousEventTLVNegElectronTClone,*fPreviousEventTLVPosElectronTClone);
4194 *fPreviousEventTLVNegElectronTClone = vCurrentTLVeNeg;
4195 *fPreviousEventTLVPosElectronTClone = vCurrentTLVePos;
a0b94e5c 4196
d7d7e825 4197 }
a0b94e5c 4198
d7d7e825 4199 /*
4200 //Photons P
4201 Double_t vtx[3];
4202 vtx[0]=0;vtx[1]=0;vtx[2]=0;
4203 for(UInt_t i=0;i<fKFReconstructedGammasChic.size();i++){
a0b94e5c 4204
d7d7e825 4205 // if(fMCGammaChicTempCut[i]->GetMother(0) < 0) continue;
a0b94e5c 4206
4207
4208
d7d7e825 4209 Int_t tempLabel = fStack->Particle(fMCGammaChicTempCut[i]->GetMother(0))->GetPdgCode();
4210 // cout<<" Label Pedro Gonzalez " <<tempLabel <<endl;
a0b94e5c 4211
d7d7e825 4212 // cout<<" Label Distance"<<fKFReconstructedGammasChic[i].GetDistanceFromVertex(vtx)<<endl;
a0b94e5c 4213
d7d7e825 4214 if( tempLabel == 10441 || tempLabel == 20443 || tempLabel == 445 )
a0b94e5c 4215
d7d7e825 4216 fHistograms->FillHistogram("ESD_PhotonsMomentum",fKFReconstructedGammasChic[i].GetMomentum());
a0b94e5c 4217
4218
d7d7e825 4219 }
a0b94e5c 4220
4221
d7d7e825 4222 */
1e7846f4 4223
4224
4225 vESDeNegTemp->Delete();
4226 vESDePosTemp->Delete();
4227 vESDxNegTemp->Delete();
4228 vESDxPosTemp->Delete();
4229 vESDeNegNoJPsi->Delete();
4230 vESDePosNoJPsi->Delete();
4231
4232 delete vESDeNegTemp;
4233 delete vESDePosTemp;
4234 delete vESDxNegTemp;
4235 delete vESDxPosTemp;
4236 delete vESDeNegNoJPsi;
4237 delete vESDePosNoJPsi;
d7d7e825 4238}
4239
6c84d371 4240/*
4241 void AliAnalysisTaskGammaConversion::FillAngle(TString histoName,vector <TLorentzVector> tlVeNeg, vector <TLorentzVector> tlVePos){
d7d7e825 4242 //see header file for documentation
4243 for( UInt_t iNeg=0; iNeg < tlVeNeg.size(); iNeg++){
6c84d371 4244 for (UInt_t iPos=0; iPos < tlVePos.size(); iPos++){
4245 fHistograms->FillHistogram(histoName.Data(),tlVeNeg[iNeg].Vect().Angle(tlVePos[iPos].Vect()));
4246 }
4247 }
4248 }
4249*/
4250void AliAnalysisTaskGammaConversion::FillAngle(TString histoName,TClonesArray const tlVeNeg, TClonesArray const tlVePos){
4251 //see header file for documentation
4252 for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++){
4253 for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){
4254 fHistograms->FillHistogram(histoName.Data(),((TLorentzVector*)(tlVeNeg.At(iNeg)))->Vect().Angle(((TLorentzVector*)(tlVePos.At(iPos)))->Vect()));
d7d7e825 4255 }
4256 }
4257}
6c84d371 4258void AliAnalysisTaskGammaConversion::FillElectronInvMass(TString histoName, TClonesArray const eNeg, TClonesArray const ePos){
d7d7e825 4259 //see header file for documentation
6c84d371 4260 for( Int_t n=0; n < eNeg.GetEntriesFast(); n++){
4261 TLorentzVector en = (*(TLorentzVector*)(eNeg.At(n)));
4262 for (Int_t p=0; p < ePos.GetEntriesFast(); p++){
4263 TLorentzVector ep = (*(TLorentzVector*)(ePos.At(p)));
d7d7e825 4264 TLorentzVector np = ep + en;
4265 fHistograms->FillHistogram(histoName.Data(),np.M());
4266 }
4267 }
d7d7e825 4268}
4269
6c84d371 4270void AliAnalysisTaskGammaConversion::FillGammaElectronInvMass(TString histoMass,TString histoDiff,TClonesArray const fKFGammas,
4271 TClonesArray const tlVeNeg,TClonesArray const tlVePos)
d7d7e825 4272{
4273 //see header file for documentation
a0b94e5c 4274
6c84d371 4275 for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++ ){
a0b94e5c 4276
6c84d371 4277 for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){
a0b94e5c 4278
6c84d371 4279 TLorentzVector xy = *((TLorentzVector *)(tlVePos.At(iPos))) + *((TLorentzVector *)(tlVeNeg.At(iNeg)));
a0b94e5c 4280
6c84d371 4281 for (Int_t iGam=0; iGam < fKFGammas.GetEntriesFast(); iGam++){
a0b94e5c 4282
6c84d371 4283 // AliKFParticle * gammaCandidate = &fKFGammas[iGam];
4284 AliKFParticle * gammaCandidate = (AliKFParticle *)(fKFGammas.At(iGam));
d7d7e825 4285 TLorentzVector g;
a0b94e5c 4286
d7d7e825 4287 g.SetXYZM(gammaCandidate->GetPx(),gammaCandidate->GetPy(),gammaCandidate->GetPz(),fGammaMass);
4288 TLorentzVector xyg = xy + g;
4289 fHistograms->FillHistogram(histoMass.Data(),xyg.M());
4290 fHistograms->FillHistogram(histoDiff.Data(),(xyg.M()-xy.M()));
4291 }
4292 }
4293 }
a0b94e5c 4294
d7d7e825 4295}
6c84d371 4296void AliAnalysisTaskGammaConversion::ElectronBackground(TString hBg, TClonesArray e)
d7d7e825 4297{
4298 // see header file for documentation
6c84d371 4299 for(Int_t i=0; i < e.GetEntriesFast(); i++)
d7d7e825 4300 {
6c84d371 4301 for (Int_t j=i+1; j < e.GetEntriesFast(); j++)
d7d7e825 4302 {
6c84d371 4303 TLorentzVector ee = (*(TLorentzVector*)(e.At(i))) + (*(TLorentzVector*)(e.At(j)));
a0b94e5c 4304
d7d7e825 4305 fHistograms->FillHistogram(hBg.Data(),ee.M());
4306 }
4307 }
4308}
4309
4310
6c84d371 4311void AliAnalysisTaskGammaConversion::CleanWithAngleCuts(TClonesArray const negativeElectrons,
4312 TClonesArray const positiveElectrons,
4313 TClonesArray const gammas){
d7d7e825 4314 // see header file for documentation
a0b94e5c 4315
6c84d371 4316 UInt_t sizeN = negativeElectrons.GetEntriesFast();
4317 UInt_t sizeP = positiveElectrons.GetEntriesFast();
4318 UInt_t sizeG = gammas.GetEntriesFast();
a0b94e5c 4319
4320
4321
d7d7e825 4322 vector <Bool_t> xNegBand(sizeN);
4323 vector <Bool_t> xPosBand(sizeP);
4324 vector <Bool_t> gammaBand(sizeG);
a0b94e5c 4325
4326
d7d7e825 4327 for(UInt_t iNeg=0; iNeg < sizeN; iNeg++) xNegBand[iNeg]=kTRUE;
4328 for(UInt_t iPos=0; iPos < sizeP; iPos++) xPosBand[iPos]=kTRUE;
4329 for(UInt_t iGam=0; iGam < sizeG; iGam++) gammaBand[iGam]=kTRUE;
d7d7e825 4330
a0b94e5c 4331
4332 for(UInt_t iPos=0; iPos < sizeP; iPos++){
4333
6c84d371 4334 Double_t aP[3];
4335 ((AliESDtrack*)(positiveElectrons.At(iPos)))->GetConstrainedPxPyPz(aP);
a0b94e5c 4336
d7d7e825 4337 TVector3 ePosV(aP[0],aP[1],aP[2]);
a0b94e5c 4338
d7d7e825 4339 for(UInt_t iNeg=0; iNeg < sizeN; iNeg++){
a0b94e5c 4340
6c84d371 4341 Double_t aN[3];
4342 ((AliESDtrack*)(negativeElectrons.At(iNeg)))->GetConstrainedPxPyPz(aN);
d7d7e825 4343 TVector3 eNegV(aN[0],aN[1],aN[2]);
a0b94e5c 4344
d7d7e825 4345 if(ePosV.Angle(eNegV) < 0.05){ //e+e- from gamma
4346 xPosBand[iPos]=kFALSE;
4347 xNegBand[iNeg]=kFALSE;
4348 }
a0b94e5c 4349
d7d7e825 4350 for(UInt_t iGam=0; iGam < sizeG; iGam++){
6c84d371 4351 AliKFParticle* gammaCandidate = (AliKFParticle*)gammas.At(iGam);
d7d7e825 4352 TVector3 gammaCandidateVector(gammaCandidate->Px(),gammaCandidate->Py(),gammaCandidate->Pz());
4353 if(ePosV.Angle(gammaCandidateVector) < 0.05 || eNegV.Angle(gammaCandidateVector) < 0.05)
4354 gammaBand[iGam]=kFALSE;
4355 }
4356 }
4357 }
a0b94e5c 4358
4359
4360
4361
d7d7e825 4362 for(UInt_t iPos=0; iPos < sizeP; iPos++){
4363 if(xPosBand[iPos]){
6c84d371 4364 new((*fCurrentEventPosElectronTClone)[fCurrentEventPosElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(positiveElectrons.At(iPos))));
4365 // fCurrentEventPosElectron.push_back(positiveElectrons[iPos]);
d7d7e825 4366 }
4367 }
4368 for(UInt_t iNeg=0;iNeg < sizeN; iNeg++){
4369 if(xNegBand[iNeg]){
6c84d371 4370 new((*fCurrentEventNegElectronTClone)[fCurrentEventNegElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(negativeElectrons.At(iNeg))));
4371 // fCurrentEventNegElectron.push_back(negativeElectrons[iNeg]);
d7d7e825 4372 }
4373 }
4374 for(UInt_t iGam=0; iGam < sizeG; iGam++){
4375 if(gammaBand[iGam]){
6c84d371 4376 new((*fKFReconstructedGammasCutTClone)[fKFReconstructedGammasCutTClone->GetEntriesFast()]) AliKFParticle((*(AliKFParticle*)(gammas.At(iGam))));
4377 //fKFReconstructedGammasCut.push_back(*(AliKFParticle*)gammas->At(iGam));
d7d7e825 4378 }
4379 }
4380}
4381
4382
77ac6f3e 4383void AliAnalysisTaskGammaConversion::GetPID(const AliESDtrack *track, Stat_t &pid, Stat_t &weight)
d7d7e825 4384{
4385 // see header file for documentation
4386 pid = -1;
4387 weight = -1;
a0b94e5c 4388
d7d7e825 4389 double wpart[5];
4390 double wpartbayes[5];
a0b94e5c 4391
d7d7e825 4392 //get probability of the diffenrent particle types
4393 track->GetESDpid(wpart);
a0b94e5c 4394
d7d7e825 4395 // Tentative particle type "concentrations"
4396 double c[5]={0.01, 0.01, 0.85, 0.10, 0.05};
a0b94e5c 4397
d7d7e825 4398 //Bayes' formula
4399 double rcc = 0.;
4400 for (int i = 0; i < 5; i++)
4401 {
4402 rcc+=(c[i] * wpart[i]);
4403 }
a0b94e5c 4404
4405
4406
d7d7e825 4407 for (int i=0; i<5; i++) {
4a6157dc 4408 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 4409 wpartbayes[i] = c[i] * wpart[i] / rcc;
4410 }
4411 }
a0b94e5c 4412
4413
4414
d7d7e825 4415 Float_t max=0.;
4416 int ipid=-1;
4417 //find most probable particle in ESD pid
4418 //0:Electron - 1:Muon - 2:Pion - 3:Kaon - 4:Proton
4419 for (int i = 0; i < 5; i++)
4420 {
4421 if (wpartbayes[i] > max)
4422 {
a0b94e5c 4423 ipid = i;
4424 max = wpartbayes[i];
d7d7e825 4425 }
4426 }
a0b94e5c 4427
d7d7e825 4428 pid = ipid;
4429 weight = max;
4430}
77ac6f3e 4431double AliAnalysisTaskGammaConversion::GetSigmaToVertex(const AliESDtrack* t)
d7d7e825 4432{
4433 // Calculates the number of sigma to the vertex.
a0b94e5c 4434
d7d7e825 4435 Float_t b[2];
4436 Float_t bRes[2];
4437 Float_t bCov[3];
4438 t->GetImpactParameters(b,bCov);
4439 if (bCov[0]<=0 || bCov[2]<=0) {
4440 AliDebug(1, "Estimated b resolution lower or equal zero!");
4441 bCov[0]=0; bCov[2]=0;
4442 }
4443 bRes[0] = TMath::Sqrt(bCov[0]);
4444 bRes[1] = TMath::Sqrt(bCov[2]);
a0b94e5c 4445
d7d7e825 4446 // -----------------------------------
4447 // How to get to a n-sigma cut?
4448 //
4449 // The accumulated statistics from 0 to d is
4450 //
4451 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
4452 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
4453 //
4454 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
4455 // Can this be expressed in a different way?
a0b94e5c 4456
d7d7e825 4457 if (bRes[0] == 0 || bRes[1] ==0)
4458 return -1;
a0b94e5c 4459
d7d7e825 4460 double d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
a0b94e5c 4461
d7d7e825 4462 // stupid rounding problem screws up everything:
4463 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
4464 if (TMath::Exp(-d * d / 2) < 1e-10)
4465 return 1000;
a0b94e5c 4466
4467
d7d7e825 4468 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
4469 return d;
4470}
6c84d371 4471
4472//vector <TLorentzVector> AliAnalysisTaskGammaConversion::GetTLorentzVector(vector <AliESDtrack*> esdTrack){
4473TClonesArray AliAnalysisTaskGammaConversion::GetTLorentzVector(TClonesArray *const esdTrack){
01b7fdcc 4474 //Return TLoresntz vector of track?
6c84d371 4475 // vector <TLorentzVector> tlVtrack(0);
4476 TClonesArray array("TLorentzVector",0);
a0b94e5c 4477
6c84d371 4478 for(Int_t itrack=0; itrack < esdTrack->GetEntriesFast(); itrack++){
4479 double p[3];
4480 //esdTrack[itrack]->GetConstrainedPxPyPz(p);
4481 ((AliESDtrack*)(esdTrack->At(itrack)))->GetConstrainedPxPyPz(p);
d7d7e825 4482 TLorentzVector currentTrack;
01b7fdcc 4483 currentTrack.SetXYZM(p[0],p[1],p[2],fElectronMass);
6c84d371 4484 new((array)[array.GetEntriesFast()]) TLorentzVector(currentTrack);
4485 // tlVtrack.push_back(currentTrack);
d7d7e825 4486 }
a0b94e5c 4487
6c84d371 4488 return array;
a0b94e5c 4489
6c84d371 4490 // return tlVtrack;
d7d7e825 4491}
77ac6f3e 4492Int_t AliAnalysisTaskGammaConversion::GetProcessType(const AliMCEvent * mcEvt) {
3c45d101 4493
4494 // Determine if the event was generated with pythia or phojet and return the process type
4495
4496 // Check if mcEvt is fine
239a56d1 4497 if (!mcEvt) { // coverty does not allow this, the check is done elsewhere
3c45d101 4498 AliFatal("NULL mc event");
239a56d1 4499 return -1;
3c45d101 4500 }
4501
4502 // Determine if it was a pythia or phojet header, and return the correct process type
4503 AliGenPythiaEventHeader * headPy = 0;
4504 AliGenDPMjetEventHeader * headPho = 0;
4505 AliGenEventHeader * htmp = mcEvt->GenEventHeader();
4506 if(!htmp) {
4507 AliFatal("Cannot Get MC Header!!");
239a56d1 4508 return -1;
3c45d101 4509 }
4510 if( TString(htmp->IsA()->GetName()) == "AliGenPythiaEventHeader") {
4511 headPy = (AliGenPythiaEventHeader*) htmp;
4512 } else if (TString(htmp->IsA()->GetName()) == "AliGenDPMjetEventHeader") {
4513 headPho = (AliGenDPMjetEventHeader*) htmp;
4514 } else {
4515 AliError("Unknown header");
4516 }
4517
4518 // Determine process type
4519 if(headPy) {
4520 if(headPy->ProcessType() == 92 || headPy->ProcessType() == 93) {
4521 // single difractive
4522 return kProcSD;
4523 } else if (headPy->ProcessType() == 94) {
4524 // double diffractive
4525 return kProcDD;
4526 }
4527 else if(headPy->ProcessType() != 92 && headPy->ProcessType() != 93 && headPy->ProcessType() != 94) {
4528 // non difractive
4529 return kProcND;
4530 }
4531 } else if (headPho) {
4532 if(headPho->ProcessType() == 5 || headPho->ProcessType() == 6 ) {
4533 // single difractive
4534 return kProcSD;
4535 } else if (headPho->ProcessType() == 7) {
4536 // double diffractive
4537 return kProcDD;
4538 } else if(headPho->ProcessType() != 5 && headPho->ProcessType() != 6 && headPho->ProcessType() != 7 ) {
4539 // non difractive
4540 return kProcND;
4541 }
4542 }
4543
4544
4545 // no process type found?
4546 AliError(Form("Unknown header: %s", htmp->IsA()->GetName()));
4547 return kProcUnknown;
4548}
6746e1e1 4549
4550
4551Int_t AliAnalysisTaskGammaConversion::CalculateMultiplicityBin(){
4552 // Get Centrality bin
4553
4554 Int_t multiplicity = 0;
4555
4556 if ( fUseMultiplicity == 1 ) {
4557
4558 if (fMultiplicity>= 0 && fMultiplicity<= 5) multiplicity=1;
4559 if (fMultiplicity>= 6 && fMultiplicity<= 9) multiplicity=2;
4560 if (fMultiplicity>=10 && fMultiplicity<=14) multiplicity=3;
4561 if (fMultiplicity>=15 && fMultiplicity<=22) multiplicity=4;
4562 if (fMultiplicity>=23 ) multiplicity=5;
4563
4564 }
4565 return multiplicity;
4566}