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