]> 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
Little change needed for the V0 studies (M. Ivanov)
[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());
e0b4c21c 2039 gamma.SetE(gammakf->GetE());
6441e967 2040 gamma.SetLabel1(fElectronv1[gammaIndex]);
2041 gamma.SetLabel2(fElectronv2[gammaIndex]);
2042 gamma.SetChi2(gammakf->Chi2());
2043 gamma.SetE(gammakf->E());
2044 gamma.SetESDEvent(dynamic_cast<AliESDEvent*>(InputEvent()));
2045 AddToAODBranch(fAODGamma, gamma);
2046 }
a0b94e5c 2047
6441e967 2048 }
d7d7e825 2049}
6272370b 2050void AliAnalysisTaskGammaConversion::ProcessGammasForOmegaMesonAnalysis(){
2051 // omega meson analysis pi0+gamma decay
2052 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
48682642 2053 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
6272370b 2054 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
48682642 2055
6272370b 2056 AliKFParticle * omegaCandidateGammaDaughter = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2057 if(fGammav1[firstPi0Index]==firstGammaIndex || fGammav2[firstPi0Index]==firstGammaIndex){
2058 continue;
2059 }
2060
48682642 2061 AliKFParticle omegaCandidate(*omegaCandidatePi0Daughter,*omegaCandidateGammaDaughter);
6272370b 2062 Double_t massOmegaCandidate = 0.;
2063 Double_t widthOmegaCandidate = 0.;
2064
48682642 2065 omegaCandidate.GetMass(massOmegaCandidate,widthOmegaCandidate);
6272370b 2066
e9aea39f 2067 if ( massOmegaCandidate > 733 && massOmegaCandidate < 833 ) {
6441e967 2068 //AddOmegaToAOD(&omegaCandidate, massOmegaCandidate, firstPi0Index, firstGammaIndex);
c59360eb 2069 }
2070
48682642 2071 fHistograms->FillHistogram("ESD_Omega_InvMass_vs_Pt",massOmegaCandidate ,omegaCandidate.GetPt());
6272370b 2072 fHistograms->FillHistogram("ESD_Omega_InvMass",massOmegaCandidate);
48682642 2073
2074 //delete omegaCandidate;
2075
2076 }// end of omega reconstruction in pi0+gamma channel
2077
2078 if(fDoJet == kTRUE){
2079 AliKFParticle* negPiKF=NULL;
2080 AliKFParticle* posPiKF=NULL;
2081
2082 // look at the pi+pi+pi0 channel
2083 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
2084 AliESDtrack* posTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
2085 if (posTrack->GetSign()<0) continue;
c0d9051e 2086 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion))>2.) continue;
48682642 2087 if (posPiKF) delete posPiKF; posPiKF=NULL;
2088 posPiKF = new AliKFParticle( *(posTrack) ,211);
2089
2090 for(Int_t jCh=0;jCh<fChargedParticles->GetEntriesFast();jCh++){
2091 AliESDtrack* negTrack = (AliESDtrack*)(fChargedParticles->At(jCh));
2092 if( negTrack->GetSign()>0) continue;
c0d9051e 2093 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion))>2.) continue;
48682642 2094 if (negPiKF) delete negPiKF; negPiKF=NULL;
2095 negPiKF = new AliKFParticle( *(negTrack) ,-211);
2096 AliKFParticle omegaCandidatePipPinPi0(*omegaCandidatePi0Daughter,*posPiKF,*negPiKF);
2097 Double_t massOmegaCandidatePipPinPi0 = 0.;
2098 Double_t widthOmegaCandidatePipPinPi0 = 0.;
2099
2100 omegaCandidatePipPinPi0.GetMass(massOmegaCandidatePipPinPi0,widthOmegaCandidatePipPinPi0);
ef2e2748 2101
2102 if ( massOmegaCandidatePipPinPi0 > 733 && massOmegaCandidatePipPinPi0 < 833 ) {
6441e967 2103 // AddOmegaToAOD(&omegaCandidatePipPinPi0, massOmegaCandidatePipPinPi0, -1, -1);
ef2e2748 2104 }
48682642 2105
2106 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass_vs_Pt",massOmegaCandidatePipPinPi0 ,omegaCandidatePipPinPi0.GetPt());
2107 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass",massOmegaCandidatePipPinPi0);
2108
2109 // delete omegaCandidatePipPinPi0;
2110 }
2111 }
48682642 2112
45936300 2113 if (posPiKF) delete posPiKF; posPiKF=NULL; if (negPiKF) delete negPiKF; negPiKF=NULL;
48682642 2114
45936300 2115 } // checking ig gammajet because in that case the chargedparticle list is created
6272370b 2116
6272370b 2117 }
48682642 2118
2119 if(fCalculateBackground){
5ce758b0 2120
2121 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2122
2123 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
2124 Int_t mbin = 0;
2125 if(fUseTrackMultiplicityForBG == kTRUE){
2126 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2127 }
2128 else{
2129 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2130 }
2131
2132 AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
2133
48682642 2134 // Background calculation for the omega
2135 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
5ce758b0 2136 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);
2137
2138 if(fMoveParticleAccordingToVertex == kTRUE){
2139 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
2140 }
48682642 2141 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2142 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
5ce758b0 2143
2144 if(fMoveParticleAccordingToVertex == kTRUE){
2145 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2146 }
2147
48682642 2148 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
2149 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
2150 AliKFParticle * omegaBckCandidate = new AliKFParticle(*omegaCandidatePi0Daughter,previousGoodV0);
2151 Double_t massOmegaBckCandidate = 0.;
2152 Double_t widthOmegaBckCandidate = 0.;
2153
2154 omegaBckCandidate->GetMass(massOmegaBckCandidate,widthOmegaBckCandidate);
ef2e2748 2155
2156
48682642 2157 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass_vs_Pt",massOmegaBckCandidate ,omegaBckCandidate->GetPt());
2158 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass",massOmegaBckCandidate);
2159
2160 delete omegaBckCandidate;
2161 }
2162 }
2163 }
2164 } // end of checking if background calculation is available
6272370b 2165}
d7d7e825 2166
04bf4381 2167
77ac6f3e 2168void AliAnalysisTaskGammaConversion::AddOmegaToAOD(const AliKFParticle * const omegakf, Double_t mass, Int_t omegaDaughter, Int_t gammaDaughter) {
04bf4381 2169 //See header file for documentation
2170 AliGammaConversionAODObject omega;
77ac6f3e 2171 omega.SetPx(omegakf->GetPx());
2172 omega.SetPy(omegakf->GetPy());
2173 omega.SetPz(omegakf->GetPz());
04bf4381 2174 omega.SetChi2(omegakf->GetChi2());
77ac6f3e 2175 omega.SetE(omegakf->GetE());
04bf4381 2176 omega.SetIMass(mass);
2177 omega.SetLabel1(omegaDaughter);
6441e967 2178 // //dynamic_cast<AliAODPWG4Particle*>(fAODBranch->At(daughter1))->SetTagged(kTRUE);
04bf4381 2179 omega.SetLabel2(gammaDaughter);
6441e967 2180 AddToAODBranch(fAODOmega, omega);
04bf4381 2181}
2182
2183
2184
d7d7e825 2185void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
2186 // see header file for documentation
2187
6c84d371 2188 // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
2189 // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
037dc2db 2190
5ce758b0 2191 fESDEvent = fV0Reader->GetESDEvent();
2192
037dc2db 2193 if(fKFReconstructedGammasTClone->GetEntriesFast()>fV0Reader->GetNumberOfV0s()){
2194 cout<<"Warning, number of entries in the tclone is bigger than number of v0s"<<endl;
2195 }
2196
6c84d371 2197 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2198 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
d7d7e825 2199
6c84d371 2200 // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
2201 // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
a0b94e5c 2202
6c84d371 2203 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2204 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
a0b94e5c 2205
d7d7e825 2206 if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
2207 continue;
2208 }
2209 if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
2210 continue;
2211 }
a0b94e5c 2212
d7d7e825 2213 AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
2214
2215 Double_t massTwoGammaCandidate = 0.;
2216 Double_t widthTwoGammaCandidate = 0.;
2217 Double_t chi2TwoGammaCandidate =10000.;
2218 twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
d707e3cf 2219 // if(twoGammaCandidate->GetNDF()>0){
2220 // chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
2221 chi2TwoGammaCandidate = twoGammaCandidate->GetChi2();
d7d7e825 2222
70ef88b5 2223 fHistograms->FillHistogram("ESD_Mother_Chi2",chi2TwoGammaCandidate);
2e2da371 2224 if((chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
d7d7e825 2225
2226 TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
2227 TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
2228
2229 Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
2230 Double_t rapidity;
96ade8ef 2231 if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() <= 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() <= 0){
2232 cout << "Error: |Pz| > E !!!! " << endl;
d7d7e825 2233 rapidity=0;
2234 }
2235 else{
2236 rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
2237 }
dc2883e4 2238
2239 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut()){
2240 delete twoGammaCandidate;
2241 continue; // rapidity cut
2242 }
2243
2244
48682642 2245 Double_t alfa=0.0;
2246 if( (twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()) != 0){
2247 alfa=TMath::Abs((twoGammaDecayCandidateDaughter0->GetE()-twoGammaDecayCandidateDaughter1->GetE())
2248 /(twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()));
2249 }
d7d7e825 2250
1e7846f4 2251 if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
2252 delete twoGammaCandidate;
2253 continue; // minimum opening angle to avoid using ghosttracks
2254 }
2255
67381a40 2256 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 2257 fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
2258 fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
2259 fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
2260 fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
2261 fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
2262 fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
2263 fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
2264 fHistograms->FillHistogram("ESD_Mother_alfa", alfa);
9da64c46 2265 if(massTwoGammaCandidate>0.1 && massTwoGammaCandidate<0.15){
2266 fHistograms->FillHistogram("ESD_Mother_alfa_Pi0", alfa);
2267 }
dcdc851f 2268 if(massTwoGammaCandidate>0.5 && massTwoGammaCandidate<0.57){
2269 fHistograms->FillHistogram("ESD_Mother_alfa_Eta", alfa);
2270 }
2271
6de3471d 2272 fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
2273 fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
2274 fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
2275 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2276 fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
9c1cb6f7 2277 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2278 }
3c45d101 2279 if(alfa<0.1){
2280 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_E_alpha",massTwoGammaCandidate ,twoGammaCandidate->GetE());
2281 }
6de3471d 2282
9640a3d1 2283
9c1cb6f7 2284 if(fCalculateBackground){
5ce758b0 2285 /* Kenneth, just for testing*/
2286 AliGammaConversionBGHandler * bgHandlerTest = fV0Reader->GetBGHandler();
2287
2288 Int_t zbin= bgHandlerTest->GetZBinIndex(fV0Reader->GetVertexZ());
2289 Int_t mbin=0;
2290 Int_t multKAA=0;
2291 if(fUseTrackMultiplicityForBG == kTRUE){
2292 multKAA=fV0Reader->CountESDTracks();
2293 mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2294 }
2295 else{// means we use #v0s for multiplicity
2296 multKAA=fV0Reader->GetNGoodV0s();
2297 mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2298 }
2299 // cout<<"Filling bin number "<<zbin<<" and "<<mbin<<endl;
2300 // cout<<"Corresponding to z = "<<fV0Reader->GetVertexZ()<<" and m = "<<multKAA<<endl;
6de3471d 2301 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2302 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass",zbin,mbin),massTwoGammaCandidate);
2303 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass_vs_Pt",zbin,mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2304 /* end Kenneth, just for testing*/
2305 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2306 }
5ce758b0 2307 }
2308 /* if(fCalculateBackground){
9c1cb6f7 2309 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2310 Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2311 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
5ce758b0 2312 }*/
ebcfaa7e 2313 // if(fDoNeutralMesonV0MCCheck){
2314 if(fDoMCTruth){
037dc2db 2315 //Kenneth: Checking the eta of the gamma to check the difference between 0.9 and 1.2
2316 Int_t indexKF1 = fKFReconstructedGammasV0Index.at(firstGammaIndex);
2317 if(indexKF1<fV0Reader->GetNumberOfV0s()){
2318 fV0Reader->GetV0(indexKF1);//updates to the correct v0
2319 Double_t eta1 = fV0Reader->GetMotherCandidateEta();
2320 Bool_t isRealPi0=kFALSE;
10e3319b 2321 Bool_t isRealEta=kFALSE;
037dc2db 2322 Int_t gamma1MotherLabel=-1;
2323 if(fV0Reader->HasSameMCMother() == kTRUE){
2324 //cout<<"This v0 is a real v0!!!!"<<endl;
2325 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2326 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
2327 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
2328 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
2329 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2330 gamma1MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
2331 }
2332 }
2333 }
2334 }
2335 Int_t indexKF2 = fKFReconstructedGammasV0Index.at(secondGammaIndex);
2336 if(indexKF1 == indexKF2){
2337 cout<<"index of the two KF particles are the same.... should not happen"<<endl;
2338 }
2339 if(indexKF2<fV0Reader->GetNumberOfV0s()){
2340 fV0Reader->GetV0(indexKF2);
2341 Double_t eta2 = fV0Reader->GetMotherCandidateEta();
2342 Int_t gamma2MotherLabel=-1;
2343 if(fV0Reader->HasSameMCMother() == kTRUE){
2344 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2345 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
2346 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
2347 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
2348 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2349 gamma2MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
2350 }
2351 }
2352 }
2353 }
2354 if(gamma1MotherLabel>=0 && gamma1MotherLabel==gamma2MotherLabel){
2355 if(fV0Reader->CheckIfPi0IsMother(gamma1MotherLabel)){
2356 isRealPi0=kTRUE;
2357 }
10e3319b 2358 if(fV0Reader->CheckIfEtaIsMother(gamma1MotherLabel)){
2359 isRealEta=kTRUE;
2360 }
2361
037dc2db 2362 }
6de3471d 2363 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2364 if(TMath::Abs(eta1)>0.9 && TMath::Abs(eta2)>0.9){
2365 // fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
2366 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2367 if(isRealPi0 || isRealEta){
2368 fHistograms->FillHistogram("ESD_TruePi0_InvMass_1212",massTwoGammaCandidate);
2369 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_1212",openingAngleTwoGammaCandidate);
2370 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2371 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2372 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2373 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2374 }
6de3471d 2375
2376 if(!isRealPi0 && !isRealEta){
2377 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
2378 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
2379 }else{
2380 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
2381 }
bd6d9fa3 2382 }
bd6d9fa3 2383
6de3471d 2384 }
2385 else if(TMath::Abs(eta1)>0.9 || TMath::Abs(eta2)>0.9){
ebcfaa7e 2386 // fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
2387 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2388
2389 if(isRealPi0 || isRealEta){
2390 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0912",massTwoGammaCandidate);
2391 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0912",openingAngleTwoGammaCandidate);
2392 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2393 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2394 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2395 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2396 }
6de3471d 2397 if(!isRealPi0 && !isRealEta){
2398 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
2399 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
2400 }else{
2401 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
2402 }
bd6d9fa3 2403 }
2404 }
6de3471d 2405 else{
ebcfaa7e 2406 // fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
2407 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2408 if(isRealPi0 || isRealEta){
2409 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0909",massTwoGammaCandidate);
2410 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0909",openingAngleTwoGammaCandidate);
2411 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2412 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2413 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2414 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
6de3471d 2415 if(gamma1MotherLabel > fV0Reader->GetMCStack()->GetNprimary()){
2416 fHistograms->FillHistogram("ESD_TruePi0Sec_InvMass",massTwoGammaCandidate);
2417 }
d707e3cf 2418 }
6de3471d 2419 if(!isRealPi0 && !isRealEta){
2420 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
2421 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
2422 }else{
2423 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
2424 }
bd6d9fa3 2425 }
037dc2db 2426 }
2427 }
2428 }
2429 }
2430 }
6de3471d 2431 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2432 if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
2433 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2434 fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
2435 }
2436
2437 if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
2438 fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
2439 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2440 }
2441 else if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 || TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
2442 fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
2443 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2444 }
2445 else{
2446 fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
2447 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2448 }
ebcfaa7e 2449
6de3471d 2450 Double_t lowMassPi0=0.1;
2451 Double_t highMassPi0=0.15;
95203575 2452 if ( ( massTwoGammaCandidate > lowMassPi0) && (massTwoGammaCandidate < highMassPi0) ){
6de3471d 2453 new((*fKFReconstructedPi0sTClone)[fKFReconstructedPi0sTClone->GetEntriesFast()]) AliKFParticle(*twoGammaCandidate);
2454 fGammav1.push_back(firstGammaIndex);
2455 fGammav2.push_back(secondGammaIndex);
95203575 2456 if( fKFCreateAOD ) {
2457 AddPionToAOD(twoGammaCandidate, massTwoGammaCandidate, firstGammaIndex, secondGammaIndex);
2458 }
6de3471d 2459 }
ebcfaa7e 2460 }
6272370b 2461
2e2da371 2462 }
6de3471d 2463 delete twoGammaCandidate;
d7d7e825 2464 }
2465 }
2466}
2467
22eae5f9 2468void AliAnalysisTaskGammaConversion::AddPionToAOD(AliKFParticle * pionkf, Double_t mass, Int_t daughter1, Int_t daughter2) {
04bf4381 2469 //See header file for documentation
22eae5f9 2470 if(fOutputAODClassName.Contains("AODObject")) {
2471 AliGammaConversionAODObject pion;
2472 pion.SetPx(pionkf->GetPx());
2473 pion.SetPy(pionkf->GetPy());
2474 pion.SetPz(pionkf->GetPz());
2475 pion.SetChi2(pionkf->GetChi2());
2476 pion.SetE(pionkf->GetE());
2477 pion.SetIMass(mass);
2478 pion.SetLabel1(daughter1);
2479 pion.SetLabel2(daughter2);
2480 AddToAODBranch(fAODPi0, pion);
2481 } else {
2482 TLorentzVector momentum(pionkf->Px(),pionkf->Py(),pionkf->Pz(), pionkf->E());
2483 AliAODConversionParticle pion = AliAODConversionParticle(momentum);
2484 pion.SetTrackLabels( daughter1, daughter2 );
2485 pion.SetChi2(pionkf->GetChi2());
2486 AddToAODBranch(fAODPi0, pion);
2487
2488 }
04bf4381 2489}
2490
2491
48682642 2492/*
77ac6f3e 2493 void AliAnalysisTaskGammaConversion::ProcessConvPHOSGammasForNeutralMesonAnalysis(){
2494
6272370b 2495 // see header file for documentation
2496 // Analyse Pi0 with one photon from Phos and 1 photon from conversions
2497
2498
2499
2500 Double_t vtx[3];
2501 vtx[0] = fV0Reader->GetPrimaryVertex()->GetX();
2502 vtx[1] = fV0Reader->GetPrimaryVertex()->GetY();
2503 vtx[2] = fV0Reader->GetPrimaryVertex()->GetZ();
2504
2505
2506 // Loop over all CaloClusters and consider only the PHOS ones:
2507 AliESDCaloCluster *clu;
2508 TLorentzVector pPHOS;
2509 TLorentzVector gammaPHOS;
2510 TLorentzVector gammaGammaConv;
2511 TLorentzVector pi0GammaConvPHOS;
2512 TLorentzVector gammaGammaConvBck;
2513 TLorentzVector pi0GammaConvPHOSBck;
2514
2515
2516 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
77ac6f3e 2517 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2518 if ( !clu->IsPHOS() || clu->E()<0.1 ) continue;
2519 clu ->GetMomentum(pPHOS ,vtx);
2520 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2521 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2522 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz(),0.);
2523 gammaPHOS.SetXYZM(pPHOS.Px(),pPHOS.Py(),pPHOS.Pz(),0.);
2524 pi0GammaConvPHOS=gammaGammaConv+gammaPHOS;
2525 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS",pi0GammaConvPHOS.M());
2526 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvPHOS",pi0GammaConvPHOS.M(),pi0GammaConvPHOS.Pt());
6272370b 2527
77ac6f3e 2528 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz());
2529 TVector3 v3D1(gammaPHOS.Px(),gammaPHOS.Py(),gammaPHOS.Pz());
2530 Double_t opanConvPHOS= v3D0.Angle(v3D1);
2531 if ( opanConvPHOS < 0.35){
2532 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanLow",pi0GammaConvPHOS.M());
2533 }else{
2534 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanHigh",pi0GammaConvPHOS.M());
2535 }
6272370b 2536
77ac6f3e 2537 }
2538
2539 // Now the LorentVector pPHOS is obtained and can be paired with the converted proton
6272370b 2540 }
2541 //==== End of the PHOS cluster selection ============
2542 TLorentzVector pEMCAL;
2543 TLorentzVector gammaEMCAL;
2544 TLorentzVector pi0GammaConvEMCAL;
2545 TLorentzVector pi0GammaConvEMCALBck;
2546
77ac6f3e 2547 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
2548 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2549 if ( !clu->IsEMCAL() || clu->E()<0.1 ) continue;
2550 if (clu->GetNCells() <= 1) continue;
2551 if ( clu->GetTOF()*1e9 < 550 || clu->GetTOF()*1e9 > 750) continue;
6272370b 2552
77ac6f3e 2553 clu ->GetMomentum(pEMCAL ,vtx);
2554 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2555 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2556 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),
2557 twoGammaDecayCandidateDaughter0->Py(),
2558 twoGammaDecayCandidateDaughter0->Pz(),0.);
2559 gammaEMCAL.SetXYZM(pEMCAL.Px(),pEMCAL.Py(),pEMCAL.Pz(),0.);
2560 pi0GammaConvEMCAL=gammaGammaConv+gammaEMCAL;
2561 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL",pi0GammaConvEMCAL.M());
2562 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL",pi0GammaConvEMCAL.M(),pi0GammaConvEMCAL.Pt());
2563 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),
2564 twoGammaDecayCandidateDaughter0->Py(),
2565 twoGammaDecayCandidateDaughter0->Pz());
2566 TVector3 v3D1(gammaEMCAL.Px(),gammaEMCAL.Py(),gammaEMCAL.Pz());
6272370b 2567
77ac6f3e 2568
2569 Double_t opanConvEMCAL= v3D0.Angle(v3D1);
2570 if ( opanConvEMCAL < 0.35){
2571 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanLow",pi0GammaConvEMCAL.M());
2572 }else{
2573 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanHigh",pi0GammaConvEMCAL.M());
2574 }
2575
2576 }
2577 if(fCalculateBackground){
2578 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2579 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
2580 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2581 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2582 gammaGammaConvBck.SetXYZM(previousGoodV0.Px(),
2583 previousGoodV0.Py(),
2584 previousGoodV0.Pz(),0.);
2585 pi0GammaConvEMCALBck=gammaGammaConvBck+gammaEMCAL;
2586 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M());
2587 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M(),
2588 pi0GammaConvEMCALBck.Pt());
2589 }
2590 }
48682642 2591
77ac6f3e 2592 // Now the LorentVector pEMCAL is obtained and can be paired with the converted proton
2593 } // end of checking if background photons are available
2594 }
6272370b 2595 //==== End of the PHOS cluster selection ============
77ac6f3e 2596
2597 }
48682642 2598*/
6272370b 2599
77ac6f3e 2600void AliAnalysisTaskGammaConversion::MoveParticleAccordingToVertex(AliKFParticle * particle,const AliGammaConversionBGHandler::GammaConversionVertex *vertex){
5ce758b0 2601 //see header file for documentation
2602
77ac6f3e 2603 Double_t dx = vertex->fX - fESDEvent->GetPrimaryVertex()->GetX();
2604 Double_t dy = vertex->fY - fESDEvent->GetPrimaryVertex()->GetY();
2605 Double_t dz = vertex->fZ - fESDEvent->GetPrimaryVertex()->GetZ();
5ce758b0 2606
2607 // cout<<"dx, dy, dz: ["<<dx<<","<<dy<<","<<dz<<"]"<<endl;
77ac6f3e 2608 particle->X() = particle->GetX() - dx;
2609 particle->Y() = particle->GetY() - dy;
2610 particle->Z() = particle->GetZ() - dz;
5ce758b0 2611}
2612
111d75df 2613void AliAnalysisTaskGammaConversion::RotateKFParticle(AliKFParticle * kfParticle,Double_t angle){
77ac6f3e 2614 // Rotate the kf particle
111d75df 2615 Double_t c = cos(angle);
2616 Double_t s = sin(angle);
2617
77ac6f3e 2618 Double_t mA[7][ 7];
111d75df 2619 for( Int_t i=0; i<7; i++ ){
2620 for( Int_t j=0; j<7; j++){
77ac6f3e 2621 mA[i][j] = 0;
111d75df 2622 }
2623 }
2624 for( int i=0; i<7; i++ ){
77ac6f3e 2625 mA[i][i] = 1;
111d75df 2626 }
77ac6f3e 2627 mA[0][0] = c; mA[0][1] = s;
2628 mA[1][0] = -s; mA[1][1] = c;
2629 mA[3][3] = c; mA[3][4] = s;
2630 mA[4][3] = -s; mA[4][4] = c;
111d75df 2631
77ac6f3e 2632 Double_t mAC[7][7];
2633 Double_t mAp[7];
111d75df 2634
2635 for( Int_t i=0; i<7; i++ ){
77ac6f3e 2636 mAp[i] = 0;
111d75df 2637 for( Int_t k=0; k<7; k++){
77ac6f3e 2638 mAp[i]+=mA[i][k] * kfParticle->GetParameter(k);
111d75df 2639 }
2640 }
2641
2642 for( Int_t i=0; i<7; i++){
77ac6f3e 2643 kfParticle->Parameter(i) = mAp[i];
111d75df 2644 }
2645
2646 for( Int_t i=0; i<7; i++ ){
2647 for( Int_t j=0; j<7; j++ ){
77ac6f3e 2648 mAC[i][j] = 0;
111d75df 2649 for( Int_t k=0; k<7; k++ ){
77ac6f3e 2650 mAC[i][j]+= mA[i][k] * kfParticle->GetCovariance(k,j);
111d75df 2651 }
2652 }
2653 }
2654
2655 for( Int_t i=0; i<7; i++ ){
2656 for( Int_t j=0; j<=i; j++ ){
2657 Double_t xx = 0;
2658 for( Int_t k=0; k<7; k++){
77ac6f3e 2659 xx+= mAC[i][k]*mA[j][k];
111d75df 2660 }
2661 kfParticle->Covariance(i,j) = xx;
2662 }
2663 }
2664}
2665
2666
d7d7e825 2667void AliAnalysisTaskGammaConversion::CalculateBackground(){
2668 // see header file for documentation
5e55d806 2669
2670
2671 TClonesArray * currentEventV0s = fV0Reader->GetCurrentEventGoodV0s();
2672
5ce758b0 2673 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2674
2675 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
2676 Int_t mbin = 0;
2677 if(fUseTrackMultiplicityForBG == kTRUE){
2678 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2679 }
2680 else{
2681 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2682 }
2683
111d75df 2684 if(fDoRotation == kTRUE){
3c45d101 2685 TRandom3 *random = new TRandom3();
111d75df 2686 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2687 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2688 for(Int_t iCurrent2=iCurrent+1;iCurrent2<currentEventV0s->GetEntriesFast();iCurrent2++){
77ac6f3e 2689 for(Int_t nRandom=0;nRandom<fNRandomEventsForBG;nRandom++){
111d75df 2690
2691 AliKFParticle currentEventGoodV02 = *(AliKFParticle *)(currentEventV0s->At(iCurrent2));
3bfbe89a 2692
2693 if(fCheckBGProbability == kTRUE){
2694 Double_t massBGprob =0.;
2695 Double_t widthBGprob = 0.;
2696 AliKFParticle *backgroundCandidateProb = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
2697 backgroundCandidateProb->GetMass(massBGprob,widthBGprob);
2698 if(massBGprob>0.1 && massBGprob<0.14){
2699 if(random->Rndm()>bgHandler->GetBGProb(zbin,mbin)){
2700 delete backgroundCandidateProb;
2701 continue;
2702 }
2703 }
2704 delete backgroundCandidateProb;
2705 }
111d75df 2706
77ac6f3e 2707 Double_t nRadiansPM = fNDegreesPMBackground*TMath::Pi()/180;
111d75df 2708
2709 Double_t rotationValue = random->Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
5ce758b0 2710
111d75df 2711 RotateKFParticle(&currentEventGoodV02,rotationValue);
5ce758b0 2712
111d75df 2713 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
5e55d806 2714
5ce758b0 2715 Double_t massBG =0.;
2716 Double_t widthBG = 0.;
2717 Double_t chi2BG =10000.;
2718 backgroundCandidate->GetMass(massBG,widthBG);
111d75df 2719
111d75df 2720 // if(backgroundCandidate->GetNDF()>0){
2721 chi2BG = backgroundCandidate->GetChi2();
2722 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
2723
2724 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2725 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2726
2727 Double_t openingAngleBG = currentEventGoodV0.GetAngle(currentEventGoodV02);
2728
2729 Double_t rapidity;
2730 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
2731 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2732
dc2883e4 2733 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2734 delete backgroundCandidate;
2735 continue; // rapidity cut
2736 }
2737
111d75df 2738
2739 Double_t alfa=0.0;
2740 if( (currentEventGoodV0.GetE()+currentEventGoodV02.GetE()) != 0){
2741 alfa=TMath::Abs((currentEventGoodV0.GetE()-currentEventGoodV02.GetE())
2742 /(currentEventGoodV0.GetE()+currentEventGoodV02.GetE()));
2743 }
2744
2745
2746 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2747 delete backgroundCandidate;
2748 continue; // minimum opening angle to avoid using ghosttracks
2749 }
2750
2751 // original
111d75df 2752 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 2753 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2754 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2755 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2756 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2757 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2758 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2759 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2760 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2761 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2762 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2763 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2764 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
111d75df 2765 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
6de3471d 2766
dcdc851f 2767 if(massBG>0.1 && massBG<0.15){
2768 fHistograms->FillHistogram("ESD_Background_alfa_Pi0", alfa);
2769 }
2770 if(massBG>0.5 && massBG<0.57){
2771 fHistograms->FillHistogram("ESD_Background_alfa_Eta", alfa);
2772 }
6de3471d 2773
2774 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
2775 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
2776 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
2777 }
2778
2779 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
2780 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
2781 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
2782 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
2783 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
2784 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
2785 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
2786 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2787 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2788 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
2789 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2790 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
2791
2792 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
2793 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2794 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
2795 }
111d75df 2796 }
3c45d101 2797 if(alfa<0.1){
2798 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
2799 }
2800
2e2da371 2801 }
111d75df 2802 //}
5ce758b0 2803 delete backgroundCandidate;
2804 }
2805 }
2806 }
111d75df 2807 delete random;
5ce758b0 2808 }
111d75df 2809 else{ // means no rotation
2810 AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
2811
2812 if(fUseTrackMultiplicityForBG){
2813 // cout<<"Using charged track multiplicity for background calculation"<<endl;
2814 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
5ce758b0 2815
111d75df 2816 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);//fV0Reader->GetBGGoodV0s(nEventsInBG);
2817
5ce758b0 2818 if(fMoveParticleAccordingToVertex == kTRUE){
2819 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
2820 }
2821
2822 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2823 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2824 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2825 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
111d75df 2826 AliKFParticle previousGoodV0test = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
5ce758b0 2827
111d75df 2828 //cout<<"Primary Vertex event: ["<<fESDEvent->GetPrimaryVertex()->GetX()<<","<<fESDEvent->GetPrimaryVertex()->GetY()<<","<<fESDEvent->GetPrimaryVertex()->GetZ()<<"]"<<endl;
2829 //cout<<"BG prim Vertex event: ["<<bgEventVertex->fX<<","<<bgEventVertex->fY<<","<<bgEventVertex->fZ<<"]"<<endl;
2830
2831 //cout<<"XYZ of particle before transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
5ce758b0 2832 if(fMoveParticleAccordingToVertex == kTRUE){
2833 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2834 }
111d75df 2835 //cout<<"XYZ of particle after transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
5ce758b0 2836
2837 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
111d75df 2838
5ce758b0 2839 Double_t massBG =0.;
2840 Double_t widthBG = 0.;
2841 Double_t chi2BG =10000.;
2842 backgroundCandidate->GetMass(massBG,widthBG);
111d75df 2843 // if(backgroundCandidate->GetNDF()>0){
2844 // chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
2845 chi2BG = backgroundCandidate->GetChi2();
2846 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
2847
2848 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2849 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2850
2851 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
2852
2853 Double_t rapidity;
2854
2855 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() <= 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() <= 0){
2856 cout << "Error: |Pz| > E !!!! " << endl;
2857 rapidity=0;
2858 } else {
2859 rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2860 }
dc2883e4 2861 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2862 delete backgroundCandidate;
2863 continue; // rapidity cut
2864 }
2865
2866
111d75df 2867 Double_t alfa=0.0;
2868 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
2869 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
2870 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
2871 }
2872
2873
2874 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2875 delete backgroundCandidate;
2876 continue; // minimum opening angle to avoid using ghosttracks
2877 }
2878
2879 // original
111d75df 2880 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 2881 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2882 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2883 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2884 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2885 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2886 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2887 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2888 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2889 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2890 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2891 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2892 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
111d75df 2893 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
6de3471d 2894
dcdc851f 2895 if(massBG>0.1 && massBG<0.15){
2896 fHistograms->FillHistogram("ESD_Background_alfa_Pi0", alfa);
2897 }
2898 if(massBG>0.5 && massBG<0.57){
2899 fHistograms->FillHistogram("ESD_Background_alfa_Eta", alfa);
2900 }
6de3471d 2901
2902 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2903 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
2904 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
2905 }
2906
2907 // test
2908 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
2909 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
2910 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
2911 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
2912 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
2913 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
2914 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
2915 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2916 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2917 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
2918 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2919 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
2920
2921 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2922 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2923 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
2924 }
2925 // }
111d75df 2926 }
3c45d101 2927 if(alfa<0.1){
2928 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
2929 }
2930
111d75df 2931 }
2932 delete backgroundCandidate;
2933 }
2934 }
2935 }
2936 }
2937 else{ // means using #V0s for multiplicity
2938
2939 // cout<<"Using the v0 multiplicity to calculate background"<<endl;
2940
2941 fHistograms->FillHistogram("ESD_Background_z_m",zbin,mbin);
2942 fHistograms->FillHistogram("ESD_Mother_multpilicityVSv0s",fV0Reader->CountESDTracks(),fV0Reader->GetNumberOfV0s());
2943
2944 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2945 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);// fV0Reader->GetBGGoodV0s(nEventsInBG);
2946 if(previousEventV0s){
2947
2948 if(fMoveParticleAccordingToVertex == kTRUE){
2949 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
2950 }
2951
2952 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2953 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2954 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2955 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2956
2957 if(fMoveParticleAccordingToVertex == kTRUE){
2958 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
5ce758b0 2959 }
111d75df 2960
2961 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
2962 Double_t massBG =0.;
2963 Double_t widthBG = 0.;
2964 Double_t chi2BG =10000.;
2965 backgroundCandidate->GetMass(massBG,widthBG);
2966 /* if(backgroundCandidate->GetNDF()>0){
2967 chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
2968 {//remember to remove
2969 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2970 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2971
2972 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
2973 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle_nochi2", openingAngleBG);
2974 }
2975 */
2e2da371 2976 chi2BG = backgroundCandidate->GetChi2();
2977 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
5ce758b0 2978 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2979 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2980
2981 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
2982
2983 Double_t rapidity;
2984 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
2985 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2986
dc2883e4 2987 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2988 delete backgroundCandidate;
2989 continue; // rapidity cut
2990 }
2991
2992
5ce758b0 2993 Double_t alfa=0.0;
2994 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
2995 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
2996 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
2997 }
2998
2999
3000 if(openingAngleBG < fMinOpeningAngleGhostCut ){
3001 delete backgroundCandidate;
3002 continue; // minimum opening angle to avoid using ghosttracks
3003 }
3004
5ce758b0 3005 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
6de3471d 3006 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
3007 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
3008 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
3009 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
3010 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
3011 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
3012 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
3013 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3014 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3015 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
3016 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
3017 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
3018
3019
5ce758b0 3020 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
3c45d101 3021
dcdc851f 3022 if(massBG>0.1 && massBG<0.15){
3023 fHistograms->FillHistogram("ESD_Background_alfa_Pi0", alfa);
3024 }
3025 if(massBG>0.5 && massBG<0.57){
3026 fHistograms->FillHistogram("ESD_Background_alfa_Eta", alfa);
3027 }
3028
6de3471d 3029 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
3030 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
3031 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
3032 }
3033
3034 if(massBG>0.5 && massBG<0.6){
3035 fHistograms->FillHistogram("ESD_Background_alfa_pt0506",momentumVectorbackgroundCandidate.Pt(),alfa);
3036 }
3037 if(massBG>0.3 && massBG<0.4){
3038 fHistograms->FillHistogram("ESD_Background_alfa_pt0304",momentumVectorbackgroundCandidate.Pt(),alfa);
3039 }
3040
3041 // test
3042 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
3043 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
3044 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
3045 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
3046 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
3047 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
3048 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
3049 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
3050 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
3051 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
3052 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3053 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
3054
3055 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
3056 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
3057 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
3058 }
5ce758b0 3059 }
5ce758b0 3060
6de3471d 3061 if(alfa<0.1){
3062 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
5ce758b0 3063 }
2e2da371 3064 // }
5ce758b0 3065 }
111d75df 3066 delete backgroundCandidate;
3067 }
5ce758b0 3068 }
3069 }
d7d7e825 3070 }
111d75df 3071 } // end else (means use #v0s as multiplicity)
3072 } // end no rotation
d7d7e825 3073}
3074
3075
d7d7e825 3076void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){
3077 //ProcessGammasForGammaJetAnalysis
a0b94e5c 3078
d7d7e825 3079 Double_t distIsoMin;
a0b94e5c 3080
d7d7e825 3081 CreateListOfChargedParticles();
a0b94e5c 3082
3083
6c84d371 3084 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
3085 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
3086 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
01b7fdcc 3087 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 3088
01b7fdcc 3089 if( momentumVectorCurrentGamma.Pt()> fMinPtForGammaJet){
d7d7e825 3090 distIsoMin=GetMinimumDistanceToCharge(gammaIndex);
a0b94e5c 3091
d7d7e825 3092 if (distIsoMin > fMinIsoConeSize && fLeadingChargedIndex>=0){
01b7fdcc 3093 CalculateJetCone(gammaIndex);
d7d7e825 3094 }
3095 }
3096 }
3097}
3098
6272370b 3099//____________________________________________________________________
77ac6f3e 3100Bool_t AliAnalysisTaskGammaConversion::IsGoodImpPar(const AliESDtrack *const track)
6272370b 3101{
3102//
3103// check whether particle has good DCAr(Pt) impact
3104// parameter. Only for TPC+ITS tracks (7*sigma cut)
3105// Origin: Andrea Dainese
3106//
3107
3108Float_t d0z0[2],covd0z0[3];
3109track->GetImpactParameters(d0z0,covd0z0);
3110Float_t sigma= 0.0050+0.0060/TMath::Power(track->Pt(),0.9);
3111Float_t d0max = 7.*sigma;
3112if(TMath::Abs(d0z0[0]) < d0max) return kTRUE;
3113
3114return kFALSE;
3115}
3116
3117
d7d7e825 3118void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){
3119 // CreateListOfChargedParticles
a0b94e5c 3120
d7d7e825 3121 fESDEvent = fV0Reader->GetESDEvent();
037dc2db 3122 Int_t numberOfESDTracks=0;
d7d7e825 3123 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
3124 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 3125
d7d7e825 3126 if(!curTrack){
3127 continue;
3128 }
d707e3cf 3129 // Not needed if Standard function used.
3130// if(!IsGoodImpPar(curTrack)){
3131// continue;
3132// }
a0b94e5c 3133
d7d7e825 3134 if(fEsdTrackCuts->AcceptTrack(curTrack) ){
6c84d371 3135 new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack);
3136 // fChargedParticles.push_back(curTrack);
d7d7e825 3137 fChargedParticlesId.push_back(iTracks);
037dc2db 3138 numberOfESDTracks++;
d7d7e825 3139 }
3140 }
6746e1e1 3141// Moved to UserExec using CountAcceptedTracks function. runjet is not needed by default
3142// fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",numberOfESDTracks);
3143// cout<<"esdtracks::"<< numberOfESDTracks<<endl;
3144// if (fV0Reader->GetNumberOfContributorsVtx()>=1){
3145// fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",numberOfESDTracks);
3146// }
d7d7e825 3147}
9c1cb6f7 3148void AliAnalysisTaskGammaConversion::RecalculateV0ForGamma(){
77ac6f3e 3149 //recalculates v0 for gamma
3150
9c1cb6f7 3151 Double_t massE=0.00051099892;
3152 TLorentzVector curElecPos;
3153 TLorentzVector curElecNeg;
3154 TLorentzVector curGamma;
3155
3156 TLorentzVector curGammaAt;
3157 TLorentzVector curElecPosAt;
3158 TLorentzVector curElecNegAt;
3159 AliKFVertex primVtxGamma(*(fESDEvent->GetPrimaryVertex()));
3160 AliKFVertex primVtxImprovedGamma = primVtxGamma;
3161
3162 const AliESDVertex *vtxT3D=fESDEvent->GetPrimaryVertex();
3163
3164 Double_t xPrimaryVertex=vtxT3D->GetXv();
3165 Double_t yPrimaryVertex=vtxT3D->GetYv();
3166 Double_t zPrimaryVertex=vtxT3D->GetZv();
70ef88b5 3167 // Float_t primvertex[3]={xPrimaryVertex,yPrimaryVertex,zPrimaryVertex};
9c1cb6f7 3168
3169 Float_t nsigmaTPCtrackPos;
3170 Float_t nsigmaTPCtrackNeg;
3171 Float_t nsigmaTPCtrackPosToPion;
3172 Float_t nsigmaTPCtrackNegToPion;
3173 AliKFParticle* negKF=NULL;
3174 AliKFParticle* posKF=NULL;
3175
3176 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
3177 AliESDtrack* posTrack = fESDEvent->GetTrack(iTracks);
3178 if(!posTrack){
3179 continue;
3180 }
3181 if (posKF) delete posKF; posKF=NULL;
3182 if(posTrack->GetSign()<0) continue;
3183 if(!(posTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
3184 if(posTrack->GetKinkIndex(0)>0 ) continue;
3185 if(posTrack->GetNcls(1)<50)continue;
3186 Double_t momPos[3];
3187 // posTrack->GetConstrainedPxPyPz(momPos);
3188 posTrack->GetPxPyPz(momPos);
3189 AliESDtrack *ptrk=fESDEvent->GetTrack(iTracks);
3190 curElecPos.SetXYZM(momPos[0],momPos[1],momPos[2],massE);
3191 if(TMath::Abs(curElecPos.Eta())<0.9) continue;
3192 posKF = new AliKFParticle( *(posTrack),-11);
3193
3194 nsigmaTPCtrackPos = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kElectron);
3195 nsigmaTPCtrackPosToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion);
3196
3197 if ( nsigmaTPCtrackPos>5.|| nsigmaTPCtrackPos<-2.){
3198 continue;
3199 }
3200
3201 if(pow((momPos[0]*momPos[0]+momPos[1]*momPos[1]+momPos[2]*momPos[2]),0.5)>0.5 && nsigmaTPCtrackPosToPion<1){
3202 continue;
3203 }
3204
3205
3206
3207 for(Int_t jTracks = 0; jTracks < fESDEvent->GetNumberOfTracks(); jTracks++){
3208 AliESDtrack* negTrack = fESDEvent->GetTrack(jTracks);
3209 if(!negTrack){
3210 continue;
3211 }
3212 if (negKF) delete negKF; negKF=NULL;
3213 if(negTrack->GetSign()>0) continue;
3214 if(!(negTrack->GetStatus() & AliESDtrack::kTPCrefit))continue;
3215 if(negTrack->GetKinkIndex(0)>0 ) continue;
3216 if(negTrack->GetNcls(1)<50)continue;
3217 Double_t momNeg[3];
3218 // negTrack->GetConstrainedPxPyPz(momNeg);
3219 negTrack->GetPxPyPz(momNeg);
3220
3221 nsigmaTPCtrackNeg = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kElectron);
3222 nsigmaTPCtrackNegToPion = fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion);
3223 if ( nsigmaTPCtrackNeg>5. || nsigmaTPCtrackNeg<-2.){
3224 continue;
3225 }
3226 if(pow((momNeg[0]*momNeg[0]+momNeg[1]*momNeg[1]+momNeg[2]*momNeg[2]),0.5)>0.5 && nsigmaTPCtrackNegToPion<1){
3227 continue;
3228 }
3229 AliESDtrack *ntrk=fESDEvent->GetTrack(jTracks);
3230 curElecNeg.SetXYZM(momNeg[0],momNeg[1],momNeg[2],massE);
3231 if(TMath::Abs(curElecNeg.Eta())<0.9) continue;
3232 negKF = new AliKFParticle( *(negTrack) ,11);
3233
3234 Double_t b=fESDEvent->GetMagneticField();
3235 Double_t xn, xp, dca=ntrk->GetDCA(ptrk,b,xn,xp);
3236 AliExternalTrackParam nt(*ntrk), pt(*ptrk);
3237 nt.PropagateTo(xn,b); pt.PropagateTo(xp,b);
3238
3239
3240 //--- Like in ITSV0Finder
3241 AliExternalTrackParam ntAt0(*ntrk), ptAt0(*ptrk);
3242 Double_t xxP,yyP,alphaP;
3243 Double_t rP[3];
3244
3245 // if (!ptAt0.GetGlobalXYZat(ptAt0->GetX(),xxP,yyP,zzP)) continue;
3246 if (!ptAt0.GetXYZAt(ptAt0.GetX(),b,rP)) continue;
3247 xxP=rP[0];
3248 yyP=rP[1];
3249 alphaP = TMath::ATan2(yyP,xxP);
3250
3251
3252 ptAt0.Propagate(alphaP,0,b);
3253 Float_t ptfacP = (1.+100.*TMath::Abs(ptAt0.GetC(b)));
3254
3255 // Double_t distP = ptAt0.GetY();
3256 Double_t normP = ptfacP*TMath::Sqrt(ptAt0.GetSigmaY2());
3257 Double_t normdist0P = TMath::Abs(ptAt0.GetY()/normP);
3258 Double_t normdist1P = TMath::Abs((ptAt0.GetZ()-zPrimaryVertex)/(ptfacP*TMath::Sqrt(ptAt0.GetSigmaZ2())));
3259 Double_t normdistP = TMath::Sqrt(normdist0P*normdist0P+normdist1P*normdist1P);
3260
3261
3262 Double_t xxN,yyN,alphaN;
3263 Double_t rN[3];
3264 // if (!ntAt0.GetGlobalXYZat(ntAt0->GetX(),xxN,yyN,zzN)) continue;
3265 if (!ntAt0.GetXYZAt(ntAt0.GetX(),b,rN)) continue;
3266 xxN=rN[0];
3267 yyN=rN[1];
3268
3269 alphaN = TMath::ATan2(yyN,xxN);
3270
3271 ntAt0.Propagate(alphaN,0,b);
3272
3273 Float_t ptfacN = (1.+100.*TMath::Abs(ntAt0.GetC(b)));
3274 // Double_t distN = ntAt0.GetY();
3275 Double_t normN = ptfacN*TMath::Sqrt(ntAt0.GetSigmaY2());
3276 Double_t normdist0N = TMath::Abs(ntAt0.GetY()/normN);
3277 Double_t normdist1N = TMath::Abs((ntAt0.GetZ()-zPrimaryVertex)/(ptfacN*TMath::Sqrt(ntAt0.GetSigmaZ2())));
3278 Double_t normdistN = TMath::Sqrt(normdist0N*normdist0N+normdist1N*normdist1N);
3279
3280 //-----------------------------
3281
3282 Double_t momNegAt[3];
3283 nt.GetPxPyPz(momNegAt);
3284 curElecNegAt.SetXYZM(momNegAt[0],momNegAt[1],momNegAt[2],massE);
3285
3286 Double_t momPosAt[3];
3287 pt.GetPxPyPz(momPosAt);
3288 curElecPosAt.SetXYZM(momPosAt[0],momPosAt[1],momPosAt[2],massE);
3289 if(dca>1){
3290 continue;
3291 }
3292
3293 // Double_t dneg= negTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
3294 // Double_t dpos= posTrack->GetD(xPrimaryVertex,yPrimaryVertex,b);
3295 AliESDv0 vertex(nt,jTracks,pt,iTracks);
3296
3297
3298 Float_t cpa=vertex.GetV0CosineOfPointingAngle(xPrimaryVertex,yPrimaryVertex,zPrimaryVertex);
3299
9c1cb6f7 3300
70ef88b5 3301
9c1cb6f7 3302 // cout<< "v0Rr::"<< v0Rr<<endl;
3303 // if (pvertex.GetRr()<0.5){
3304 // continue;
3305 //}
3306 // cout<<"vertex.GetChi2V0()"<<vertex.GetChi2V0()<<endl;
3307 if(cpa<0.9)continue;
3308 // if (vertex.GetChi2V0() > 30) continue;
3309 // cout<<"xp+xn::"<<xp<<" "<<xn<<endl;
3310 if ((xn+xp) < 0.4) continue;
3311 if (TMath::Abs(ntrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05)
3312 if (TMath::Abs(ptrk->GetD(xPrimaryVertex,yPrimaryVertex,b))<0.05) continue;
3313
3314 //cout<<"pass"<<endl;
3315
3316 AliKFParticle v0GammaC;
3317 v0GammaC+=(*negKF);
3318 v0GammaC+=(*posKF);
3319 v0GammaC.SetMassConstraint(0,0.001);
3320 primVtxImprovedGamma+=v0GammaC;
3321 v0GammaC.SetProductionVertex(primVtxImprovedGamma);
3322
3323
3324 curGamma=curElecNeg+curElecPos;
3325 curGammaAt=curElecNegAt+curElecPosAt;
3326
3327 // invariant mass versus pt of K0short
3328
3329 Double_t chi2V0GammaC=100000.;
3330 if( v0GammaC.GetNDF() != 0) {
3331 chi2V0GammaC = v0GammaC.GetChi2()/v0GammaC.GetNDF();
3332 }else{
3333 cout<< "ERROR::v0K0C.GetNDF()" << endl;
3334 }
3335
3336 if(chi2V0GammaC<200 &&chi2V0GammaC>0 ){
3337 if(fHistograms != NULL){
3338 fHistograms->FillHistogram("ESD_RecalculateV0_InvMass",v0GammaC.GetMass());
3339 fHistograms->FillHistogram("ESD_RecalculateV0_Pt",v0GammaC.GetPt());
3340 fHistograms->FillHistogram("ESD_RecalculateV0_E_dEdxP",curElecNegAt.P(),negTrack->GetTPCsignal());
3341 fHistograms->FillHistogram("ESD_RecalculateV0_P_dEdxP",curElecPosAt.P(),posTrack->GetTPCsignal());
3342 fHistograms->FillHistogram("ESD_RecalculateV0_cpa",cpa);
3343 fHistograms->FillHistogram("ESD_RecalculateV0_dca",dca);
9c1cb6f7 3344 fHistograms->FillHistogram("ESD_RecalculateV0_normdistP",normdistP);
3345 fHistograms->FillHistogram("ESD_RecalculateV0_normdistN",normdistN);
3346
3347 new((*fKFRecalculatedGammasTClone)[fKFRecalculatedGammasTClone->GetEntriesFast()]) AliKFParticle(v0GammaC);
3348 fElectronRecalculatedv1.push_back(iTracks);
3349 fElectronRecalculatedv2.push_back(jTracks);
3350 }
3351 }
3352 }
3353 }
3354
3355 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();firstGammaIndex++){
3356 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFRecalculatedGammasTClone->GetEntriesFast();secondGammaIndex++){
3357 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(firstGammaIndex);
3358 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFRecalculatedGammasTClone->At(secondGammaIndex);
3359
3360 if(fElectronRecalculatedv1[firstGammaIndex]==fElectronRecalculatedv1[secondGammaIndex]){
3361 continue;
3362 }
3363 if( fElectronRecalculatedv2[firstGammaIndex]==fElectronRecalculatedv2[secondGammaIndex]){
3364 continue;
3365 }
3366
3367 AliKFParticle twoGammaCandidate(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
3368 fHistograms->FillHistogram("ESD_RecalculateGG_InvMass",twoGammaCandidate.GetMass());
3369 fHistograms->FillHistogram("ESD_RecalculateGG_InvMass_vs_Pt",twoGammaCandidate.GetMass(),twoGammaCandidate.GetPt());
3370 }
3371 }
3372}
01b7fdcc 3373void AliAnalysisTaskGammaConversion::CalculateJetCone(Int_t gammaIndex){
6c84d371 3374 // CaculateJetCone
a0b94e5c 3375
d7d7e825 3376 Double_t cone;
3377 Double_t coneSize=0.3;
3378 Double_t ptJet=0;
a0b94e5c 3379
6c84d371 3380 // AliKFParticle * currentGamma = &fKFReconstructedGammas[gammaIndex];
3381 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
98778c17 3382
01b7fdcc 3383 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 3384
6c84d371 3385 AliESDtrack* leadingCharged = (AliESDtrack*)(fChargedParticles->At(fLeadingChargedIndex));
98778c17 3386
d7d7e825 3387 Double_t momLeadingCharged[3];
3388 leadingCharged->GetConstrainedPxPyPz(momLeadingCharged);
a0b94e5c 3389
01b7fdcc 3390 TVector3 momentumVectorLeadingCharged(momLeadingCharged[0],momLeadingCharged[1],momLeadingCharged[2]);
a0b94e5c 3391
01b7fdcc 3392 Double_t phi1=momentumVectorLeadingCharged.Phi();
3393 Double_t eta1=momentumVectorLeadingCharged.Eta();
3394 Double_t phi3=momentumVectorCurrentGamma.Phi();
a0b94e5c 3395
6c84d371 3396 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
3397 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 3398 Int_t chId = fChargedParticlesId[iCh];
3399 if(fLeadingChargedIndex==chId || fLeadingChargedIndex==chId) continue;
3400 Double_t mom[3];
3401 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 3402 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
3403 Double_t phi2=momentumVectorChargedParticle.Phi();
3404 Double_t eta2=momentumVectorChargedParticle.Eta();
a0b94e5c 3405
3406
d7d7e825 3407 cone=100.;
3408 if( TMath::Abs(phi2 - phi1) <= ( TMath::TwoPi()-coneSize) ){
3409 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-phi1),2) );
3410 }else{
3411 if( (phi2 - phi1)> TMath::TwoPi()-coneSize ){
3412 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2-TMath::TwoPi()-phi1),2) );
3413 }
3414 if( (phi2 - phi1)< -(TMath::TwoPi()-coneSize) ){
3415 cone = TMath::Sqrt( TMath::Power((eta2-eta1),2)+ TMath::Power((phi2+TMath::TwoPi()-phi1),2) );
3416 }
3417 }
a0b94e5c 3418
01b7fdcc 3419 if(cone <coneSize&& momentumVectorChargedParticle.Pt()>fMinPtJetCone ){
3420 ptJet+= momentumVectorChargedParticle.Pt();
3421 Double_t ffzHdrGam = momentumVectorChargedParticle.Pt()/momentumVectorCurrentGamma.Pt();
3422 Double_t imbalanceHdrGam=-momentumVectorChargedParticle.Dot(momentumVectorCurrentGamma)/momentumVectorCurrentGamma.Mag2();
d7d7e825 3423 fHistograms->FillHistogram("ESD_FFzHdrGam",ffzHdrGam);
3424 fHistograms->FillHistogram("ESD_ImbalanceHdrGam",imbalanceHdrGam);
a0b94e5c 3425
d7d7e825 3426 }
a0b94e5c 3427
d7d7e825 3428 Double_t dphiHdrGam=phi3-phi2;
3429 if ( dphiHdrGam < (-TMath::PiOver2())){
3430 dphiHdrGam+=(TMath::TwoPi());
3431 }
a0b94e5c 3432
d7d7e825 3433 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
3434 dphiHdrGam-=(TMath::TwoPi());
3435 }
a0b94e5c 3436
01b7fdcc 3437 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 3438 fHistograms->FillHistogram("ESD_dphiHdrGamIsolated",dphiHdrGam);
3439 }
3440 }//track loop
a0b94e5c 3441
3442
d7d7e825 3443}
3444
3445Double_t AliAnalysisTaskGammaConversion::GetMinimumDistanceToCharge(Int_t indexHighestPtGamma){
3446 // GetMinimumDistanceToCharge
a0b94e5c 3447
d7d7e825 3448 Double_t fIsoMin=100.;
3449 Double_t ptLeadingCharged=-1.;
98778c17 3450
3451 fLeadingChargedIndex=-1;
a0b94e5c 3452
6c84d371 3453 AliKFParticle * gammaHighestPt = (AliKFParticle*)fKFReconstructedGammasTClone->At(indexHighestPtGamma);
01b7fdcc 3454 TVector3 momentumVectorgammaHighestPt(gammaHighestPt->GetPx(),gammaHighestPt->GetPy(),gammaHighestPt->GetPz());
a0b94e5c 3455
01b7fdcc 3456 Double_t phi1=momentumVectorgammaHighestPt.Phi();
3457 Double_t eta1=momentumVectorgammaHighestPt.Eta();
a0b94e5c 3458
6c84d371 3459 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
3460 AliESDtrack* curTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
d7d7e825 3461 Int_t chId = fChargedParticlesId[iCh];
3462 if(fElectronv1[indexHighestPtGamma]==chId || fElectronv2[indexHighestPtGamma]==chId) continue;
3463 Double_t mom[3];
3464 curTrack->GetConstrainedPxPyPz(mom);
01b7fdcc 3465 TVector3 momentumVectorChargedParticle(mom[0],mom[1],mom[2]);
3466 Double_t phi2=momentumVectorChargedParticle.Phi();
3467 Double_t eta2=momentumVectorChargedParticle.Eta();
d7d7e825 3468 Double_t iso=pow( (pow( (eta1-eta2),2)+ pow((phi1-phi2),2)),0.5 );
a0b94e5c 3469
01b7fdcc 3470 if(momentumVectorChargedParticle.Pt()>fMinPtIsoCone ){
d7d7e825 3471 if (iso<fIsoMin){
3472 fIsoMin=iso;
3473 }
3474 }
a0b94e5c 3475
d7d7e825 3476 Double_t dphiHdrGam=phi1-phi2;
3477 if ( dphiHdrGam < (-TMath::PiOver2())){
3478 dphiHdrGam+=(TMath::TwoPi());
3479 }
a0b94e5c 3480
d7d7e825 3481 if ( dphiHdrGam > (3.*TMath::PiOver2()) ){
3482 dphiHdrGam-=(TMath::TwoPi());
3483 }
01b7fdcc 3484 if (momentumVectorChargedParticle.Pt()>fMinPtGamChargedCorr){
d7d7e825 3485 fHistograms->FillHistogram("ESD_dphiHdrGam",dphiHdrGam);
3486 }
a0b94e5c 3487
d7d7e825 3488 if (dphiHdrGam>0.9*TMath::Pi() && dphiHdrGam<1.1*TMath::Pi()){
a0b94e5c 3489 if (momentumVectorChargedParticle.Pt()> ptLeadingCharged && momentumVectorChargedParticle.Pt()>0.1*momentumVectorgammaHighestPt.Pt()){
01b7fdcc 3490 ptLeadingCharged=momentumVectorChargedParticle.Pt();
d7d7e825 3491 fLeadingChargedIndex=iCh;
3492 }
3493 }
a0b94e5c 3494
d7d7e825 3495 }//track loop
3496 fHistograms->FillHistogram("ESD_MinimumIsoDistance",fIsoMin);
3497 return fIsoMin;
a0b94e5c 3498
d7d7e825 3499}
3500
a0b94e5c 3501Int_t AliAnalysisTaskGammaConversion::GetIndexHighestPtGamma(){
3502 //GetIndexHighestPtGamma
3503
d7d7e825 3504 Int_t indexHighestPtGamma=-1;
01b7fdcc 3505 //Double_t
3506 fGammaPtHighest = -100.;
a0b94e5c 3507
6c84d371 3508 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
3509 AliKFParticle * gammaHighestPtCandidate = (AliKFParticle*)fKFReconstructedGammasTClone->At(firstGammaIndex);
a0b94e5c 3510 TVector3 momentumVectorgammaHighestPtCandidate(gammaHighestPtCandidate->GetPx(),gammaHighestPtCandidate->GetPy(),gammaHighestPtCandidate->GetPz());
3511 if (momentumVectorgammaHighestPtCandidate.Pt() > fGammaPtHighest){
3512 fGammaPtHighest=momentumVectorgammaHighestPtCandidate.Pt();
3513 //gammaHighestPt = gammaHighestPtCandidate;
3514 indexHighestPtGamma=firstGammaIndex;
3515 }
d7d7e825 3516 }
a0b94e5c 3517
d7d7e825 3518 return indexHighestPtGamma;
a0b94e5c 3519
d7d7e825 3520}
3521
3522
3523void AliAnalysisTaskGammaConversion::Terminate(Option_t */*option*/)
3524{
3525 // Terminate analysis
3526 //
3527 AliDebug(1,"Do nothing in Terminate");
3528}
3529
3530void AliAnalysisTaskGammaConversion::UserCreateOutputObjects()
3531{
332f1f44 3532
63477d82 3533 if(fKFCreateAOD) {
3534
3535 //AOD
3536 if(!fAODGamma) fAODGamma = new TClonesArray(fOutputAODClassName.Data(), 0);
3537 else fAODGamma->Delete();
3538 fAODGamma->SetName(Form("%s_gamma", fAODBranchName.Data()));
3539
3540 if(!fAODPi0) fAODPi0 = new TClonesArray(fOutputAODClassName.Data(), 0);
3541 else fAODPi0->Delete();
3542 fAODPi0->SetName(Form("%s_Pi0", fAODBranchName.Data()));
3543
3544 if(!fAODOmega) fAODOmega = new TClonesArray(fOutputAODClassName.Data(), 0);
3545 else fAODOmega->Delete();
3546 fAODOmega->SetName(Form("%s_Omega", fAODBranchName.Data()));
3547
3548 //If delta AOD file name set, add in separate file. Else add in standard aod file.
3549 if(GetDeltaAODFileName().Length() > 0) {
3550 AddAODBranch("TClonesArray", &fAODGamma, GetDeltaAODFileName().Data());
3551 AddAODBranch("TClonesArray", &fAODPi0, GetDeltaAODFileName().Data());
3552 AddAODBranch("TClonesArray", &fAODOmega, GetDeltaAODFileName().Data());
3553 AliAnalysisManager::GetAnalysisManager()->RegisterExtraFile(GetDeltaAODFileName().Data());
3554 } else {
3555 AddAODBranch("TClonesArray", &fAODGamma);
3556 AddAODBranch("TClonesArray", &fAODPi0);
3557 AddAODBranch("TClonesArray", &fAODOmega);
3558 }
332f1f44 3559 }
04bf4381 3560
d7d7e825 3561 // Create the output container
3562 if(fOutputContainer != NULL){
3563 delete fOutputContainer;
3564 fOutputContainer = NULL;
3565 }
3566 if(fOutputContainer == NULL){
3567 fOutputContainer = new TList();
b58d3c74 3568 fOutputContainer->SetOwner(kTRUE);
d7d7e825 3569 }
3570
3571 //Adding the histograms to the output container
3572 fHistograms->GetOutputContainer(fOutputContainer);
3573
3574
3575 if(fWriteNtuple){
3576 if(fGammaNtuple == NULL){
3577 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);
3578 }
3579 if(fNeutralMesonNtuple == NULL){
3580 fNeutralMesonNtuple = new TNtuple("NeutralMesonNtuple","NeutralMesonNtuple","test");
3581 }
3582 TList * ntupleTList = new TList();
b58d3c74 3583 ntupleTList->SetOwner(kTRUE);
d7d7e825 3584 ntupleTList->SetName("Ntuple");
3585 ntupleTList->Add((TNtuple*)fGammaNtuple);
3586 fOutputContainer->Add(ntupleTList);
3587 }
3588
3589 fOutputContainer->SetName(GetName());
3590}
3591
77ac6f3e 3592Double_t AliAnalysisTaskGammaConversion::GetMCOpeningAngle(const TParticle* const daughter0, const TParticle* const daughter1) const{
d7d7e825 3593 //helper function
3594 TVector3 v3D0(daughter0->Px(),daughter0->Py(),daughter0->Pz());
3595 TVector3 v3D1(daughter1->Px(),daughter1->Py(),daughter1->Pz());
3596 return v3D0.Angle(v3D1);
3597}
3598
3599void AliAnalysisTaskGammaConversion::CheckV0Efficiency(){
3600 // see header file for documentation
5e55d806 3601
d7d7e825 3602 vector<Int_t> indexOfGammaParticle;
a0b94e5c 3603
d7d7e825 3604 fStack = fV0Reader->GetMCStack();
a0b94e5c 3605
d7d7e825 3606 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
3607 return; // aborts if the primary vertex does not have contributors.
3608 }
a0b94e5c 3609
d7d7e825 3610 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
3611 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
3612 if(particle->GetPdgCode()==22){ //Gamma
3613 if(particle->GetNDaughters() >= 2){
3614 TParticle* electron=NULL;
3615 TParticle* positron=NULL;
3616 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
3617 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
3618 if(tmpDaughter->GetUniqueID() == 5){
3619 if(tmpDaughter->GetPdgCode() == 11){
3620 electron = tmpDaughter;
3621 }
3622 else if(tmpDaughter->GetPdgCode() == -11){
3623 positron = tmpDaughter;
3624 }
3625 }
3626 }
3627 if(electron!=NULL && positron!=0){
3628 if(electron->R()<160){
3629 indexOfGammaParticle.push_back(iTracks);
3630 }
3631 }
3632 }
3633 }
3634 }
a0b94e5c 3635
d7d7e825 3636 Int_t nFoundGammas=0;
3637 Int_t nNotFoundGammas=0;
a0b94e5c 3638
d7d7e825 3639 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
3640 for(Int_t i=0;i<numberOfV0s;i++){
3641 fV0Reader->GetV0(i);
a0b94e5c 3642
d7d7e825 3643 if(fV0Reader->HasSameMCMother() == kFALSE){
3644 continue;
3645 }
a0b94e5c 3646
d7d7e825 3647 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
3648 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 3649
d7d7e825 3650 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
3651 continue;
3652 }
3653 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
3654 continue;
3655 }
a0b94e5c 3656
d7d7e825 3657 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
3658 //TParticle * v0Gamma = fV0Reader->GetMotherMCParticle();
3659 for(UInt_t mcIndex=0;mcIndex<indexOfGammaParticle.size();mcIndex++){
3660 if(negativeMC->GetFirstMother()==indexOfGammaParticle[mcIndex]){
3661 nFoundGammas++;
3662 }
3663 else{
3664 nNotFoundGammas++;
3665 }
3666 }
3667 }
3668 }
d7d7e825 3669}
3670
3671
3672void AliAnalysisTaskGammaConversion::ProcessGammaElectronsForChicAnalysis(){
3673 // see header file for documantation
a0b94e5c 3674
d7d7e825 3675 fESDEvent = fV0Reader->GetESDEvent();
a0b94e5c 3676
3677
3678 TClonesArray * vESDeNegTemp = new TClonesArray("AliESDtrack",0);
6c84d371 3679 TClonesArray * vESDePosTemp = new TClonesArray("AliESDtrack",0);
3680 TClonesArray * vESDxNegTemp = new TClonesArray("AliESDtrack",0);
3681 TClonesArray * vESDxPosTemp = new TClonesArray("AliESDtrack",0);
3682 TClonesArray * vESDeNegNoJPsi = new TClonesArray("AliESDtrack",0);
3683 TClonesArray * vESDePosNoJPsi = new TClonesArray("AliESDtrack",0);
a0b94e5c 3684
6c84d371 3685 /*
3686 vector <AliESDtrack*> vESDeNegTemp(0);
3687 vector <AliESDtrack*> vESDePosTemp(0);
3688 vector <AliESDtrack*> vESDxNegTemp(0);
3689 vector <AliESDtrack*> vESDxPosTemp(0);
3690 vector <AliESDtrack*> vESDeNegNoJPsi(0);
3691 vector <AliESDtrack*> vESDePosNoJPsi(0);
3692 */
a0b94e5c 3693
3694
d7d7e825 3695 fHistograms->FillTable("Table_Electrons",0);//Count number of Events
a0b94e5c 3696
d7d7e825 3697 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
3698 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 3699
d7d7e825 3700 if(!curTrack){
3701 //print warning here
3702 continue;
3703 }
a0b94e5c 3704
d7d7e825 3705 double p[3];if(!curTrack->GetConstrainedPxPyPz(p))continue;
3706 double r[3];curTrack->GetConstrainedXYZ(r);
a0b94e5c 3707
d7d7e825 3708 TVector3 rXYZ(r);
a0b94e5c 3709
d7d7e825 3710 fHistograms->FillTable("Table_Electrons",4);//Count number of ESD tracks
a0b94e5c 3711
d7d7e825 3712 Bool_t flagKink = kTRUE;
3713 Bool_t flagTPCrefit = kTRUE;
3714 Bool_t flagTRDrefit = kTRUE;
3715 Bool_t flagITSrefit = kTRUE;
3716 Bool_t flagTRDout = kTRUE;
3717 Bool_t flagVertex = kTRUE;
a0b94e5c 3718
3719
d7d7e825 3720 //Cuts ---------------------------------------------------------------
a0b94e5c 3721
d7d7e825 3722 if(curTrack->GetKinkIndex(0) > 0){
3723 fHistograms->FillHistogram("Table_Electrons",5);//Count kink
3724 flagKink = kFALSE;
3725 }
a0b94e5c 3726
d7d7e825 3727 ULong_t trkStatus = curTrack->GetStatus();
a0b94e5c 3728
d7d7e825 3729 ULong_t tpcRefit = (trkStatus & AliESDtrack::kTPCrefit);
a0b94e5c 3730
d7d7e825 3731 if(!tpcRefit){
3732 fHistograms->FillHistogram("Table_Electrons",9);//Count not TPCrefit
3733 flagTPCrefit = kFALSE;
3734 }
a0b94e5c 3735
d7d7e825 3736 ULong_t itsRefit = (trkStatus & AliESDtrack::kITSrefit);
3737 if(!itsRefit){
3738 fHistograms->FillHistogram("Table_Electrons",10);//Count not ITSrefit
3739 flagITSrefit = kFALSE;
3740 }
a0b94e5c 3741
d7d7e825 3742 ULong_t trdRefit = (trkStatus & AliESDtrack::kTRDrefit);
a0b94e5c 3743
d7d7e825 3744 if(!trdRefit){
3745 fHistograms->FillHistogram("Table_Electrons",8); //Count not TRDrefit
3746 flagTRDrefit = kFALSE;
3747 }
a0b94e5c 3748
d7d7e825 3749 ULong_t trdOut = (trkStatus & AliESDtrack::kTRDout);
a0b94e5c 3750
d7d7e825 3751 if(!trdOut) {
3752 fHistograms->FillHistogram("Table_Electrons",7); //Count not TRDout
3753 flagTRDout = kFALSE;
3754 }
a0b94e5c 3755
d7d7e825 3756 double nsigmaToVxt = GetSigmaToVertex(curTrack);
a0b94e5c 3757
d7d7e825 3758 if(nsigmaToVxt > 3){
3759 fHistograms->FillHistogram("Table_Electrons",6); //Count Tracks with number of sigmas > 3
3760 flagVertex = kFALSE;
3761 }
a0b94e5c 3762
d7d7e825 3763 if(! (flagKink && flagTPCrefit && flagITSrefit && flagTRDrefit && flagTRDout && flagVertex ) ) continue;
3764 fHistograms->FillHistogram("Table_Electrons",11);//Count Tracks passed Cuts
a0b94e5c 3765
3766
d7d7e825 3767 Stat_t pid, weight;
3768 GetPID(curTrack, pid, weight);
a0b94e5c 3769
d7d7e825 3770 if(pid!=0){
3771 fHistograms->FillHistogram("Table_Electrons",12); //Count Tracks with pid != 0
3772 }
a0b94e5c 3773
d7d7e825 3774 if(pid == 0){
3775 fHistograms->FillHistogram("Table_Electrons",13); //Count Tracks with pid != 0
3776 }
a0b94e5c 3777
3778
3779
3780
a0b94e5c 3781
3782
d7d7e825 3783 TLorentzVector curElec;
3784 curElec.SetXYZM(p[0],p[1],p[2],fElectronMass);
a0b94e5c 3785
3786
113d8432 3787 if(fDoMCTruth){
3788 Int_t labelMC = TMath::Abs(curTrack->GetLabel());
3789 TParticle* curParticle = fStack->Particle(labelMC);
3790 if(curTrack->GetSign() > 0){
3791 if( pid == 0){
3792 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
3793 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
3794 }
3795 else{
3796 fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
3797 fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
3798 }
3799 }
3800 }
a0b94e5c 3801
3802
d7d7e825 3803 if(curTrack->GetSign() > 0){
a0b94e5c 3804
6c84d371 3805 // vESDxPosTemp.push_back(curTrack);
3806 new((*vESDxPosTemp)[vESDxPosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 3807
d7d7e825 3808 if( pid == 0){
a0b94e5c 3809
d7d7e825 3810 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
3811 fHistograms->FillHistogram("ESD_ElectronPosPt",curElec.Pt());
113d8432 3812 // fHistograms->FillHistogram("MC_ElectronPosNegPt",curParticle->Pt());
d7d7e825 3813 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
113d8432 3814 // fHistograms->FillHistogram("MC_ElectronPosNegEta",curParticle->Eta());
6c84d371 3815 // vESDePosTemp.push_back(curTrack);
3816 new((*vESDePosTemp)[vESDePosTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 3817
d7d7e825 3818 }
a0b94e5c 3819
d7d7e825 3820 }
3821 else {
5e55d806 3822
6c84d371 3823 new((*vESDxNegTemp)[vESDxNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 3824
d7d7e825 3825 if( pid == 0){
a0b94e5c 3826
d7d7e825 3827 fHistograms->FillHistogram("ESD_ElectronPosNegPt",curElec.Pt());
3828 fHistograms->FillHistogram("ESD_ElectronNegPt",curElec.Pt());
d7d7e825 3829 fHistograms->FillHistogram("ESD_ElectronPosNegEta",curElec.Eta());
6c84d371 3830 new((*vESDeNegTemp)[vESDeNegTemp->GetEntriesFast()]) AliESDtrack(*curTrack);
a0b94e5c 3831
d7d7e825 3832 }
a0b94e5c 3833
d7d7e825 3834 }
a0b94e5c 3835
d7d7e825 3836 }
a0b94e5c 3837
3838
d7d7e825 3839 Bool_t ePosJPsi = kFALSE;
3840 Bool_t eNegJPsi = kFALSE;
3841 Bool_t ePosPi0 = kFALSE;
3842 Bool_t eNegPi0 = kFALSE;
3843
3844 UInt_t iePosJPsi=0,ieNegJPsi=0,iePosPi0=0,ieNegPi0=0;
a0b94e5c 3845
6c84d371 3846 for(Int_t iNeg=0; iNeg < vESDeNegTemp->GetEntriesFast(); iNeg++){
3847 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetPdgCode() == 11)
3848 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0) > -1){
3849 Int_t labelMother = fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(iNeg)))->GetLabel()))->GetMother(0);
d7d7e825 3850 TParticle* partMother = fStack ->Particle(labelMother);
3851 if (partMother->GetPdgCode() == 111){
3852 ieNegPi0 = iNeg;
3853 eNegPi0 = kTRUE;
3854 }
3855 if(partMother->GetPdgCode() == 443){ //Mother JPsi
3856 fHistograms->FillTable("Table_Electrons",14);
3857 ieNegJPsi = iNeg;
3858 eNegJPsi = kTRUE;
3859 }
3860 else{
6c84d371 3861 // vESDeNegNoJPsi.push_back(vESDeNegTemp[iNeg]);
3862 new((*vESDeNegNoJPsi)[vESDeNegNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDeNegTemp->At(iNeg)));
d7d7e825 3863 // cout<<"ESD No Positivo JPsi "<<endl;
3864 }
a0b94e5c 3865
d7d7e825 3866 }
3867 }
a0b94e5c 3868
6c84d371 3869 for(Int_t iPos=0; iPos < vESDePosTemp->GetEntriesFast(); iPos++){
3870 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetPdgCode() == -11)
3871 if(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0) > -1){
3872 Int_t labelMother = fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iPos)))->GetLabel()))->GetMother(0);
d7d7e825 3873 TParticle* partMother = fStack ->Particle(labelMother);
3874 if (partMother->GetPdgCode() == 111){
3875 iePosPi0 = iPos;
3876 ePosPi0 = kTRUE;
3877 }
3878 if(partMother->GetPdgCode() == 443){ //Mother JPsi
3879 fHistograms->FillTable("Table_Electrons",15);
3880 iePosJPsi = iPos;
3881 ePosJPsi = kTRUE;
3882 }
3883 else{
6c84d371 3884 // vESDePosNoJPsi.push_back(vESDePosTemp[iPos]);
3885 new((*vESDePosNoJPsi)[vESDePosNoJPsi->GetEntriesFast()]) AliESDtrack(*(AliESDtrack*)(vESDePosTemp->At(iPos)));
d7d7e825 3886 // cout<<"ESD No Negativo JPsi "<<endl;
3887 }
a0b94e5c 3888
d7d7e825 3889 }
3890 }
3891
3892 if( eNegJPsi && ePosJPsi ){
3893 TVector3 tempeNegV,tempePosV;
6c84d371 3894 tempeNegV.SetXYZ(((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Px(),((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Py(),((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->Pz());
3895 tempePosV.SetXYZ(((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Px(),((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Py(),((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->Pz());
d7d7e825 3896 fHistograms->FillTable("Table_Electrons",16);
3897 fHistograms->FillHistogram("ESD_ElectronPosNegJPsiAngle",tempeNegV.Angle(tempePosV));
6c84d371 3898 fHistograms->FillHistogram("MC_ElectronPosNegJPsiAngle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegJPsi)))->GetLabel())),
3899 fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosJPsi)))->GetLabel()))));
d7d7e825 3900 }
3901
3902 if( eNegPi0 && ePosPi0 ){
3903 TVector3 tempeNegV,tempePosV;
6c84d371 3904 tempeNegV.SetXYZ(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Px(),((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Py(),((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->Pz());
3905 tempePosV.SetXYZ(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Px(),((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Py(),((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->Pz());
d7d7e825 3906 fHistograms->FillHistogram("ESD_ElectronPosNegPi0Angle",tempeNegV.Angle(tempePosV));
6c84d371 3907 fHistograms->FillHistogram("MC_ElectronPosNegPi0Angle",GetMCOpeningAngle(fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDeNegTemp->At(ieNegPi0)))->GetLabel())),
3908 fStack->Particle(TMath::Abs(((AliESDtrack*)(vESDePosTemp->At(iePosPi0)))->GetLabel()))));
d7d7e825 3909 }
a0b94e5c 3910
3911
d7d7e825 3912 FillAngle("ESD_eNegePosAngleBeforeCut",GetTLorentzVector(vESDeNegTemp),GetTLorentzVector(vESDePosTemp));
a0b94e5c 3913
6c84d371 3914 CleanWithAngleCuts(*vESDeNegTemp,*vESDePosTemp,*fKFReconstructedGammasTClone);
d7d7e825 3915
6c84d371 3916 // vector <TLorentzVector> vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectron);
3917 // vector <TLorentzVector> vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectron);
a0b94e5c 3918
6c84d371 3919 TClonesArray vCurrentTLVeNeg = GetTLorentzVector(fCurrentEventNegElectronTClone);
3920 TClonesArray vCurrentTLVePos = GetTLorentzVector(fCurrentEventPosElectronTClone);
a0b94e5c 3921
3922
d7d7e825 3923 FillAngle("ESD_eNegePosAngleAfterCut",vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 3924
3925
3926
3927
d7d7e825 3928 //FillAngle("ESD_eNegePosAngleAfterCut",CurrentTLVeNeg,CurrentTLVePos);
a0b94e5c 3929
3930
d7d7e825 3931 FillElectronInvMass("ESD_InvMass_ePluseMinus",vCurrentTLVeNeg,vCurrentTLVePos);
3932 FillElectronInvMass("ESD_InvMass_xPlusxMinus",GetTLorentzVector(vESDxNegTemp),GetTLorentzVector(vESDxPosTemp));
a0b94e5c 3933
3934
3935
6c84d371 3936 FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusChiC","ESD_InvMass_GammaePluseMinusChiCDiff",*fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 3937
d7d7e825 3938 FillGammaElectronInvMass("ESD_InvMass_GammaePluseMinusPi0","ESD_InvMass_GammaePluseMinusPi0Diff",
6c84d371 3939 *fKFReconstructedGammasCutTClone,vCurrentTLVeNeg,vCurrentTLVePos);
a0b94e5c 3940
d7d7e825 3941 //BackGround
a0b94e5c 3942
d7d7e825 3943 //Like Sign e+e-
3944 ElectronBackground("ESD_ENegBackground",vCurrentTLVeNeg);
3945 ElectronBackground("ESD_EPosBackground",vCurrentTLVePos);
3946 ElectronBackground("ESD_EPosENegBackground",vCurrentTLVeNeg);
3947 ElectronBackground("ESD_EPosENegBackground",vCurrentTLVePos);
a0b94e5c 3948
d7d7e825 3949 // Like Sign e+e- no JPsi
3950 ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDeNegNoJPsi));
3951 ElectronBackground("ESD_EPosENegNoJPsiBG",GetTLorentzVector(vESDePosNoJPsi));
a0b94e5c 3952
d7d7e825 3953 //Mixed Event
a0b94e5c 3954
6c84d371 3955 if( fCurrentEventPosElectronTClone->GetEntriesFast() > 0 && fCurrentEventNegElectronTClone->GetEntriesFast() > 0 && fKFReconstructedGammasCutTClone->GetEntriesFast() > 0 ){
d7d7e825 3956 FillGammaElectronInvMass("ESD_EPosENegGammaBackgroundMX","ESD_EPosENegGammaBackgroundMXDiff",
6c84d371 3957 *fKFReconstructedGammasCutTClone,*fPreviousEventTLVNegElectronTClone,*fPreviousEventTLVPosElectronTClone);
3958 *fPreviousEventTLVNegElectronTClone = vCurrentTLVeNeg;
3959 *fPreviousEventTLVPosElectronTClone = vCurrentTLVePos;
a0b94e5c 3960
d7d7e825 3961 }
a0b94e5c 3962
d7d7e825 3963 /*
3964 //Photons P
3965 Double_t vtx[3];
3966 vtx[0]=0;vtx[1]=0;vtx[2]=0;
3967 for(UInt_t i=0;i<fKFReconstructedGammasChic.size();i++){
a0b94e5c 3968
d7d7e825 3969 // if(fMCGammaChicTempCut[i]->GetMother(0) < 0) continue;
a0b94e5c 3970
3971
3972
d7d7e825 3973 Int_t tempLabel = fStack->Particle(fMCGammaChicTempCut[i]->GetMother(0))->GetPdgCode();
3974 // cout<<" Label Pedro Gonzalez " <<tempLabel <<endl;
a0b94e5c 3975
d7d7e825 3976 // cout<<" Label Distance"<<fKFReconstructedGammasChic[i].GetDistanceFromVertex(vtx)<<endl;
a0b94e5c 3977
d7d7e825 3978 if( tempLabel == 10441 || tempLabel == 20443 || tempLabel == 445 )
a0b94e5c 3979
d7d7e825 3980 fHistograms->FillHistogram("ESD_PhotonsMomentum",fKFReconstructedGammasChic[i].GetMomentum());
a0b94e5c 3981
3982
d7d7e825 3983 }
a0b94e5c 3984
3985
d7d7e825 3986 */
1e7846f4 3987
3988
3989 vESDeNegTemp->Delete();
3990 vESDePosTemp->Delete();
3991 vESDxNegTemp->Delete();
3992 vESDxPosTemp->Delete();
3993 vESDeNegNoJPsi->Delete();
3994 vESDePosNoJPsi->Delete();
3995
3996 delete vESDeNegTemp;
3997 delete vESDePosTemp;
3998 delete vESDxNegTemp;
3999 delete vESDxPosTemp;
4000 delete vESDeNegNoJPsi;
4001 delete vESDePosNoJPsi;
d7d7e825 4002}
4003
6c84d371 4004/*
4005 void AliAnalysisTaskGammaConversion::FillAngle(TString histoName,vector <TLorentzVector> tlVeNeg, vector <TLorentzVector> tlVePos){
d7d7e825 4006 //see header file for documentation
4007 for( UInt_t iNeg=0; iNeg < tlVeNeg.size(); iNeg++){
6c84d371 4008 for (UInt_t iPos=0; iPos < tlVePos.size(); iPos++){
4009 fHistograms->FillHistogram(histoName.Data(),tlVeNeg[iNeg].Vect().Angle(tlVePos[iPos].Vect()));
4010 }
4011 }
4012 }
4013*/
4014void AliAnalysisTaskGammaConversion::FillAngle(TString histoName,TClonesArray const tlVeNeg, TClonesArray const tlVePos){
4015 //see header file for documentation
4016 for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++){
4017 for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){
4018 fHistograms->FillHistogram(histoName.Data(),((TLorentzVector*)(tlVeNeg.At(iNeg)))->Vect().Angle(((TLorentzVector*)(tlVePos.At(iPos)))->Vect()));
d7d7e825 4019 }
4020 }
4021}
6c84d371 4022void AliAnalysisTaskGammaConversion::FillElectronInvMass(TString histoName, TClonesArray const eNeg, TClonesArray const ePos){
d7d7e825 4023 //see header file for documentation
6c84d371 4024 for( Int_t n=0; n < eNeg.GetEntriesFast(); n++){
4025 TLorentzVector en = (*(TLorentzVector*)(eNeg.At(n)));
4026 for (Int_t p=0; p < ePos.GetEntriesFast(); p++){
4027 TLorentzVector ep = (*(TLorentzVector*)(ePos.At(p)));
d7d7e825 4028 TLorentzVector np = ep + en;
4029 fHistograms->FillHistogram(histoName.Data(),np.M());
4030 }
4031 }
d7d7e825 4032}
4033
6c84d371 4034void AliAnalysisTaskGammaConversion::FillGammaElectronInvMass(TString histoMass,TString histoDiff,TClonesArray const fKFGammas,
4035 TClonesArray const tlVeNeg,TClonesArray const tlVePos)
d7d7e825 4036{
4037 //see header file for documentation
a0b94e5c 4038
6c84d371 4039 for( Int_t iNeg=0; iNeg < tlVeNeg.GetEntriesFast(); iNeg++ ){
a0b94e5c 4040
6c84d371 4041 for (Int_t iPos=0; iPos < tlVePos.GetEntriesFast(); iPos++){
a0b94e5c 4042
6c84d371 4043 TLorentzVector xy = *((TLorentzVector *)(tlVePos.At(iPos))) + *((TLorentzVector *)(tlVeNeg.At(iNeg)));
a0b94e5c 4044
6c84d371 4045 for (Int_t iGam=0; iGam < fKFGammas.GetEntriesFast(); iGam++){
a0b94e5c 4046
6c84d371 4047 // AliKFParticle * gammaCandidate = &fKFGammas[iGam];
4048 AliKFParticle * gammaCandidate = (AliKFParticle *)(fKFGammas.At(iGam));
d7d7e825 4049 TLorentzVector g;
a0b94e5c 4050
d7d7e825 4051 g.SetXYZM(gammaCandidate->GetPx(),gammaCandidate->GetPy(),gammaCandidate->GetPz(),fGammaMass);
4052 TLorentzVector xyg = xy + g;
4053 fHistograms->FillHistogram(histoMass.Data(),xyg.M());
4054 fHistograms->FillHistogram(histoDiff.Data(),(xyg.M()-xy.M()));
4055 }
4056 }
4057 }
a0b94e5c 4058
d7d7e825 4059}
6c84d371 4060void AliAnalysisTaskGammaConversion::ElectronBackground(TString hBg, TClonesArray e)
d7d7e825 4061{
4062 // see header file for documentation
6c84d371 4063 for(Int_t i=0; i < e.GetEntriesFast(); i++)
d7d7e825 4064 {
6c84d371 4065 for (Int_t j=i+1; j < e.GetEntriesFast(); j++)
d7d7e825 4066 {
6c84d371 4067 TLorentzVector ee = (*(TLorentzVector*)(e.At(i))) + (*(TLorentzVector*)(e.At(j)));
a0b94e5c 4068
d7d7e825 4069 fHistograms->FillHistogram(hBg.Data(),ee.M());
4070 }
4071 }
4072}
4073
4074
6c84d371 4075void AliAnalysisTaskGammaConversion::CleanWithAngleCuts(TClonesArray const negativeElectrons,
4076 TClonesArray const positiveElectrons,
4077 TClonesArray const gammas){
d7d7e825 4078 // see header file for documentation
a0b94e5c 4079
6c84d371 4080 UInt_t sizeN = negativeElectrons.GetEntriesFast();
4081 UInt_t sizeP = positiveElectrons.GetEntriesFast();
4082 UInt_t sizeG = gammas.GetEntriesFast();
a0b94e5c 4083
4084
4085
d7d7e825 4086 vector <Bool_t> xNegBand(sizeN);
4087 vector <Bool_t> xPosBand(sizeP);
4088 vector <Bool_t> gammaBand(sizeG);
a0b94e5c 4089
4090
d7d7e825 4091 for(UInt_t iNeg=0; iNeg < sizeN; iNeg++) xNegBand[iNeg]=kTRUE;
4092 for(UInt_t iPos=0; iPos < sizeP; iPos++) xPosBand[iPos]=kTRUE;
4093 for(UInt_t iGam=0; iGam < sizeG; iGam++) gammaBand[iGam]=kTRUE;
d7d7e825 4094
a0b94e5c 4095
4096 for(UInt_t iPos=0; iPos < sizeP; iPos++){
4097
6c84d371 4098 Double_t aP[3];
4099 ((AliESDtrack*)(positiveElectrons.At(iPos)))->GetConstrainedPxPyPz(aP);
a0b94e5c 4100
d7d7e825 4101 TVector3 ePosV(aP[0],aP[1],aP[2]);
a0b94e5c 4102
d7d7e825 4103 for(UInt_t iNeg=0; iNeg < sizeN; iNeg++){
a0b94e5c 4104
6c84d371 4105 Double_t aN[3];
4106 ((AliESDtrack*)(negativeElectrons.At(iNeg)))->GetConstrainedPxPyPz(aN);
d7d7e825 4107 TVector3 eNegV(aN[0],aN[1],aN[2]);
a0b94e5c 4108
d7d7e825 4109 if(ePosV.Angle(eNegV) < 0.05){ //e+e- from gamma
4110 xPosBand[iPos]=kFALSE;
4111 xNegBand[iNeg]=kFALSE;
4112 }
a0b94e5c 4113
d7d7e825 4114 for(UInt_t iGam=0; iGam < sizeG; iGam++){
6c84d371 4115 AliKFParticle* gammaCandidate = (AliKFParticle*)gammas.At(iGam);
d7d7e825 4116 TVector3 gammaCandidateVector(gammaCandidate->Px(),gammaCandidate->Py(),gammaCandidate->Pz());
4117 if(ePosV.Angle(gammaCandidateVector) < 0.05 || eNegV.Angle(gammaCandidateVector) < 0.05)
4118 gammaBand[iGam]=kFALSE;
4119 }
4120 }
4121 }
a0b94e5c 4122
4123
4124
4125
d7d7e825 4126 for(UInt_t iPos=0; iPos < sizeP; iPos++){
4127 if(xPosBand[iPos]){
6c84d371 4128 new((*fCurrentEventPosElectronTClone)[fCurrentEventPosElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(positiveElectrons.At(iPos))));
4129 // fCurrentEventPosElectron.push_back(positiveElectrons[iPos]);
d7d7e825 4130 }
4131 }
4132 for(UInt_t iNeg=0;iNeg < sizeN; iNeg++){
4133 if(xNegBand[iNeg]){
6c84d371 4134 new((*fCurrentEventNegElectronTClone)[fCurrentEventNegElectronTClone->GetEntriesFast()]) AliESDtrack((*(AliESDtrack*)(negativeElectrons.At(iNeg))));
4135 // fCurrentEventNegElectron.push_back(negativeElectrons[iNeg]);
d7d7e825 4136 }
4137 }
4138 for(UInt_t iGam=0; iGam < sizeG; iGam++){
4139 if(gammaBand[iGam]){
6c84d371 4140 new((*fKFReconstructedGammasCutTClone)[fKFReconstructedGammasCutTClone->GetEntriesFast()]) AliKFParticle((*(AliKFParticle*)(gammas.At(iGam))));
4141 //fKFReconstructedGammasCut.push_back(*(AliKFParticle*)gammas->At(iGam));
d7d7e825 4142 }
4143 }
4144}
4145
4146
77ac6f3e 4147void AliAnalysisTaskGammaConversion::GetPID(const AliESDtrack *track, Stat_t &pid, Stat_t &weight)
d7d7e825 4148{
4149 // see header file for documentation
4150 pid = -1;
4151 weight = -1;
a0b94e5c 4152
d7d7e825 4153 double wpart[5];
4154 double wpartbayes[5];
a0b94e5c 4155
d7d7e825 4156 //get probability of the diffenrent particle types
4157 track->GetESDpid(wpart);
a0b94e5c 4158
d7d7e825 4159 // Tentative particle type "concentrations"
4160 double c[5]={0.01, 0.01, 0.85, 0.10, 0.05};
a0b94e5c 4161
d7d7e825 4162 //Bayes' formula
4163 double rcc = 0.;
4164 for (int i = 0; i < 5; i++)
4165 {
4166 rcc+=(c[i] * wpart[i]);
4167 }
a0b94e5c 4168
4169
4170
d7d7e825 4171 for (int i=0; i<5; i++) {
4a6157dc 4172 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 4173 wpartbayes[i] = c[i] * wpart[i] / rcc;
4174 }
4175 }
a0b94e5c 4176
4177
4178
d7d7e825 4179 Float_t max=0.;
4180 int ipid=-1;
4181 //find most probable particle in ESD pid
4182 //0:Electron - 1:Muon - 2:Pion - 3:Kaon - 4:Proton
4183 for (int i = 0; i < 5; i++)
4184 {
4185 if (wpartbayes[i] > max)
4186 {
a0b94e5c 4187 ipid = i;
4188 max = wpartbayes[i];
d7d7e825 4189 }
4190 }
a0b94e5c 4191
d7d7e825 4192 pid = ipid;
4193 weight = max;
4194}
77ac6f3e 4195double AliAnalysisTaskGammaConversion::GetSigmaToVertex(const AliESDtrack* t)
d7d7e825 4196{
4197 // Calculates the number of sigma to the vertex.
a0b94e5c 4198
d7d7e825 4199 Float_t b[2];
4200 Float_t bRes[2];
4201 Float_t bCov[3];
4202 t->GetImpactParameters(b,bCov);
4203 if (bCov[0]<=0 || bCov[2]<=0) {
4204 AliDebug(1, "Estimated b resolution lower or equal zero!");
4205 bCov[0]=0; bCov[2]=0;
4206 }
4207 bRes[0] = TMath::Sqrt(bCov[0]);
4208 bRes[1] = TMath::Sqrt(bCov[2]);
a0b94e5c 4209
d7d7e825 4210 // -----------------------------------
4211 // How to get to a n-sigma cut?
4212 //
4213 // The accumulated statistics from 0 to d is
4214 //
4215 // -> Erf(d/Sqrt(2)) for a 1-dim gauss (d = n_sigma)
4216 // -> 1 - Exp(-d**2) for a 2-dim gauss (d*d = dx*dx + dy*dy != n_sigma)
4217 //
4218 // It means that for a 2-dim gauss: n_sigma(d) = Sqrt(2)*ErfInv(1 - Exp((-x**2)/2)
4219 // Can this be expressed in a different way?
a0b94e5c 4220
d7d7e825 4221 if (bRes[0] == 0 || bRes[1] ==0)
4222 return -1;
a0b94e5c 4223
d7d7e825 4224 double d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
a0b94e5c 4225
d7d7e825 4226 // stupid rounding problem screws up everything:
4227 // if d is too big, TMath::Exp(...) gets 0, and TMath::ErfInverse(1) that should be infinite, gets 0 :(
4228 if (TMath::Exp(-d * d / 2) < 1e-10)
4229 return 1000;
a0b94e5c 4230
4231
d7d7e825 4232 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
4233 return d;
4234}
6c84d371 4235
4236//vector <TLorentzVector> AliAnalysisTaskGammaConversion::GetTLorentzVector(vector <AliESDtrack*> esdTrack){
4237TClonesArray AliAnalysisTaskGammaConversion::GetTLorentzVector(TClonesArray *const esdTrack){
01b7fdcc 4238 //Return TLoresntz vector of track?
6c84d371 4239 // vector <TLorentzVector> tlVtrack(0);
4240 TClonesArray array("TLorentzVector",0);
a0b94e5c 4241
6c84d371 4242 for(Int_t itrack=0; itrack < esdTrack->GetEntriesFast(); itrack++){
4243 double p[3];
4244 //esdTrack[itrack]->GetConstrainedPxPyPz(p);
4245 ((AliESDtrack*)(esdTrack->At(itrack)))->GetConstrainedPxPyPz(p);
d7d7e825 4246 TLorentzVector currentTrack;
01b7fdcc 4247 currentTrack.SetXYZM(p[0],p[1],p[2],fElectronMass);
6c84d371 4248 new((array)[array.GetEntriesFast()]) TLorentzVector(currentTrack);
4249 // tlVtrack.push_back(currentTrack);
d7d7e825 4250 }
a0b94e5c 4251
6c84d371 4252 return array;
a0b94e5c 4253
6c84d371 4254 // return tlVtrack;
d7d7e825 4255}
77ac6f3e 4256Int_t AliAnalysisTaskGammaConversion::GetProcessType(const AliMCEvent * mcEvt) {
3c45d101 4257
4258 // Determine if the event was generated with pythia or phojet and return the process type
4259
4260 // Check if mcEvt is fine
4261 if (!mcEvt) {
4262 AliFatal("NULL mc event");
4263 }
4264
4265 // Determine if it was a pythia or phojet header, and return the correct process type
4266 AliGenPythiaEventHeader * headPy = 0;
4267 AliGenDPMjetEventHeader * headPho = 0;
4268 AliGenEventHeader * htmp = mcEvt->GenEventHeader();
4269 if(!htmp) {
4270 AliFatal("Cannot Get MC Header!!");
4271 }
4272 if( TString(htmp->IsA()->GetName()) == "AliGenPythiaEventHeader") {
4273 headPy = (AliGenPythiaEventHeader*) htmp;
4274 } else if (TString(htmp->IsA()->GetName()) == "AliGenDPMjetEventHeader") {
4275 headPho = (AliGenDPMjetEventHeader*) htmp;
4276 } else {
4277 AliError("Unknown header");
4278 }
4279
4280 // Determine process type
4281 if(headPy) {
4282 if(headPy->ProcessType() == 92 || headPy->ProcessType() == 93) {
4283 // single difractive
4284 return kProcSD;
4285 } else if (headPy->ProcessType() == 94) {
4286 // double diffractive
4287 return kProcDD;
4288 }
4289 else if(headPy->ProcessType() != 92 && headPy->ProcessType() != 93 && headPy->ProcessType() != 94) {
4290 // non difractive
4291 return kProcND;
4292 }
4293 } else if (headPho) {
4294 if(headPho->ProcessType() == 5 || headPho->ProcessType() == 6 ) {
4295 // single difractive
4296 return kProcSD;
4297 } else if (headPho->ProcessType() == 7) {
4298 // double diffractive
4299 return kProcDD;
4300 } else if(headPho->ProcessType() != 5 && headPho->ProcessType() != 6 && headPho->ProcessType() != 7 ) {
4301 // non difractive
4302 return kProcND;
4303 }
4304 }
4305
4306
4307 // no process type found?
4308 AliError(Form("Unknown header: %s", htmp->IsA()->GetName()));
4309 return kProcUnknown;
4310}
6746e1e1 4311
4312
4313Int_t AliAnalysisTaskGammaConversion::CalculateMultiplicityBin(){
4314 // Get Centrality bin
4315
4316 Int_t multiplicity = 0;
4317
4318 if ( fUseMultiplicity == 1 ) {
4319
4320 if (fMultiplicity>= 0 && fMultiplicity<= 5) multiplicity=1;
4321 if (fMultiplicity>= 6 && fMultiplicity<= 9) multiplicity=2;
4322 if (fMultiplicity>=10 && fMultiplicity<=14) multiplicity=3;
4323 if (fMultiplicity>=15 && fMultiplicity<=22) multiplicity=4;
4324 if (fMultiplicity>=23 ) multiplicity=5;
4325
4326 }
4327 return multiplicity;
4328}