]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/GammaConv/AliAnalysisTaskGammaConversion.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
New macro to plot SDD offline calibration quantities
[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"
4a6157dc 35#include "AliGammaConversionAODObject.h"
5e55d806 36#include "AliGammaConversionBGHandler.h"
6272370b 37#include "AliESDCaloCluster.h" // for combining PHOS and GammaConv
9c1cb6f7 38#include "AliKFVertex.h"
3c45d101 39#include "AliGenPythiaEventHeader.h"
40#include "AliGenDPMjetEventHeader.h"
41#include "AliGenEventHeader.h"
e60f3265 42#include <AliMCEventHandler.h>
d707e3cf 43class AliESDTrackCuts;
4a6157dc 44class AliCFContainer;
45class AliCFManager;
d7d7e825 46class AliKFVertex;
47class AliAODHandler;
48class AliAODEvent;
49class ALiESDEvent;
50class AliMCEvent;
51class AliMCEventHandler;
52class AliESDInputHandler;
53class AliAnalysisManager;
54class Riostream;
55class TFile;
56class TInterpreter;
57class TSystem;
58class TROOT;
59
60ClassImp(AliAnalysisTaskGammaConversion)
61
62
63AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion():
64AliAnalysisTaskSE(),
65 fV0Reader(NULL),
66 fStack(NULL),
a0b94e5c 67 fMCTruth(NULL), // for CF
4a6157dc 68 fGCMCEvent(NULL), // for CF
d7d7e825 69 fESDEvent(NULL),
70 fOutputContainer(NULL),
a0b94e5c 71 fCFManager(0x0), // for CF
d7d7e825 72 fHistograms(NULL),
b5832f95 73 fTriggerCINT1B(kFALSE),
d7d7e825 74 fDoMCTruth(kFALSE),
75 fDoNeutralMeson(kFALSE),
6272370b 76 fDoOmegaMeson(kFALSE),
d7d7e825 77 fDoJet(kFALSE),
78 fDoChic(kFALSE),
9c1cb6f7 79 fRecalculateV0ForGamma(kFALSE),
6c84d371 80 fKFReconstructedGammasTClone(NULL),
6272370b 81 fKFReconstructedPi0sTClone(NULL),
9c1cb6f7 82 fKFRecalculatedGammasTClone(NULL),
6c84d371 83 fCurrentEventPosElectronTClone(NULL),
84 fCurrentEventNegElectronTClone(NULL),
85 fKFReconstructedGammasCutTClone(NULL),
86 fPreviousEventTLVNegElectronTClone(NULL),
87 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 88 fElectronv1(),
89 fElectronv2(),
6272370b 90 fGammav1(),
91 fGammav2(),
9c1cb6f7 92 fElectronRecalculatedv1(),
93 fElectronRecalculatedv2(),
d7d7e825 94 fElectronMass(-1),
95 fGammaMass(-1),
96 fPi0Mass(-1),
97 fEtaMass(-1),
98 fGammaWidth(-1),
99 fPi0Width(-1),
100 fEtaWidth(-1),
101 fMinOpeningAngleGhostCut(0.),
9640a3d1 102 fEsdTrackCuts(NULL),
d7d7e825 103 fCalculateBackground(kFALSE),
104 fWriteNtuple(kFALSE),
105 fGammaNtuple(NULL),
106 fNeutralMesonNtuple(NULL),
107 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 108 fChargedParticles(NULL),
d7d7e825 109 fChargedParticlesId(),
110 fGammaPtHighest(0.),
111 fMinPtForGammaJet(1.),
112 fMinIsoConeSize(0.2),
113 fMinPtIsoCone(0.7),
114 fMinPtGamChargedCorr(0.5),
115 fMinPtJetCone(0.5),
116 fLeadingChargedIndex(-1),
9640a3d1 117 fLowPtMapping(1.),
118 fHighPtMapping(3.),
1e7846f4 119 fDoCF(kFALSE),
04bf4381 120 fAODGamma(NULL),
121 fAODPi0(NULL),
122 fAODOmega(NULL),
037dc2db 123 fAODBranchName("GammaConv"),
d765d400 124 fKFForceAOD(kFALSE),
332f1f44 125 fKFDeltaAODFileName(""),
037dc2db 126 fDoNeutralMesonV0MCCheck(kFALSE),
5ce758b0 127 fUseTrackMultiplicityForBG(kTRUE),
128 fMoveParticleAccordingToVertex(kFALSE),
2e2da371 129 fApplyChi2Cut(kFALSE),
111d75df 130 nRandomEventsForBG(15),
131 nDegreesPMBackground(15),
132 fDoRotation(kTRUE),
133 fCheckBGProbability(kTRUE),
037dc2db 134 fKFReconstructedGammasV0Index()
d7d7e825 135{
136 // Default constructor
5a34881d 137
138 /* Kenneth: the default constructor should not have any define input/output or the call to SetESDtrackCuts
d7d7e825 139 // Common I/O in slot 0
140 DefineInput (0, TChain::Class());
141 DefineOutput(0, TTree::Class());
142
143 // Your private output
144 DefineOutput(1, TList::Class());
a0b94e5c 145
d7d7e825 146 // Define standard ESD track cuts for Gamma-hadron correlation
147 SetESDtrackCuts();
5a34881d 148 */
d7d7e825 149}
150
151AliAnalysisTaskGammaConversion::AliAnalysisTaskGammaConversion(const char* name):
152 AliAnalysisTaskSE(name),
153 fV0Reader(NULL),
154 fStack(NULL),
a0b94e5c 155 fMCTruth(NULL), // for CF
4a6157dc 156 fGCMCEvent(NULL), // for CF
d7d7e825 157 fESDEvent(NULL),
158 fOutputContainer(0x0),
a0b94e5c 159 fCFManager(0x0), // for CF
d7d7e825 160 fHistograms(NULL),
b5832f95 161 fTriggerCINT1B(kFALSE),
d7d7e825 162 fDoMCTruth(kFALSE),
163 fDoNeutralMeson(kFALSE),
6272370b 164 fDoOmegaMeson(kFALSE),
d7d7e825 165 fDoJet(kFALSE),
166 fDoChic(kFALSE),
9c1cb6f7 167 fRecalculateV0ForGamma(kFALSE),
6c84d371 168 fKFReconstructedGammasTClone(NULL),
6272370b 169 fKFReconstructedPi0sTClone(NULL),
9c1cb6f7 170 fKFRecalculatedGammasTClone(NULL),
6c84d371 171 fCurrentEventPosElectronTClone(NULL),
172 fCurrentEventNegElectronTClone(NULL),
173 fKFReconstructedGammasCutTClone(NULL),
174 fPreviousEventTLVNegElectronTClone(NULL),
175 fPreviousEventTLVPosElectronTClone(NULL),
d7d7e825 176 fElectronv1(),
177 fElectronv2(),
6272370b 178 fGammav1(),
179 fGammav2(),
9c1cb6f7 180 fElectronRecalculatedv1(),
181 fElectronRecalculatedv2(),
d7d7e825 182 fElectronMass(-1),
183 fGammaMass(-1),
184 fPi0Mass(-1),
185 fEtaMass(-1),
186 fGammaWidth(-1),
187 fPi0Width(-1),
188 fEtaWidth(-1),
189 fMinOpeningAngleGhostCut(0.),
9640a3d1 190 fEsdTrackCuts(NULL),
d7d7e825 191 fCalculateBackground(kFALSE),
192 fWriteNtuple(kFALSE),
193 fGammaNtuple(NULL),
194 fNeutralMesonNtuple(NULL),
195 fTotalNumberOfAddedNtupleEntries(0),
6c84d371 196 fChargedParticles(NULL),
d7d7e825 197 fChargedParticlesId(),
198 fGammaPtHighest(0.),
199 fMinPtForGammaJet(1.),
200 fMinIsoConeSize(0.2),
201 fMinPtIsoCone(0.7),
202 fMinPtGamChargedCorr(0.5),
203 fMinPtJetCone(0.5),
204 fLeadingChargedIndex(-1),
9640a3d1 205 fLowPtMapping(1.),
206 fHighPtMapping(3.),
c59360eb 207 fDoCF(kFALSE),
208 fAODGamma(NULL),
209 fAODPi0(NULL),
210 fAODOmega(NULL),
211 fAODBranchName("GammaConv"),
212 fKFForceAOD(kFALSE),
213 fKFDeltaAODFileName(""),
214 fDoNeutralMesonV0MCCheck(kFALSE),
5ce758b0 215 fUseTrackMultiplicityForBG(kTRUE),
216 fMoveParticleAccordingToVertex(kFALSE),
2e2da371 217 fApplyChi2Cut(kFALSE),
111d75df 218 nRandomEventsForBG(15),
219 nDegreesPMBackground(15),
220 fDoRotation(kTRUE),
221 fCheckBGProbability(kTRUE),
037dc2db 222 fKFReconstructedGammasV0Index()
d7d7e825 223{
224 // Common I/O in slot 0
225 DefineInput (0, TChain::Class());
226 DefineOutput(0, TTree::Class());
227
228 // Your private output
229 DefineOutput(1, TList::Class());
a0b94e5c 230 DefineOutput(2, AliCFContainer::Class()); // for CF
231
232
d7d7e825 233 // Define standard ESD track cuts for Gamma-hadron correlation
234 SetESDtrackCuts();
9c1cb6f7 235
d7d7e825 236}
237
238AliAnalysisTaskGammaConversion::~AliAnalysisTaskGammaConversion()
239{
240 // Remove all pointers
241
242 if(fOutputContainer){
243 fOutputContainer->Clear() ;
244 delete fOutputContainer ;
245 }
246 if(fHistograms){
247 delete fHistograms;
248 }
249 if(fV0Reader){
250 delete fV0Reader;
251 }
a0b94e5c 252
253 // for CF
254 if(fCFManager){
255 delete fCFManager;
256 }
9640a3d1 257
258 if(fEsdTrackCuts){
259 delete fEsdTrackCuts;
260 }
261
04bf4381 262 //Delete AODs
263 if (fAODGamma) {
264 fAODGamma->Clear();
265 delete fAODGamma;
d7d7e825 266 }
04bf4381 267 fAODGamma = NULL;
268
269 if (fAODPi0) {
270 fAODPi0->Clear();
271 delete fAODPi0;
272 }
273 fAODPi0 = NULL;
d7d7e825 274
04bf4381 275 if (fAODOmega) {
276 fAODOmega->Clear();
277 delete fAODOmega;
278 }
279 fAODOmega = NULL;
280
04bf4381 281}
d7d7e825 282
12464034 283
d7d7e825 284void AliAnalysisTaskGammaConversion::Init()
285{
286 // Initialization
4a6157dc 287 // AliLog::SetGlobalLogLevel(AliLog::kError);
d7d7e825 288}
289void AliAnalysisTaskGammaConversion::SetESDtrackCuts()
290{
291 // SetESDtrackCuts
9640a3d1 292 if (fEsdTrackCuts!=NULL){
293 delete fEsdTrackCuts;
294 }
d7d7e825 295 fEsdTrackCuts = new AliESDtrackCuts("AliESDtrackCuts");
a0b94e5c 296 //standard cuts from:
297 //http://aliceinfo.cern.ch/alicvs/viewvc/PWG0/dNdEta/CreateCuts.C?revision=1.4&view=markup
ebcfaa7e 298
299 // Cuts used up to 3rd of March
300
301 // fEsdTrackCuts->SetMinNClustersTPC(50);
302// fEsdTrackCuts->SetMaxChi2PerClusterTPC(3.5);
303// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
304// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
305// fEsdTrackCuts->SetMaxNsigmaToVertex(3);
306// fEsdTrackCuts->SetRequireSigmaToVertex(kTRUE);
307
308 //------- To be tested-----------
d707e3cf 309 // Cuts used up to 26th of Agost
310// Int_t minNClustersTPC = 70;
311// Double_t maxChi2PerClusterTPC = 4.0;
312// Double_t maxDCAtoVertexXY = 2.4; // cm
313// Double_t maxDCAtoVertexZ = 3.2; // cm
314// fEsdTrackCuts->SetRequireSigmaToVertex(kFALSE);
315// fEsdTrackCuts->SetRequireTPCRefit(kTRUE);
316// fEsdTrackCuts->SetRequireITSRefit(kTRUE);
317// // fEsdTrackCuts->SetRequireTPCStandAlone(kTRUE);
318// fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
319// fEsdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
320// fEsdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
321// fEsdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
322// fEsdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
323// fEsdTrackCuts->SetDCAToVertex2D(kTRUE);
324// fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
325// fEsdTrackCuts->SetPtRange(0.15);
326
327// fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kFirst);
328
329
330// Using standard function for setting Cuts
331 Bool_t selectPrimaries=kTRUE;
332 fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selectPrimaries);
333 fEsdTrackCuts->SetEtaRange(-0.8, 0.8);
334 fEsdTrackCuts->SetPtRange(0.15);
ebcfaa7e 335
336 //----- From Jacek 10.03.03 ------------------/
337// minNClustersTPC = 70;
338// maxChi2PerClusterTPC = 4.0;
339// maxDCAtoVertexXY = 2.4; // cm
340// maxDCAtoVertexZ = 3.2; // cm
341
342// esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
343// esdTrackCuts->SetRequireTPCRefit(kFALSE);
344// esdTrackCuts->SetRequireTPCStandAlone(kTRUE);
345// esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
346// esdTrackCuts->SetMinNClustersTPC(minNClustersTPC);
347// esdTrackCuts->SetMaxChi2PerClusterTPC(maxChi2PerClusterTPC);
348// esdTrackCuts->SetMaxDCAToVertexXY(maxDCAtoVertexXY);
349// esdTrackCuts->SetMaxDCAToVertexZ(maxDCAtoVertexZ);
350// esdTrackCuts->SetDCAToVertex2D(kTRUE);
351
352
353
d7d7e825 354 // fEsdTrackCuts->SetAcceptKinkDaughters(kFALSE);
037dc2db 355 // fV0Reader->SetESDtrackCuts(fEsdTrackCuts);
d7d7e825 356}
357
c00009fb 358void AliAnalysisTaskGammaConversion::UserExec(Option_t */*option*/)
d7d7e825 359{
360 // Execute analysis for current event
6272370b 361
9c1cb6f7 362 // Load the esdpid from the esdhandler if exists (tender was applied) otherwise set the Bethe Bloch parameters
3c45d101 363 Int_t eventQuality=-1;
9c1cb6f7 364
365 AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager();
366 AliESDInputHandler *esdHandler=0x0;
367 if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){
368 AliV0Reader::SetESDpid(esdHandler->GetESDpid());
369 } else {
370 //load esd pid bethe bloch parameters depending on the existance of the MC handler
371 // yes: MC parameters
372 // no: data parameters
373 if (!AliV0Reader::GetESDpid()){
374 if (fMCEvent ) {
375 AliV0Reader::InitESDpid();
376 } else {
377 AliV0Reader::InitESDpid(1);
378 }
379 }
380 }
381
d765d400 382 //Must set fForceAOD to true for the AOD to get filled. Should only be done when running independent chain / train.
383 if(fKFForceAOD) {
384 if (!AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()) AliFatal("Cannot run ESD filter without an output event handler");
385 AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler()->SetFillAOD(kTRUE);
386 }
48682642 387
388 if(fV0Reader == NULL){
389 // Write warning here cuts and so on are default if this ever happens
6272370b 390 }
48682642 391
e60f3265 392 if (fMCEvent ) {
393 // To avoid crashes due to unzip errors. Sometimes the trees are not there.
394
395 AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
3c45d101 396 if (!mcHandler){
397 AliError("Could not retrive MC event handler!");
398 return;
399
400 eventQuality=0;
401 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
402 }
403 if (!mcHandler->InitOk() ){
404 eventQuality=0;
405 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
406 return;
407 }
408 if (!mcHandler->TreeK() ){
409 eventQuality=0;
410 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
411 return;
412 }
413 if (!mcHandler->TreeTR() ) {
414 eventQuality=0;
415 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
416 return;
417 }
e60f3265 418 }
419
48682642 420 fV0Reader->SetInputAndMCEvent(InputEvent(), MCEvent());
421
422 fV0Reader->Initialize();
423 fDoMCTruth = fV0Reader->GetDoMCTruth();
424
12464034 425 if(fAODGamma) fAODGamma->Delete();
426 if(fAODPi0) fAODPi0->Delete();
427 if(fAODOmega) fAODOmega->Delete();
a0b94e5c 428
6c84d371 429 if(fKFReconstructedGammasTClone == NULL){
430 fKFReconstructedGammasTClone = new TClonesArray("AliKFParticle",0);
431 }
432 if(fCurrentEventPosElectronTClone == NULL){
433 fCurrentEventPosElectronTClone = new TClonesArray("AliESDtrack",0);
434 }
435 if(fCurrentEventNegElectronTClone == NULL){
436 fCurrentEventNegElectronTClone = new TClonesArray("AliESDtrack",0);
437 }
438 if(fKFReconstructedGammasCutTClone == NULL){
439 fKFReconstructedGammasCutTClone = new TClonesArray("AliKFParticle",0);
440 }
441 if(fPreviousEventTLVNegElectronTClone == NULL){
442 fPreviousEventTLVNegElectronTClone = new TClonesArray("TLorentzVector",0);
443 }
444 if(fPreviousEventTLVPosElectronTClone == NULL){
445 fPreviousEventTLVPosElectronTClone = new TClonesArray("TLorentzVector",0);
446 }
447 if(fChargedParticles == NULL){
448 fChargedParticles = new TClonesArray("AliESDtrack",0);
449 }
6272370b 450
451 if(fKFReconstructedPi0sTClone == NULL){
452 fKFReconstructedPi0sTClone = new TClonesArray("AliKFParticle",0);
453 }
9c1cb6f7 454
455 if(fKFRecalculatedGammasTClone == NULL){
456 fKFRecalculatedGammasTClone = new TClonesArray("AliKFParticle",0);
457 }
458
a0b94e5c 459
6c84d371 460 //clear TClones
1e7846f4 461 fKFReconstructedGammasTClone->Delete();
462 fCurrentEventPosElectronTClone->Delete();
463 fCurrentEventNegElectronTClone->Delete();
464 fKFReconstructedGammasCutTClone->Delete();
465 fPreviousEventTLVNegElectronTClone->Delete();
466 fPreviousEventTLVPosElectronTClone->Delete();
6272370b 467 fKFReconstructedPi0sTClone->Delete();
9c1cb6f7 468 fKFRecalculatedGammasTClone->Delete();
5e55d806 469
d7d7e825 470 //clear vectors
d7d7e825 471 fElectronv1.clear();
472 fElectronv2.clear();
6272370b 473
474 fGammav1.clear();
475 fGammav2.clear();
476
9c1cb6f7 477 fElectronRecalculatedv1.clear();
478 fElectronRecalculatedv2.clear();
479
d7d7e825 480
1e7846f4 481 fChargedParticles->Delete();
5e55d806 482
d7d7e825 483 fChargedParticlesId.clear();
037dc2db 484
485 fKFReconstructedGammasV0Index.clear();
a0b94e5c 486
d7d7e825 487 //Clear the data in the v0Reader
5e55d806 488 // fV0Reader->UpdateEventByEventData();
b5832f95 489
490 //Take Only events with proper trigger
5e55d806 491 /*
b5832f95 492 if(fTriggerCINT1B){
493 if(!fV0Reader->GetESDEvent()->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL")) return;
494 }
5e55d806 495 */
10e3319b 496
c8206114 497 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
498 // cout<< "Event not taken"<< endl;
3c45d101 499 eventQuality=1;
500 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
501 if(fDoMCTruth){
502 CheckMesonProcessTypeEventQuality(eventQuality);
503 }
c8206114 504 return; // aborts if the primary vertex does not have contributors.
505 }
506
507
885114d1 508 if(!fV0Reader->CheckForPrimaryVertexZ() ){
3c45d101 509 eventQuality=2;
510 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
511 if(fDoMCTruth){
512 CheckMesonProcessTypeEventQuality(eventQuality);
513 }
885114d1 514 return;
515 }
3c45d101 516 eventQuality=3;
517 fHistograms->FillHistogram("ESD_EventQuality",eventQuality);
885114d1 518
d7d7e825 519 // Process the MC information
520 if(fDoMCTruth){
521 ProcessMCData();
522 }
a0b94e5c 523
d7d7e825 524 //Process the v0 information with no cuts
525 ProcessV0sNoCut();
cb90a330 526
d7d7e825 527 // Process the v0 information
528 ProcessV0s();
cb90a330 529
5e55d806 530
d7d7e825 531 //Fill Gamma AOD
532 FillAODWithConversionGammas() ;
a0b94e5c 533
d7d7e825 534 // Process reconstructed gammas
535 if(fDoNeutralMeson == kTRUE){
a0b94e5c 536 ProcessGammasForNeutralMesonAnalysis();
77880bd8 537
48682642 538 }
539
77880bd8 540 if(fDoMCTruth == kTRUE){
541 CheckV0Efficiency();
542 }
d7d7e825 543 //Process reconstructed gammas electrons for Chi_c Analysis
6c84d371 544 if(fDoChic == kTRUE){
d7d7e825 545 ProcessGammaElectronsForChicAnalysis();
546 }
547 // Process reconstructed gammas for gamma Jet/hadron correlations
548 if(fDoJet == kTRUE){
549 ProcessGammasForGammaJetAnalysis();
550 }
a0b94e5c 551
5e55d806 552 //calculate background if flag is set
553 if(fCalculateBackground){
554 CalculateBackground();
555 }
48682642 556
557 if(fDoNeutralMeson == kTRUE){
558// ProcessConvPHOSGammasForNeutralMesonAnalysis();
559 if(fDoOmegaMeson == kTRUE){
560 ProcessGammasForOmegaMesonAnalysis();
561 }
562 }
563
5e55d806 564 //Clear the data in the v0Reader
565 fV0Reader->UpdateEventByEventData();
9c1cb6f7 566 if(fRecalculateV0ForGamma==kTRUE){
567 RecalculateV0ForGamma();
568 }
d7d7e825 569 PostData(1, fOutputContainer);
a0b94e5c 570 PostData(2, fCFManager->GetParticleContainer()); // for CF
d7d7e825 571
572}
573
48682642 574// void AliAnalysisTaskGammaConversion::ConnectInputData(Option_t *option){
575// // see header file for documentation
576// // printf(" ConnectInputData %s\n", GetName());
8a685cf3 577
48682642 578// AliAnalysisTaskSE::ConnectInputData(option);
8a685cf3 579
48682642 580// if(fV0Reader == NULL){
581// // Write warning here cuts and so on are default if this ever happens
582// }
583// fV0Reader->Initialize();
584// fDoMCTruth = fV0Reader->GetDoMCTruth();
585// }
d7d7e825 586
3c45d101 587void AliAnalysisTaskGammaConversion::CheckMesonProcessTypeEventQuality(Int_t evtQ){
588 fStack= MCEvent()->Stack();
589 fGCMCEvent=MCEvent();
590
591 for (Int_t iTracks = 0; iTracks < fStack->GetNprimary(); iTracks++) {
592 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
593 if (!particle) {
594 //print warning here
595 continue;
596 }
597 if(particle->GetPdgCode()!=111){ //Pi0
598 continue;
599 }
600 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
601 if(evtQ==1){
602 switch(GetProcessType(fGCMCEvent)){
603 case kProcSD:
604 fHistograms->FillHistogram("MC_SD_EvtQ1_Pi0_Pt", particle->Pt());
605 break;
606 case kProcDD:
607 fHistograms->FillHistogram("MC_DD_EvtQ1_Pi0_Pt", particle->Pt());
608 break;
609 case kProcND:
610 fHistograms->FillHistogram("MC_ND_EvtQ1_Pi0_Pt", particle->Pt());
611 break;
612 default:
613 AliError("Unknown Process");
614 }
615 }
616 if(evtQ==2){
617 switch(GetProcessType(fGCMCEvent)){
618 case kProcSD:
619 fHistograms->FillHistogram("MC_SD_EvtQ2_Pi0_Pt", particle->Pt());
620 break;
621 case kProcDD:
622 fHistograms->FillHistogram("MC_DD_EvtQ2_Pi0_Pt", particle->Pt());
623 break;
624 case kProcND:
625 fHistograms->FillHistogram("MC_ND_EvtQ2_Pi0_Pt", particle->Pt());
626 break;
627 default:
628 AliError("Unknown Process");
629 }
630 }
d7d7e825 631
632
3c45d101 633 }
634
635}
636
d7d7e825 637void AliAnalysisTaskGammaConversion::ProcessMCData(){
638 // see header file for documentation
48682642 639 //InputEvent(), MCEvent());
640 /* TestAnaMarin
d7d7e825 641 fStack = fV0Reader->GetMCStack();
a0b94e5c 642 fMCTruth = fV0Reader->GetMCTruth(); // for CF
4a6157dc 643 fGCMCEvent = fV0Reader->GetMCEvent(); // for CF
48682642 644 */
645 fStack= MCEvent()->Stack();
646 fGCMCEvent=MCEvent();
a0b94e5c 647
648 // for CF
1e7846f4 649 Double_t containerInput[3];
650 if(fDoCF){
651 if(!fGCMCEvent) cout << "NO MC INFO FOUND" << endl;
652 fCFManager->SetEventInfo(fGCMCEvent);
653 }
a0b94e5c 654 // end for CF
655
656
d7d7e825 657 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
658 return; // aborts if the primary vertex does not have contributors.
659 }
bd6d9fa3 660
d7d7e825 661 for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++) {
662 TParticle* particle = (TParticle *)fStack->Particle(iTracks);
bd6d9fa3 663
111d75df 664
bd6d9fa3 665
d7d7e825 666 if (!particle) {
667 //print warning here
668 continue;
669 }
a0b94e5c 670
d7d7e825 671 ///////////////////////Begin Chic Analysis/////////////////////////////
a0b94e5c 672
a68437fb 673 if(particle->GetPdgCode() == 443){//Is JPsi
d7d7e825 674 if(particle->GetNDaughters()==2){
675 if(TMath::Abs(fStack->Particle(particle->GetFirstDaughter())->GetPdgCode()) == 11 &&
676 TMath::Abs(fStack->Particle(particle->GetLastDaughter())->GetPdgCode()) == 11){
63e16c52 677
d7d7e825 678 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
679 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
680 if(TMath::Abs(daug0->Eta()) < 0.9 && TMath::Abs(daug1->Eta()) < 0.9)
681 fHistograms->FillTable("Table_Electrons",3);//e+ e- from J/Psi inside acceptance
a0b94e5c 682
d7d7e825 683 if( TMath::Abs(daug0->Eta()) < 0.9){
684 if(daug0->GetPdgCode() == -11)
685 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
686 else
687 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
a0b94e5c 688
d7d7e825 689 }
690 if(TMath::Abs(daug1->Eta()) < 0.9){
691 if(daug1->GetPdgCode() == -11)
692 fHistograms->FillTable("Table_Electrons",1);//e+ from J/Psi inside acceptance
693 else
694 fHistograms->FillTable("Table_Electrons",2);//e- from J/Psi inside acceptance
695 }
696 }
697 }
698 }
699 // const int CHI_C0 = 10441;
700 // const int CHI_C1 = 20443;
701 // const int CHI_C2 = 445
702 if(particle->GetPdgCode() == 22){//gamma from JPsi
703 if(particle->GetMother(0) > -1){
704 if(fStack->Particle(particle->GetMother(0))->GetPdgCode() == 10441 ||
705 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 20443 ||
706 fStack->Particle(particle->GetMother(0))->GetPdgCode() == 445){
707 if(TMath::Abs(particle->Eta()) < 1.2)
708 fHistograms->FillTable("Table_Electrons",17);// gamma from chic inside accptance
709 }
710 }
711 }
712 if(particle->GetPdgCode() == 10441 || particle->GetPdgCode() == 20443 || particle->GetPdgCode() == 445){
713 if( particle->GetNDaughters() == 2){
714 TParticle* daug0 = fStack->Particle(particle->GetFirstDaughter());
715 TParticle* daug1 = fStack->Particle(particle->GetLastDaughter());
a0b94e5c 716
d7d7e825 717 if( (daug0->GetPdgCode() == 443 || daug0->GetPdgCode() == 22) && (daug1->GetPdgCode() == 443 || daug1->GetPdgCode() == 22) ){
718 if( daug0->GetPdgCode() == 443){
719 TParticle* daugE0 = fStack->Particle(daug0->GetFirstDaughter());
720 TParticle* daugE1 = fStack->Particle(daug0->GetLastDaughter());
721 if( TMath::Abs(daug1->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
722 fHistograms->FillTable("Table_Electrons",18);
a0b94e5c 723
d7d7e825 724 }//if
725 else if (daug1->GetPdgCode() == 443){
726 TParticle* daugE0 = fStack->Particle(daug1->GetFirstDaughter());
727 TParticle* daugE1 = fStack->Particle(daug1->GetLastDaughter());
728 if( TMath::Abs(daug0->Eta()) < 1.2 && TMath::Abs(daugE0->Eta()) < 0.9 && TMath::Abs(daugE1->Eta()) < 0.9 )
729 fHistograms->FillTable("Table_Electrons",18);
730 }//else if
731 }//gamma o Jpsi
732 }//GetNDaughters
733 }
a0b94e5c 734
735
d7d7e825 736 /////////////////////End Chic Analysis////////////////////////////
a0b94e5c 737
738
dc2883e4 739 // if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
a0b94e5c 740
d7d7e825 741 if(particle->R()>fV0Reader->GetMaxRCut()) continue; // cuts on distance from collision point
742
743 Double_t tmpPhi=particle->Phi();
744
745 if(particle->Phi()> TMath::Pi()){
746 tmpPhi = particle->Phi()-(2*TMath::Pi());
747 }
748
749 Double_t rapidity;
750 if(particle->Energy() - particle->Pz() == 0 || particle->Energy() + particle->Pz() == 0){
751 rapidity=0;
752 }
753 else{
754 rapidity = 0.5*(TMath::Log((particle->Energy()+particle->Pz()) / (particle->Energy()-particle->Pz())));
755 }
756
111d75df 757
758
759 if(iTracks<=fStack->GetNprimary() ){
760 if ( particle->GetPdgCode()== -211 || particle->GetPdgCode()== 211 ||
761 particle->GetPdgCode()== 2212 || particle->GetPdgCode()==-2212 ||
762 particle->GetPdgCode()== 321 || particle->GetPdgCode()==-321 ){
dc2883e4 763 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
3c45d101 764 fHistograms->FillHistogram("MC_PhysicalPrimaryCharged_Pt", particle->Pt());
111d75df 765 }
bd6d9fa3 766 }
111d75df 767
bd6d9fa3 768
769
d7d7e825 770 //process the gammas
771 if (particle->GetPdgCode() == 22){
dc2883e4 772 if(TMath::Abs(particle->Eta())> fV0Reader->GetEtaCut() ) continue;
a607c218 773
d7d7e825 774 if(particle->GetMother(0) >-1 && fStack->Particle(particle->GetMother(0))->GetPdgCode() == 22){
775 continue; // no photon as mothers!
776 }
a607c218 777
d7d7e825 778 if(particle->GetMother(0) >= fStack->GetNprimary()){
779 continue; // the gamma has a mother, and it is not a primary particle
780 }
f8017b04 781
782 if(particle->GetMother(0) >-1){
783 fHistograms->FillHistogram("MC_DecayAllGamma_Pt", particle->Pt()); // All
784 switch(fStack->Particle(particle->GetMother(0))->GetPdgCode()){
785 case 111: // Pi0
786 fHistograms->FillHistogram("MC_DecayPi0Gamma_Pt", particle->Pt());
787 break;
788 case 113: // Rho0
789 fHistograms->FillHistogram("MC_DecayRho0Gamma_Pt", particle->Pt());
790 break;
791 case 221: // Eta
792 fHistograms->FillHistogram("MC_DecayEtaGamma_Pt", particle->Pt());
793 break;
794 case 223: // Omega
795 fHistograms->FillHistogram("MC_DecayOmegaGamma_Pt", particle->Pt());
796 break;
797 case 310: // K_s0
798 fHistograms->FillHistogram("MC_DecayK0sGamma_Pt", particle->Pt());
799 break;
800 case 331: // Eta'
801 fHistograms->FillHistogram("MC_DecayEtapGamma_Pt", particle->Pt());
802 break;
803 }
804 }
805
d7d7e825 806 fHistograms->FillHistogram("MC_allGamma_Energy", particle->Energy());
807 fHistograms->FillHistogram("MC_allGamma_Pt", particle->Pt());
808 fHistograms->FillHistogram("MC_allGamma_Eta", particle->Eta());
809 fHistograms->FillHistogram("MC_allGamma_Phi", tmpPhi);
810 fHistograms->FillHistogram("MC_allGamma_Rapid", rapidity);
811
a0b94e5c 812 // for CF
1e7846f4 813 if(fDoCF){
814 containerInput[0] = particle->Pt();
815 containerInput[1] = particle->Eta();
816 if(particle->GetMother(0) >=0){
817 containerInput[2] = fStack->Particle(particle->GetMother(0))->GetMass();
818 }
819 else{
820 containerInput[2]=-1;
821 }
a607c218 822
1e7846f4 823 fCFManager->GetParticleContainer()->Fill(containerInput,kStepGenerated); // generated gamma
824 }
d7d7e825 825 if(particle->GetMother(0) < 0){ // direct gamma
826 fHistograms->FillHistogram("MC_allDirectGamma_Energy",particle->Energy());
827 fHistograms->FillHistogram("MC_allDirectGamma_Pt", particle->Pt());
828 fHistograms->FillHistogram("MC_allDirectGamma_Eta", particle->Eta());
829 fHistograms->FillHistogram("MC_allDirectGamma_Phi", tmpPhi);
830 fHistograms->FillHistogram("MC_allDirectGamma_Rapid", rapidity);
831 }
832
d7d7e825 833 // looking for conversion (electron + positron from pairbuilding (= 5) )
834 TParticle* ePos = NULL;
835 TParticle* eNeg = NULL;
836
837 if(particle->GetNDaughters() >= 2){
838 for(Int_t daughterIndex=particle->GetFirstDaughter();daughterIndex<=particle->GetLastDaughter();daughterIndex++){
839 TParticle *tmpDaughter = fStack->Particle(daughterIndex);
840 if(tmpDaughter->GetUniqueID() == 5){
841 if(tmpDaughter->GetPdgCode() == 11){
842 eNeg = tmpDaughter;
843 }
844 else if(tmpDaughter->GetPdgCode() == -11){
845 ePos = tmpDaughter;
846 }
847 }
848 }
849 }
850
851
852 if(ePos == NULL || eNeg == NULL){ // means we do not have two daughters from pair production
853 continue;
854 }
855
856
857 Double_t ePosPhi = ePos->Phi();
858 if(ePos->Phi()> TMath::Pi()) ePosPhi = ePos->Phi()-(2*TMath::Pi());
859
860 Double_t eNegPhi = eNeg->Phi();
861 if(eNeg->Phi()> TMath::Pi()) eNegPhi = eNeg->Phi()-(2*TMath::Pi());
862
863
864 if(ePos->Pt()<fV0Reader->GetPtCut() || eNeg->Pt()<fV0Reader->GetPtCut()){
865 continue; // no reconstruction below the Pt cut
866 }
a0b94e5c 867
d7d7e825 868 if(TMath::Abs(ePos->Eta())> fV0Reader->GetEtaCut() || TMath::Abs(eNeg->Eta())> fV0Reader->GetEtaCut()){
869 continue;
870 }
a0b94e5c 871
d7d7e825 872 if(ePos->R()>fV0Reader->GetMaxRCut()){
873 continue; // cuts on distance from collision point
874 }
a0b94e5c 875
876 if(TMath::Abs(ePos->Vz()) > fV0Reader->GetMaxZCut()){
d7d7e825 877 continue; // outside material
878 }
879
a0b94e5c 880
d7d7e825 881 if((TMath::Abs(ePos->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() > ePos->R()){
882 continue; // line cut to exclude regions where we do not reconstruct
883 }
884
a0b94e5c 885
886 // for CF
1e7846f4 887 if(fDoCF){
888 fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructable); // reconstructable gamma
889 }
d7d7e825 890 fHistograms->FillHistogram("MC_ConvGamma_Energy", particle->Energy());
891 fHistograms->FillHistogram("MC_ConvGamma_Pt", particle->Pt());
892 fHistograms->FillHistogram("MC_ConvGamma_Eta", particle->Eta());
893 fHistograms->FillHistogram("MC_ConvGamma_Phi", tmpPhi);
894 fHistograms->FillHistogram("MC_ConvGamma_Rapid", rapidity);
895 fHistograms->FillHistogram("MC_ConvGamma_Pt_Eta", particle->Pt(),particle->Eta());
896
897 fHistograms->FillHistogram("MC_E_Energy", eNeg->Energy());
898 fHistograms->FillHistogram("MC_E_Pt", eNeg->Pt());
899 fHistograms->FillHistogram("MC_E_Eta", eNeg->Eta());
900 fHistograms->FillHistogram("MC_E_Phi", eNegPhi);
901
902 fHistograms->FillHistogram("MC_P_Energy", ePos->Energy());
903 fHistograms->FillHistogram("MC_P_Pt", ePos->Pt());
904 fHistograms->FillHistogram("MC_P_Eta", ePos->Eta());
905 fHistograms->FillHistogram("MC_P_Phi", ePosPhi);
906
907
d7d7e825 908 // begin Mapping
909 Int_t rBin = fHistograms->GetRBin(ePos->R());
9640a3d1 910 Int_t zBin = fHistograms->GetZBin(ePos->Vz());
d7d7e825 911 Int_t phiBin = fHistograms->GetPhiBin(particle->Phi());
9c1cb6f7 912 Double_t rFMD=30;
9640a3d1 913
914 TVector3 vtxPos(ePos->Vx(),ePos->Vy(),ePos->Vz());
915
d7d7e825 916 TString nameMCMappingPhiR="";
e158cbc3 917 nameMCMappingPhiR.Form("MC_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 918 // fHistograms->FillHistogram(nameMCMappingPhiR, ePos->Vz(), particle->Eta());
d7d7e825 919
920 TString nameMCMappingPhi="";
e158cbc3 921 nameMCMappingPhi.Form("MC_Conversion_Mapping_Phi%02d",phiBin);
9640a3d1 922 // fHistograms->FillHistogram(nameMCMappingPhi, particle->Eta());
ebcfaa7e 923 //fHistograms->FillHistogram(nameMCMappingPhi, ePos->Vz(), particle->Eta());
d7d7e825 924
925 TString nameMCMappingR="";
e158cbc3 926 nameMCMappingR.Form("MC_Conversion_Mapping_R%02d",rBin);
9640a3d1 927 // fHistograms->FillHistogram(nameMCMappingR, particle->Eta());
ebcfaa7e 928 //fHistograms->FillHistogram(nameMCMappingR,ePos->Vz(), particle->Eta());
d7d7e825 929
930 TString nameMCMappingPhiInR="";
e158cbc3 931 nameMCMappingPhiInR.Form("MC_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 932 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
933 fHistograms->FillHistogram(nameMCMappingPhiInR, vtxPos.Phi());
934
935 TString nameMCMappingZInR="";
936 nameMCMappingZInR.Form("MC_Conversion_Mapping_Z_in_R_%02d",rBin);
937 fHistograms->FillHistogram(nameMCMappingZInR,ePos->Vz() );
938
939
940 TString nameMCMappingPhiInZ="";
941 nameMCMappingPhiInZ.Form("MC_Conversion_Mapping_Phi_in_Z_%02d",zBin);
942 // fHistograms->FillHistogram(nameMCMappingPhiInR, tmpPhi);
943 fHistograms->FillHistogram(nameMCMappingPhiInZ, vtxPos.Phi());
944
9c1cb6f7 945
946 if(ePos->R()<rFMD){
947 TString nameMCMappingFMDPhiInZ="";
948 nameMCMappingFMDPhiInZ.Form("MC_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
949 fHistograms->FillHistogram(nameMCMappingFMDPhiInZ, vtxPos.Phi());
950 }
951
9640a3d1 952 TString nameMCMappingRInZ="";
953 nameMCMappingRInZ.Form("MC_Conversion_Mapping_R_in_Z_%02d",zBin);
954 fHistograms->FillHistogram(nameMCMappingRInZ,ePos->R() );
955
956 if(particle->Pt() > fLowPtMapping && particle->Pt()< fHighPtMapping){
957 TString nameMCMappingMidPtPhiInR="";
958 nameMCMappingMidPtPhiInR.Form("MC_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
959 fHistograms->FillHistogram(nameMCMappingMidPtPhiInR, vtxPos.Phi());
960
961 TString nameMCMappingMidPtZInR="";
962 nameMCMappingMidPtZInR.Form("MC_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
963 fHistograms->FillHistogram(nameMCMappingMidPtZInR,ePos->Vz() );
964
965
966 TString nameMCMappingMidPtPhiInZ="";
967 nameMCMappingMidPtPhiInZ.Form("MC_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
968 fHistograms->FillHistogram(nameMCMappingMidPtPhiInZ, vtxPos.Phi());
9c1cb6f7 969
970
971 if(ePos->R()<rFMD){
972 TString nameMCMappingMidPtFMDPhiInZ="";
973 nameMCMappingMidPtFMDPhiInZ.Form("MC_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
974 fHistograms->FillHistogram(nameMCMappingMidPtFMDPhiInZ, vtxPos.Phi());
975 }
976
9640a3d1 977
978 TString nameMCMappingMidPtRInZ="";
979 nameMCMappingMidPtRInZ.Form("MC_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
980 fHistograms->FillHistogram(nameMCMappingMidPtRInZ,ePos->R() );
981
982 }
983
d7d7e825 984 //end mapping
985
986 fHistograms->FillHistogram("MC_Conversion_R",ePos->R());
987 fHistograms->FillHistogram("MC_Conversion_ZR",ePos->Vz(),ePos->R());
988 fHistograms->FillHistogram("MC_Conversion_XY",ePos->Vx(),ePos->Vy());
989 fHistograms->FillHistogram("MC_Conversion_OpeningAngle",GetMCOpeningAngle(ePos, eNeg));
9640a3d1 990 fHistograms->FillHistogram("MC_ConvGamma_E_AsymmetryP",particle->P(),eNeg->P()/particle->P());
991 fHistograms->FillHistogram("MC_ConvGamma_P_AsymmetryP",particle->P(),ePos->P()/particle->P());
992
993
d7d7e825 994 if(particle->GetMother(0) < 0){ // no mother = direct gamma, still inside converted
995 fHistograms->FillHistogram("MC_ConvDirectGamma_Energy",particle->Energy());
996 fHistograms->FillHistogram("MC_ConvDirectGamma_Pt", particle->Pt());
997 fHistograms->FillHistogram("MC_ConvDirectGamma_Eta", particle->Eta());
998 fHistograms->FillHistogram("MC_ConvDirectGamma_Phi", tmpPhi);
999 fHistograms->FillHistogram("MC_ConvDirectGamma_Rapid", rapidity);
1000
1001 } // end direct gamma
1002 else{ // mother exits
6c84d371 1003 /* if( fStack->Particle(particle->GetMother(0))->GetPdgCode()==10441 ||//chic0
a0b94e5c 1004 fStack->Particle(particle->GetMother(0))->GetPdgCode()==20443 ||//psi2S
1005 fStack->Particle(particle->GetMother(0))->GetPdgCode()==445 //chic2
1006 ){
1007 fMCGammaChic.push_back(particle);
1008 }
6c84d371 1009 */
d7d7e825 1010 } // end if mother exits
1011 } // end if particle is a photon
1012
1013
1014
1015 // process motherparticles (2 gammas as daughters)
1016 // the motherparticle had already to pass the R and the eta cut, but no line cut.
1017 // the line cut is just valid for the conversions!
1018
1019 if(particle->GetNDaughters() == 2){
1020
1021 TParticle* daughter0 = (TParticle*)fStack->Particle(particle->GetFirstDaughter());
1022 TParticle* daughter1 = (TParticle*)fStack->Particle(particle->GetLastDaughter());
1023
1024 if(daughter0->GetPdgCode() != 22 || daughter1->GetPdgCode() != 22) continue; //check for gamma gamma daughters
dc2883e4 1025
1026 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ) continue;
1027
d7d7e825 1028 // Check the acceptance for both gammas
01b7fdcc 1029 Bool_t gammaEtaCut = kTRUE;
1030 if(TMath::Abs(daughter0->Eta()) > fV0Reader->GetEtaCut() || TMath::Abs(daughter1->Eta()) > fV0Reader->GetEtaCut() ) gammaEtaCut = kFALSE;
1031 Bool_t gammaRCut = kTRUE;
1032 if(daughter0->R() > fV0Reader->GetMaxRCut() || daughter1->R() > fV0Reader->GetMaxRCut() ) gammaRCut = kFALSE;
d7d7e825 1033
1034
1035
1036 // check for conversions now -> have to pass eta, R and line cut!
1037 Bool_t daughter0Electron = kFALSE;
1038 Bool_t daughter0Positron = kFALSE;
1039 Bool_t daughter1Electron = kFALSE;
1040 Bool_t daughter1Positron = kFALSE;
1041
1042 if(daughter0->GetNDaughters() >= 2){ // first gamma
1043 for(Int_t TrackIndex=daughter0->GetFirstDaughter();TrackIndex<=daughter0->GetLastDaughter();TrackIndex++){
1044 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
1045 if(tmpDaughter->GetUniqueID() == 5){
1046 if(tmpDaughter->GetPdgCode() == 11){
1047 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1048 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1049 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1050 daughter0Electron = kTRUE;
1051 }
1052 }
1053
1054 }
1055 }
1056 else if(tmpDaughter->GetPdgCode() == -11){
1057 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1058 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1059 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1060 daughter0Positron = kTRUE;
1061 }
a68437fb 1062 }
1063 }
d7d7e825 1064 }
1065 }
1066 }
1067 }
1068
1069
1070 if(daughter1->GetNDaughters() >= 2){ // second gamma
1071 for(Int_t TrackIndex=daughter1->GetFirstDaughter();TrackIndex<=daughter1->GetLastDaughter();TrackIndex++){
1072 TParticle *tmpDaughter = fStack->Particle(TrackIndex);
1073 if(tmpDaughter->GetUniqueID() == 5){
1074 if(tmpDaughter->GetPdgCode() == 11){
1075 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1076 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1077 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1078 daughter1Electron = kTRUE;
1079 }
1080 }
1081
1082 }
1083 }
1084 else if(tmpDaughter->GetPdgCode() == -11){
1085 if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() ){
1086 if( ( TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue() < tmpDaughter->R() ){
1087 if(tmpDaughter->R()< fV0Reader->GetMaxRCut()){
1088 daughter1Positron = kTRUE;
1089 }
1090 }
1091
1092 }
1093
1094 }
1095 }
1096 }
1097 }
1098
a0b94e5c 1099
d7d7e825 1100 if(particle->GetPdgCode()==111){ //Pi0
1101 if( iTracks >= fStack->GetNprimary()){
1102 fHistograms->FillHistogram("MC_Pi0_Secondaries_Eta", particle->Eta());
1103 fHistograms->FillHistogram("MC_Pi0_Secondaries_Rapid", rapidity);
1104 fHistograms->FillHistogram("MC_Pi0_Secondaries_Phi", tmpPhi);
1105 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt", particle->Pt());
1106 fHistograms->FillHistogram("MC_Pi0_Secondaries_Energy", particle->Energy());
1107 fHistograms->FillHistogram("MC_Pi0_Secondaries_R", particle->R());
1108 fHistograms->FillHistogram("MC_Pi0_Secondaries_ZR", particle->Vz(),particle->R());
1109 fHistograms->FillHistogram("MC_Pi0_Secondaries_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1110 fHistograms->FillHistogram("MC_Pi0_Secondaries_XY", particle->Vx(),particle->Vy());//only fill from one daughter to avoid multiple filling
1111
01b7fdcc 1112 if(gammaEtaCut && gammaRCut){
d7d7e825 1113 //if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
1114 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1115 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1116 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1117 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1118 fHistograms->FillHistogram("MC_Pi0_Secondaries_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1119 }
1120 }
1121 }
1122 else{
1123 fHistograms->FillHistogram("MC_Pi0_Eta", particle->Eta());
1124 fHistograms->FillHistogram("MC_Pi0_Rapid", rapidity);
1125 fHistograms->FillHistogram("MC_Pi0_Phi", tmpPhi);
1126 fHistograms->FillHistogram("MC_Pi0_Pt", particle->Pt());
dc2883e4 1127 fHistograms->FillHistogram("MC_Pi0_Pt_vs_Rapid", particle->Pt(),rapidity);
d7d7e825 1128 fHistograms->FillHistogram("MC_Pi0_Energy", particle->Energy());
1129 fHistograms->FillHistogram("MC_Pi0_R", particle->R());
1130 fHistograms->FillHistogram("MC_Pi0_ZR", particle->Vz(),particle->R());
1131 fHistograms->FillHistogram("MC_Pi0_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1132 fHistograms->FillHistogram("MC_Pi0_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
63e16c52 1133 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_Fiducial", particle->Pt());
3c45d101 1134
1135 switch(GetProcessType(fGCMCEvent)){
1136 case kProcSD:
1137 fHistograms->FillHistogram("MC_SD_EvtQ3_Pi0_Pt", particle->Pt());
1138 break;
1139 case kProcDD:
1140 fHistograms->FillHistogram("MC_DD_EvtQ3_Pi0_Pt", particle->Pt());
1141 break;
1142 case kProcND:
1143 fHistograms->FillHistogram("MC_ND_EvtQ3_Pi0_Pt", particle->Pt());
1144 break;
1145 default:
1146 AliError("Unknown Process");
1147 }
63e16c52 1148
01b7fdcc 1149 if(gammaEtaCut && gammaRCut){
d7d7e825 1150 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
1151 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1152 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
63e16c52 1153 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_withinAcceptance_Fiducial", particle->Pt());
1154
d7d7e825 1155 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1156 fHistograms->FillHistogram("MC_Pi0_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1157 fHistograms->FillHistogram("MC_Pi0_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1158 fHistograms->FillHistogram("MC_Pi0_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
10e3319b 1159 fHistograms->FillHistogram("MC_Pi0_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
1160 fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
1161 fHistograms->FillHistogram("MC_Pi0_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
1162
48682642 1163 Double_t alfa=0.;
1164 if((daughter0->Energy()+daughter1->Energy())!= 0.){
1165 alfa= TMath::Abs((daughter0->Energy()-daughter1->Energy())/(daughter0->Energy()+daughter1->Energy()));
1166 }
1167 fHistograms->FillHistogram("MC_Pi0_alpha",alfa);
1168 if(TMath::Abs(particle->Eta())<0.9)fHistograms->FillHistogram("MC_Pi0_Pt_ConvGamma_withinAcceptance_Fiducial", particle->Pt());
1169
d7d7e825 1170 }
1171 }
1172 }
1173 }
1174
1175 if(particle->GetPdgCode()==221){ //Eta
1176 fHistograms->FillHistogram("MC_Eta_Eta", particle->Eta());
1177 fHistograms->FillHistogram("MC_Eta_Rapid", rapidity);
1178 fHistograms->FillHistogram("MC_Eta_Phi",tmpPhi);
1179 fHistograms->FillHistogram("MC_Eta_Pt", particle->Pt());
dc2883e4 1180 fHistograms->FillHistogram("MC_Eta_Pt_vs_Rapid", particle->Pt(),rapidity);
d7d7e825 1181 fHistograms->FillHistogram("MC_Eta_Energy", particle->Energy());
1182 fHistograms->FillHistogram("MC_Eta_R", particle->R());
1183 fHistograms->FillHistogram("MC_Eta_ZR", particle->Vz(),particle->R());
1184 fHistograms->FillHistogram("MC_Eta_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1185 fHistograms->FillHistogram("MC_Eta_XY", particle->Vx(), particle->Vy());//only fill from one daughter to avoid multiple filling
1186
01b7fdcc 1187 if(gammaEtaCut && gammaRCut){
a0b94e5c 1188 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 1189 fHistograms->FillHistogram("MC_Eta_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1190 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1191 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1192 fHistograms->FillHistogram("MC_Eta_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1193 fHistograms->FillHistogram("MC_Eta_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1194 fHistograms->FillHistogram("MC_Eta_ZR_ConvGamma_withinAcceptance", particle->Vz(),particle->R());
10e3319b 1195 fHistograms->FillHistogram("MC_Eta_ConvGamma_OpeningAngle_Pt", particle->Pt(),GetMCOpeningAngle(daughter0,daughter1));
1196 fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter0->Pt());
1197 fHistograms->FillHistogram("MC_Eta_ConvGamma_PtGamma_Pt", particle->Pt(),daughter1->Pt());
1198
d7d7e825 1199 }
1200
1201 }
1202
1203 }
1204
1205 // all motherparticles with 2 gammas as daughters
1206 fHistograms->FillHistogram("MC_Mother_R", particle->R());
1207 fHistograms->FillHistogram("MC_Mother_ZR", particle->Vz(),particle->R());
1208 fHistograms->FillHistogram("MC_Mother_XY", particle->Vx(),particle->Vy());
1209 fHistograms->FillHistogram("MC_Mother_Mass", particle->GetCalcMass());
1210 fHistograms->FillHistogram("MC_Mother_GammaDaughter_OpeningAngle", GetMCOpeningAngle(daughter0,daughter1));
1211 fHistograms->FillHistogram("MC_Mother_Energy", particle->Energy());
1212 fHistograms->FillHistogram("MC_Mother_Pt", particle->Pt());
1213 fHistograms->FillHistogram("MC_Mother_Eta", particle->Eta());
1214 fHistograms->FillHistogram("MC_Mother_Rapid", rapidity);
1215 fHistograms->FillHistogram("MC_Mother_Phi",tmpPhi);
1216 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt",particle->GetMass(),particle->Pt());
a0b94e5c 1217
01b7fdcc 1218 if(gammaEtaCut && gammaRCut){
a0b94e5c 1219 // if(TMath::Abs(daughter0->Eta()) <= fV0Reader->GetEtaCut() && TMath::Abs(daughter1->Eta()) <= fV0Reader->GetEtaCut() ){
d7d7e825 1220 fHistograms->FillHistogram("MC_Mother_Pt_Eta_withinAcceptance", particle->Pt(),particle->Eta());
1221 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_withinAcceptance", particle->Pt(),rapidity);
1222 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_withinAcceptance",particle->GetMass(),particle->Pt());
1223 if(daughter0Electron && daughter0Positron && daughter1Electron && daughter1Positron){
1224 fHistograms->FillHistogram("MC_Mother_Pt_Eta_ConvGamma_withinAcceptance", particle->Pt(),particle->Eta());
1225 fHistograms->FillHistogram("MC_Mother_Pt_Rapid_ConvGamma_withinAcceptance", particle->Pt(),rapidity);
1226 fHistograms->FillHistogram("MC_Mother_InvMass_vs_Pt_ConvGamma_withinAcceptance",particle->GetMass(),particle->Pt());
a0b94e5c 1227
d7d7e825 1228 }
1229
1230
1231 } // end passed R and eta cut
a0b94e5c 1232
d7d7e825 1233 } // end if(particle->GetNDaughters() == 2)
1234
1235 }// end for (Int_t iTracks = 0; iTracks < fStack->GetNtrack(); iTracks++)
a0b94e5c 1236
d7d7e825 1237} // end ProcessMCData
1238
1239
1240
1241void AliAnalysisTaskGammaConversion::FillNtuple(){
1242 //Fills the ntuple with the different values
a0b94e5c 1243
d7d7e825 1244 if(fGammaNtuple == NULL){
1245 return;
1246 }
1247 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1248 for(Int_t i=0;i<numberOfV0s;i++){
1249 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};
1250 AliESDv0 * cV0 = fV0Reader->GetV0(i);
1251 Double_t negPID=0;
1252 Double_t posPID=0;
1253 fV0Reader->GetPIDProbability(negPID,posPID);
1254 values[0]=cV0->GetOnFlyStatus();
1255 values[1]=fV0Reader->CheckForPrimaryVertex();
1256 values[2]=negPID;
1257 values[3]=posPID;
1258 values[4]=fV0Reader->GetX();
1259 values[5]=fV0Reader->GetY();
1260 values[6]=fV0Reader->GetZ();
1261 values[7]=fV0Reader->GetXYRadius();
1262 values[8]=fV0Reader->GetMotherCandidateNDF();
1263 values[9]=fV0Reader->GetMotherCandidateChi2();
1264 values[10]=fV0Reader->GetMotherCandidateEnergy();
1265 values[11]=fV0Reader->GetMotherCandidateEta();
1266 values[12]=fV0Reader->GetMotherCandidatePt();
1267 values[13]=fV0Reader->GetMotherCandidateMass();
1268 values[14]=fV0Reader->GetMotherCandidateWidth();
1269 // values[15]=fV0Reader->GetMotherMCParticle()->Pt(); MOVED TO THE END, HAS TO BE CALLED AFTER HasSameMother NB: still has the same entry in the array
1270 values[16]=fV0Reader->GetOpeningAngle();
1271 values[17]=fV0Reader->GetNegativeTrackEnergy();
1272 values[18]=fV0Reader->GetNegativeTrackPt();
1273 values[19]=fV0Reader->GetNegativeTrackEta();
1274 values[20]=fV0Reader->GetNegativeTrackPhi();
1275 values[21]=fV0Reader->GetPositiveTrackEnergy();
1276 values[22]=fV0Reader->GetPositiveTrackPt();
1277 values[23]=fV0Reader->GetPositiveTrackEta();
1278 values[24]=fV0Reader->GetPositiveTrackPhi();
1279 values[25]=fV0Reader->HasSameMCMother();
1280 if(values[25] != 0){
1281 values[26]=fV0Reader->GetMotherMCParticlePDGCode();
1282 values[15]=fV0Reader->GetMotherMCParticle()->Pt();
1283 }
1284 fTotalNumberOfAddedNtupleEntries++;
1285 fGammaNtuple->Fill(values);
1286 }
1287 fV0Reader->ResetV0IndexNumber();
1288
1289}
1290
1291void AliAnalysisTaskGammaConversion::ProcessV0sNoCut(){
1292 // Process all the V0's without applying any cuts to it
a0b94e5c 1293
d7d7e825 1294 Int_t numberOfV0s = fV0Reader->GetNumberOfV0s();
1295 for(Int_t i=0;i<numberOfV0s;i++){
1296 /*AliESDv0 * cV0 = */fV0Reader->GetV0(i);
a0b94e5c 1297
d7d7e825 1298 if(fV0Reader->CheckForPrimaryVertex() == kFALSE){
cb90a330 1299 continue;
d7d7e825 1300 }
9640a3d1 1301
77880bd8 1302 // if( !fV0Reader->GetV0(i)->GetOnFlyStatus()){
1303 if( !fV0Reader->CheckV0FinderStatus(i)){
cb90a330 1304 continue;
9640a3d1 1305 }
1306
cb90a330 1307
1308 if( !((fV0Reader->GetNegativeESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ||
1309 !((fV0Reader->GetPositiveESDTrack())->GetStatus() & AliESDtrack::kTPCrefit) ){
1310 continue;
9640a3d1 1311 }
cb90a330 1312
ebcfaa7e 1313 if( fV0Reader->GetNegativeESDTrack()->GetSign()== fV0Reader->GetPositiveESDTrack()->GetSign()){
1314 continue;
1315 }
1316
1317 if( fV0Reader->GetNegativeESDTrack()->GetKinkIndex(0) > 0 ||
1318 fV0Reader->GetPositiveESDTrack()->GetKinkIndex(0) > 0) {
1319 continue;
1320 }
9c1cb6f7 1321 if(TMath::Abs(fV0Reader->GetMotherCandidateEta())> fV0Reader->GetEtaCut()){
1322 continue;
1323 }
1324 if(TMath::Abs(fV0Reader->GetPositiveTrackEta())> fV0Reader->GetEtaCut()){
1325 continue;
1326 }
1327 if(TMath::Abs(fV0Reader->GetNegativeTrackEta())> fV0Reader->GetEtaCut()){
1328 continue;
1329 }
1330 if((TMath::Abs(fV0Reader->GetZ())*fV0Reader->GetLineCutZRSlope())-fV0Reader->GetLineCutZValue() > fV0Reader->GetXYRadius() ){ // cuts out regions where we do not reconstruct
1331 continue;
1332 }
d7d7e825 1333 if(fDoMCTruth){
a0b94e5c 1334
d7d7e825 1335 if(fV0Reader->HasSameMCMother() == kFALSE){
1336 continue;
1337 }
a0b94e5c 1338
d7d7e825 1339 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
1340 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
a0b94e5c 1341
d7d7e825 1342 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1343 continue;
1344 }
26923b22 1345 if(negativeMC->GetPdgCode() == positiveMC->GetPdgCode()){
d7d7e825 1346 continue;
1347 }
a68437fb 1348
26923b22 1349 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){ // id 5 is conversion
a68437fb 1350 continue;
1351 }
a0b94e5c 1352
d7d7e825 1353 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a0b94e5c 1354
d7d7e825 1355 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1356 fHistograms->FillHistogram("ESD_NoCutConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1357 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1358 fHistograms->FillHistogram("ESD_NoCutConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1359 fHistograms->FillHistogram("ESD_NoCutConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1360 fHistograms->FillHistogram("ESD_NoCutConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1361 fHistograms->FillHistogram("ESD_NoCutConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1362 fHistograms->FillHistogram("ESD_NoCutConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1363 fHistograms->FillHistogram("ESD_NoCutConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1364 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
a0b94e5c 1365
d7d7e825 1366 fHistograms->FillHistogram("ESD_NoCutConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1367 fHistograms->FillHistogram("ESD_NoCutConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
1368
1369 fHistograms->FillHistogram("ESD_NoCutConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1370 fHistograms->FillHistogram("ESD_NoCutConversion_R", fV0Reader->GetXYRadius());
1371 fHistograms->FillHistogram("ESD_NoCutConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1372 fHistograms->FillHistogram("ESD_NoCutConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1373 fHistograms->FillHistogram("ESD_NoCutConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1374 fHistograms->FillHistogram("ESD_NoCutConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1375 fHistograms->FillHistogram("ESD_NoCutConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1376 fHistograms->FillHistogram("ESD_NoCutConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1377
1378 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1379 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1380 fHistograms->FillHistogram("ESD_NoCutConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1381 fHistograms->FillHistogram("ESD_NoCutConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1382
d7d7e825 1383 //store MCTruth properties
1384 fHistograms->FillHistogram("ESD_NoCutConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1385 fHistograms->FillHistogram("ESD_NoCutConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1386 fHistograms->FillHistogram("ESD_NoCutConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
d7d7e825 1387 }
1388 }
1389 }
1390 fV0Reader->ResetV0IndexNumber();
1391}
1392
1393void AliAnalysisTaskGammaConversion::ProcessV0s(){
1394 // see header file for documentation
037dc2db 1395
1396
d7d7e825 1397 if(fWriteNtuple == kTRUE){
1398 FillNtuple();
1399 }
1400
1401 Int_t nSurvivingV0s=0;
5ce758b0 1402 fV0Reader->ResetNGoodV0s();
d7d7e825 1403 while(fV0Reader->NextV0()){
1404 nSurvivingV0s++;
1405
e60f3265 1406
1407 TVector3 vtxConv(fV0Reader->GetX(),fV0Reader->GetY(), fV0Reader->GetZ());
1408
d7d7e825 1409 //-------------------------- filling v0 information -------------------------------------
1410 fHistograms->FillHistogram("ESD_Conversion_R", fV0Reader->GetXYRadius());
1411 fHistograms->FillHistogram("ESD_Conversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1412 fHistograms->FillHistogram("ESD_Conversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1413 fHistograms->FillHistogram("ESD_Conversion_OpeningAngle", fV0Reader->GetOpeningAngle());
e60f3265 1414
1415 // Specific histograms for beam pipe studies
1416 if( TMath::Abs(fV0Reader->GetZ()) < fV0Reader->GetLineCutZValue() ){
1417 fHistograms->FillHistogram("ESD_Conversion_XY_BeamPipe", fV0Reader->GetX(),fV0Reader->GetY());
1418 fHistograms->FillHistogram("ESD_Conversion_RPhi_BeamPipe", vtxConv.Phi(),fV0Reader->GetXYRadius());
1419 }
1420
1421
d7d7e825 1422 fHistograms->FillHistogram("ESD_E_Energy", fV0Reader->GetNegativeTrackEnergy());
1423 fHistograms->FillHistogram("ESD_E_Pt", fV0Reader->GetNegativeTrackPt());
1424 fHistograms->FillHistogram("ESD_E_Eta", fV0Reader->GetNegativeTrackEta());
1425 fHistograms->FillHistogram("ESD_E_Phi", fV0Reader->GetNegativeTrackPhi());
ebcfaa7e 1426 fHistograms->FillHistogram("ESD_E_nTPCClusters", fV0Reader->GetNegativeTracknTPCClusters());
1427 fHistograms->FillHistogram("ESD_E_nITSClusters", fV0Reader->GetNegativeTracknITSClusters());
5ce758b0 1428 if(fV0Reader->GetNegativeTracknTPCFClusters()!=0 && fV0Reader->GetNegativeTracknTPCClusters()!=0 ){
1429 Double_t EclsToF= (Double_t)fV0Reader->GetNegativeTracknTPCClusters()/(Double_t)fV0Reader->GetNegativeTracknTPCFClusters();
1430 fHistograms->FillHistogram("ESD_E_nTPCClustersToFP", fV0Reader->GetNegativeTrackP(),EclsToF );
1431 fHistograms->FillHistogram("ESD_E_TPCchi2", fV0Reader->GetNegativeTrackTPCchi2()/(Double_t)fV0Reader->GetNegativeTracknTPCClusters());
1432 }
1433
1434
1435
d7d7e825 1436 fHistograms->FillHistogram("ESD_P_Energy", fV0Reader->GetPositiveTrackEnergy());
1437 fHistograms->FillHistogram("ESD_P_Pt", fV0Reader->GetPositiveTrackPt());
1438 fHistograms->FillHistogram("ESD_P_Eta", fV0Reader->GetPositiveTrackEta());
1439 fHistograms->FillHistogram("ESD_P_Phi", fV0Reader->GetPositiveTrackPhi());
ebcfaa7e 1440 fHistograms->FillHistogram("ESD_P_nTPCClusters", fV0Reader->GetPositiveTracknTPCClusters());
1441 fHistograms->FillHistogram("ESD_P_nITSClusters", fV0Reader->GetPositiveTracknITSClusters());
5ce758b0 1442 if(fV0Reader->GetPositiveTracknTPCFClusters()!=0 && (Double_t)fV0Reader->GetPositiveTracknTPCClusters()!=0 ){
1443 Double_t PclsToF= (Double_t)fV0Reader->GetPositiveTracknTPCClusters()/(Double_t)fV0Reader->GetPositiveTracknTPCFClusters();
1444 fHistograms->FillHistogram("ESD_P_nTPCClustersToFP",fV0Reader->GetPositiveTrackP(), PclsToF);
1445 fHistograms->FillHistogram("ESD_P_TPCchi2", fV0Reader->GetPositiveTrackTPCchi2()/(Double_t)fV0Reader->GetPositiveTracknTPCClusters());
1446 }
1447
1448
d7d7e825 1449 fHistograms->FillHistogram("ESD_ConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1450 fHistograms->FillHistogram("ESD_ConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1451 fHistograms->FillHistogram("ESD_ConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1452 fHistograms->FillHistogram("ESD_ConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1453 fHistograms->FillHistogram("ESD_ConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1454 fHistograms->FillHistogram("ESD_ConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1455 fHistograms->FillHistogram("ESD_ConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1456 fHistograms->FillHistogram("ESD_ConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1457 fHistograms->FillHistogram("ESD_ConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1458 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1459
1460 fHistograms->FillHistogram("ESD_ConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1461 fHistograms->FillHistogram("ESD_ConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
9640a3d1 1462
1463 fHistograms->FillHistogram("ESD_ConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1464 fHistograms->FillHistogram("ESD_ConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1465 fHistograms->FillHistogram("ESD_ConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1466 fHistograms->FillHistogram("ESD_ConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
1467
1468 fHistograms->FillHistogram("ESD_ConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1469 fHistograms->FillHistogram("ESD_ConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1470 fHistograms->FillHistogram("ESD_ConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1471 fHistograms->FillHistogram("ESD_ConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
1472
9c1cb6f7 1473 Double_t negPID=0;
1474 Double_t posPID=0;
1475 fV0Reader->GetPIDProbability(negPID,posPID);
1476 fHistograms->FillHistogram("ESD_ConvGamma_E_EProbP",fV0Reader->GetNegativeTrackP(),negPID);
1477 fHistograms->FillHistogram("ESD_ConvGamma_P_EProbP",fV0Reader->GetPositiveTrackP(),posPID);
1478
1479 Double_t negPIDmupi=0;
1480 Double_t posPIDmupi=0;
1481 fV0Reader->GetPIDProbabilityMuonPion(negPIDmupi,posPIDmupi);
1482 fHistograms->FillHistogram("ESD_ConvGamma_E_mupiProbP",fV0Reader->GetNegativeTrackP(),negPIDmupi);
1483 fHistograms->FillHistogram("ESD_ConvGamma_P_mupiProbP",fV0Reader->GetPositiveTrackP(),posPIDmupi);
1484
48682642 1485 Double_t armenterosQtAlfa[2];
1486 fV0Reader->GetArmenterosQtAlfa(fV0Reader-> GetNegativeKFParticle(),
1487 fV0Reader-> GetPositiveKFParticle(),
1488 fV0Reader->GetMotherCandidateKFCombination(),
1489 armenterosQtAlfa);
1490
1491 fHistograms->FillHistogram("ESD_ConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
1492
037dc2db 1493
d7d7e825 1494 // begin mapping
1495 Int_t rBin = fHistograms->GetRBin(fV0Reader->GetXYRadius());
9640a3d1 1496 Int_t zBin = fHistograms->GetZBin(fV0Reader->GetZ());
d7d7e825 1497 Int_t phiBin = fHistograms->GetPhiBin(fV0Reader->GetNegativeTrackPhi());
9c1cb6f7 1498 Double_t rFMD=30;
1499
e60f3265 1500
9640a3d1 1501
ebcfaa7e 1502 // Double_t motherCandidateEta= fV0Reader->GetMotherCandidateEta();
d7d7e825 1503
1504 TString nameESDMappingPhiR="";
e158cbc3 1505 nameESDMappingPhiR.Form("ESD_Conversion_Mapping_Phi%02d_R%02d",phiBin,rBin);
ebcfaa7e 1506 //fHistograms->FillHistogram(nameESDMappingPhiR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1507
1508 TString nameESDMappingPhi="";
e158cbc3 1509 nameESDMappingPhi.Form("ESD_Conversion_Mapping_Phi%02d",phiBin);
ebcfaa7e 1510 //fHistograms->FillHistogram(nameESDMappingPhi, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1511
1512 TString nameESDMappingR="";
e158cbc3 1513 nameESDMappingR.Form("ESD_Conversion_Mapping_R%02d",rBin);
ebcfaa7e 1514 //fHistograms->FillHistogram(nameESDMappingR, fV0Reader->GetZ(), motherCandidateEta);
d7d7e825 1515
1516 TString nameESDMappingPhiInR="";
e158cbc3 1517 nameESDMappingPhiInR.Form("ESD_Conversion_Mapping_Phi_in_R_%02d",rBin);
9640a3d1 1518 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1519 fHistograms->FillHistogram(nameESDMappingPhiInR, vtxConv.Phi());
1520
1521 TString nameESDMappingZInR="";
1522 nameESDMappingZInR.Form("ESD_Conversion_Mapping_Z_in_R_%02d",rBin);
1523 fHistograms->FillHistogram(nameESDMappingZInR, fV0Reader->GetZ());
1524
1525 TString nameESDMappingPhiInZ="";
1526 nameESDMappingPhiInZ.Form("ESD_Conversion_Mapping_Phi_in_Z_%02d",zBin);
1527 // fHistograms->FillHistogram(nameESDMappingPhiInR, fV0Reader->GetMotherCandidatePhi());
1528 fHistograms->FillHistogram(nameESDMappingPhiInZ, vtxConv.Phi());
1529
9c1cb6f7 1530 if(fV0Reader->GetXYRadius()<rFMD){
1531 TString nameESDMappingFMDPhiInZ="";
1532 nameESDMappingFMDPhiInZ.Form("ESD_Conversion_Mapping_FMD_Phi_in_Z_%02d",zBin);
1533 fHistograms->FillHistogram(nameESDMappingFMDPhiInZ, vtxConv.Phi());
1534 }
1535
1536
9640a3d1 1537 TString nameESDMappingRInZ="";
1538 nameESDMappingRInZ.Form("ESD_Conversion_Mapping_R_in_Z_%02d",zBin);
1539 fHistograms->FillHistogram(nameESDMappingRInZ, fV0Reader->GetXYRadius());
1540
1541 if(fV0Reader->GetMotherCandidatePt() > fLowPtMapping && fV0Reader->GetMotherCandidatePt()< fHighPtMapping){
1542 TString nameESDMappingMidPtPhiInR="";
1543 nameESDMappingMidPtPhiInR.Form("ESD_Conversion_Mapping_MidPt_Phi_in_R_%02d",rBin);
1544 fHistograms->FillHistogram(nameESDMappingMidPtPhiInR, vtxConv.Phi());
1545
1546 TString nameESDMappingMidPtZInR="";
1547 nameESDMappingMidPtZInR.Form("ESD_Conversion_Mapping_MidPt_Z_in_R_%02d",rBin);
1548 fHistograms->FillHistogram(nameESDMappingMidPtZInR, fV0Reader->GetZ());
1549
1550 TString nameESDMappingMidPtPhiInZ="";
1551 nameESDMappingMidPtPhiInZ.Form("ESD_Conversion_Mapping_MidPt_Phi_in_Z_%02d",zBin);
1552 fHistograms->FillHistogram(nameESDMappingMidPtPhiInZ, vtxConv.Phi());
9c1cb6f7 1553 if(fV0Reader->GetXYRadius()<rFMD){
1554 TString nameESDMappingMidPtFMDPhiInZ="";
1555 nameESDMappingMidPtFMDPhiInZ.Form("ESD_Conversion_Mapping_MidPt_FMD_Phi_in_Z_%02d",zBin);
1556 fHistograms->FillHistogram(nameESDMappingMidPtFMDPhiInZ, vtxConv.Phi());
1557 }
1558
9640a3d1 1559
1560 TString nameESDMappingMidPtRInZ="";
1561 nameESDMappingMidPtRInZ.Form("ESD_Conversion_Mapping_MidPt_R_in_Z_%02d",zBin);
1562 fHistograms->FillHistogram(nameESDMappingMidPtRInZ, fV0Reader->GetXYRadius());
1563
1564 }
1565
1566
d7d7e825 1567 // end mapping
1568
6c84d371 1569 new((*fKFReconstructedGammasTClone)[fKFReconstructedGammasTClone->GetEntriesFast()]) AliKFParticle(*fV0Reader->GetMotherCandidateKFCombination());
037dc2db 1570 fKFReconstructedGammasV0Index.push_back(fV0Reader->GetCurrentV0IndexNumber()-1);
6c84d371 1571 // fKFReconstructedGammas.push_back(*fV0Reader->GetMotherCandidateKFCombination());
d7d7e825 1572 fElectronv1.push_back(fV0Reader->GetCurrentV0()->GetPindex());
1573 fElectronv2.push_back(fV0Reader->GetCurrentV0()->GetNindex());
a0b94e5c 1574
d7d7e825 1575 //----------------------------------- checking for "real" conversions (MC match) --------------------------------------
1576 if(fDoMCTruth){
1577
a68437fb 1578 if(fV0Reader->HasSameMCMother() == kFALSE){
1579 continue;
1580 }
d7d7e825 1581 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
a68437fb 1582 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
1583
1584 if(TMath::Abs(negativeMC->GetPdgCode())!=11 || TMath::Abs(positiveMC->GetPdgCode())!=11){
1585 continue;
1586 }
1587
1588 if(negativeMC->GetPdgCode()==positiveMC->GetPdgCode()){
1589 continue;
1590 }
bd6d9fa3 1591 if( (negativeMC->GetUniqueID() == 4 && positiveMC->GetUniqueID() ==4) ||
1592 (negativeMC->GetUniqueID() == 0 && positiveMC->GetUniqueID() ==0) ){// fill r distribution for Dalitz decays
037dc2db 1593 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 111){ //pi0
1594 fHistograms->FillHistogram("ESD_TrueDalitzContamination_R", fV0Reader->GetXYRadius());
1595 }
1596 }
a68437fb 1597
1598 if(negativeMC->GetUniqueID() != 5 || positiveMC->GetUniqueID() !=5){// check if the daughters come from a conversion
1599 continue;
1600 }
a0b94e5c 1601
d7d7e825 1602 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
a68437fb 1603
26923b22 1604 if(fDoCF){
1605 Double_t containerInput[3];
1606 containerInput[0] = fV0Reader->GetMotherCandidatePt();
1607 containerInput[1] = fV0Reader->GetMotherCandidateEta();
1608 containerInput[2] = fV0Reader->GetMotherCandidateMass();
1609 fCFManager->GetParticleContainer()->Fill(containerInput,kStepTrueGamma); // for CF
1610 }
1611
d7d7e825 1612 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt", fV0Reader->GetMotherCandidatePt());
1613 fHistograms->FillHistogram("ESD_TrueConvGamma_Energy", fV0Reader->GetMotherCandidateEnergy());
1614 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta", fV0Reader->GetMotherCandidateEta());
1615 fHistograms->FillHistogram("ESD_TrueConvGamma_Phi", fV0Reader->GetMotherCandidatePhi());
1616 fHistograms->FillHistogram("ESD_TrueConvGamma_Mass", fV0Reader->GetMotherCandidateMass());
1617 fHistograms->FillHistogram("ESD_TrueConvGamma_Width", fV0Reader->GetMotherCandidateWidth());
1618 fHistograms->FillHistogram("ESD_TrueConvGamma_Chi2", fV0Reader->GetMotherCandidateChi2());
1619 fHistograms->FillHistogram("ESD_TrueConvGamma_NDF", fV0Reader->GetMotherCandidateNDF());
1620 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Eta", fV0Reader->GetMotherCandidatePt(),fV0Reader->GetMotherCandidateEta());
1621 fHistograms->FillHistogram("ESD_TrueConvGamma_Rapid", fV0Reader->GetMotherCandidateRapidity());
1622 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters());
1623 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLength", /*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters());
1624 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetNegativeTrackLength()*/fV0Reader->GetNegativeNTPCClusters(),fV0Reader->GetMotherCandidateMass());
1625 fHistograms->FillHistogram("ESD_TrueConvGamma_TrackLengthVSInvMass",/*fV0Reader->GetPositiveTrackLength()*/fV0Reader->GetPositiveNTPCClusters(),fV0Reader->GetMotherCandidateMass());
a0b94e5c 1626
d7d7e825 1627 fHistograms->FillHistogram("ESD_TrueConvGamma_Pt_Chi2", fV0Reader->GetMotherCandidatePt(), fV0Reader->GetMotherCandidateChi2());
1628 fHistograms->FillHistogram("ESD_TrueConvGamma_Eta_Chi2", fV0Reader->GetMotherCandidateEta(), fV0Reader->GetMotherCandidateChi2());
a0b94e5c 1629
d7d7e825 1630
1631 fHistograms->FillHistogram("ESD_TrueConversion_XY", fV0Reader->GetX(),fV0Reader->GetY());
1632 fHistograms->FillHistogram("ESD_TrueConversion_R", fV0Reader->GetXYRadius());
1633 fHistograms->FillHistogram("ESD_TrueConversion_ZR", fV0Reader->GetZ(),fV0Reader->GetXYRadius());
1634 fHistograms->FillHistogram("ESD_TrueConversion_OpeningAngle", fV0Reader->GetOpeningAngle());
9640a3d1 1635
1636 fHistograms->FillHistogram("ESD_TrueConvGamma_CosPointingAngle", fV0Reader->GetCosPointingAngle());
1637 fHistograms->FillHistogram("ESD_TrueConvGamma_DcaDaughters", fV0Reader->GetDcaDaughters());
1638 fHistograms->FillHistogram("ESD_TrueConvGamma_NormDcaDistDaughters", fV0Reader->GetNormDcaDistDaughters());
1639 fHistograms->FillHistogram("ESD_TrueConvGamma_LikelihoodAP", fV0Reader->GetLikelihoodAP());
96ade8ef 1640 if (fV0Reader->GetMotherCandidateP() != 0) {
1641 fHistograms->FillHistogram("ESD_TrueConvGamma_E_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetNegativeTrackP()/fV0Reader->GetMotherCandidateP());
1642 fHistograms->FillHistogram("ESD_TrueConvGamma_P_AsymmetryP",fV0Reader->GetMotherCandidateP(),fV0Reader->GetPositiveTrackP()/fV0Reader->GetMotherCandidateP());
1643 } else { cout << "Error::fV0Reader->GetNegativeTrackP() == 0 !!!" << endl; }
1644
9640a3d1 1645 fHistograms->FillHistogram("ESD_TrueConvGamma_E_dEdxP",fV0Reader->GetNegativeTrackP(),fV0Reader->GetNegativeTrackTPCdEdx());
1646 fHistograms->FillHistogram("ESD_TrueConvGamma_P_dEdxP",fV0Reader->GetPositiveTrackP(),fV0Reader->GetPositiveTrackTPCdEdx());
70ef88b5 1647 fHistograms->FillHistogram("ESD_TrueConvGamma_alfa_qt",armenterosQtAlfa[1],armenterosQtAlfa[0]);
1648
9640a3d1 1649
a0b94e5c 1650
d7d7e825 1651 //store MCTruth properties
1652 fHistograms->FillHistogram("ESD_TrueConvGamma_MC_Pt_Eta", fV0Reader->GetMotherMCParticle()->Pt(),fV0Reader->GetMotherMCParticle()->Eta());
1653 fHistograms->FillHistogram("ESD_TrueConversion_MC_ZR", negativeMC->Vz(),negativeMC->R());
1654 fHistograms->FillHistogram("ESD_TrueConversion_MC_XY", negativeMC->Vx(),negativeMC->Vy());
48682642 1655
d7d7e825 1656 //resolution
1657 Double_t mcpt = fV0Reader->GetMotherMCParticle()->Pt();
1658 Double_t esdpt = fV0Reader->GetMotherCandidatePt();
9c1cb6f7 1659 Double_t resdPt = 0.;
4a6157dc 1660 if(mcpt > 0){
9c1cb6f7 1661 resdPt = ((esdpt - mcpt)/mcpt)*100.;
d7d7e825 1662 }
4a6157dc 1663 else if(mcpt < 0){
1664 cout<<"Pt of MC particle is negative, this will cause wrong calculation of resPt"<<endl;
1665 }
d7d7e825 1666
ca6d4600 1667 fHistograms->FillHistogram("Resolution_Gamma_dPt_Pt", mcpt, resdPt);
d7d7e825 1668 fHistograms->FillHistogram("Resolution_MC_Pt", mcpt);
1669 fHistograms->FillHistogram("Resolution_ESD_Pt", esdpt);
ca6d4600 1670 fHistograms->FillHistogram("Resolution_Gamma_dPt_Phi", fV0Reader->GetMotherCandidatePhi(), resdPt);
d7d7e825 1671
9c1cb6f7 1672 Double_t resdZ = 0.;
d7d7e825 1673 if(fV0Reader->GetNegativeMCParticle()->Vz() != 0){
9c1cb6f7 1674 resdZ = ((fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz())/fV0Reader->GetNegativeMCParticle()->Vz())*100.;
d7d7e825 1675 }
9c1cb6f7 1676 Double_t resdZAbs = 0.;
037dc2db 1677 resdZAbs = (fV0Reader->GetZ() -fV0Reader->GetNegativeMCParticle()->Vz());
1678
1679 fHistograms->FillHistogram("Resolution_dZAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdZAbs);
d7d7e825 1680 fHistograms->FillHistogram("Resolution_dZ", fV0Reader->GetNegativeMCParticle()->Vz(), resdZ);
1681 fHistograms->FillHistogram("Resolution_MC_Z", fV0Reader->GetNegativeMCParticle()->Vz());
1682 fHistograms->FillHistogram("Resolution_ESD_Z", fV0Reader->GetZ());
ca6d4600 1683
1684 // new for dPt_Pt-histograms for Electron and Positron
1685 Double_t mcEpt = fV0Reader->GetNegativeMCParticle()->Pt();
9c1cb6f7 1686 Double_t resEdPt = 0.;
ca6d4600 1687 if (mcEpt > 0){
9c1cb6f7 1688 resEdPt = ((fV0Reader->GetNegativeTrackPt()-mcEpt)/mcEpt)*100.;
ca6d4600 1689 }
1690 UInt_t statusN = fV0Reader->GetNegativeESDTrack()->GetStatus();
1691 // AliESDtrack * negTrk = fV0Reader->GetNegativeESDTrack();
1692 UInt_t kTRDoutN = (statusN & AliESDtrack::kTRDout);
1693
1694 Int_t ITSclsE= fV0Reader->GetNegativeTracknITSClusters();
1695 // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
1696 switch(ITSclsE){
1697 case 0: // 0 ITS clusters
1698 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS0", mcEpt, resEdPt);
1699 break;
1700 case 1: // 1 ITS cluster
1701 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS1", mcEpt, resEdPt);
1702 break;
1703 case 2: // 2 ITS clusters
1704 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS2", mcEpt, resEdPt);
1705 break;
1706 case 3: // 3 ITS clusters
1707 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS3", mcEpt, resEdPt);
1708 break;
1709 case 4: // 4 ITS clusters
1710 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS4", mcEpt, resEdPt);
1711 break;
1712 case 5: // 5 ITS clusters
1713 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS5", mcEpt, resEdPt);
1714 break;
1715 case 6: // 6 ITS clusters
1716 fHistograms->FillHistogram("Resolution_E_dPt_Pt_ITS6", mcEpt, resEdPt);
1717 break;
1718 }
1719 //Filling histograms with respect to Electron resolution
1720 fHistograms->FillHistogram("Resolution_E_dPt_Pt", mcEpt, resEdPt);
1721 fHistograms->FillHistogram("Resolution_E_dPt_Phi", fV0Reader->GetNegativeTrackPhi(), resEdPt);
1722 if(kTRDoutN){
1723 fHistograms->FillHistogram("Resolution_E_nTRDtracklets_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
1724 fHistograms->FillHistogram("Resolution_E_nTRDtracklets_MCPt", mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDntracklets());
1725 fHistograms->FillHistogram("Resolution_E_nTRDclusters_ESDPt",fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDncls());
1726 fHistograms->FillHistogram("Resolution_E_nTRDclusters_MCPt",mcEpt, fV0Reader->GetNegativeESDTrack()->GetTRDncls());
1727 fHistograms->FillHistogram("Resolution_E_TRDsignal_ESDPt", fV0Reader->GetNegativeTrackPt(), fV0Reader->GetNegativeESDTrack()->GetTRDsignal());
1728 }
1729
1730 Double_t mcPpt = fV0Reader->GetPositiveMCParticle()->Pt();
1731 Double_t resPdPt = 0;
1732 if (mcPpt > 0){
9c1cb6f7 1733 resPdPt = ((fV0Reader->GetPositiveTrackPt()-mcPpt)/mcPpt)*100.;
ca6d4600 1734 }
1735
1736 UInt_t statusP = fV0Reader->GetPositiveESDTrack()->GetStatus();
1737// AliESDtrack * posTr= fV0Reader->GetPositiveESDTrack();
1738 UInt_t kTRDoutP = (statusP & AliESDtrack::kTRDout);
1739
1740 Int_t ITSclsP = fV0Reader->GetPositiveTracknITSClusters();
1741 // filling Resolution_Pt_dPt with respect to the Number of ITS clusters for Positrons
1742 switch(ITSclsP){
1743 case 0: // 0 ITS clusters
1744 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS0", mcPpt, resPdPt);
1745 break;
1746 case 1: // 1 ITS cluster
1747 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS1", mcPpt, resPdPt);
1748 break;
1749 case 2: // 2 ITS clusters
1750 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS2", mcPpt, resPdPt);
1751 break;
1752 case 3: // 3 ITS clusters
1753 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS3", mcPpt, resPdPt);
1754 break;
1755 case 4: // 4 ITS clusters
1756 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS4", mcPpt, resPdPt);
1757 break;
1758 case 5: // 5 ITS clusters
1759 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS5", mcPpt, resPdPt);
1760 break;
1761 case 6: // 6 ITS clusters
1762 fHistograms->FillHistogram("Resolution_P_dPt_Pt_ITS6", mcPpt, resPdPt);
1763 break;
1764 }
1765 //Filling histograms with respect to Positron resolution
1766 fHistograms->FillHistogram("Resolution_P_dPt_Pt", mcPpt, resPdPt);
1767 fHistograms->FillHistogram("Resolution_P_dPt_Phi", fV0Reader->GetPositiveTrackPhi(), resPdPt);
1768 if(kTRDoutP){
1769 fHistograms->FillHistogram("Resolution_P_nTRDtracklets_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
1770 fHistograms->FillHistogram("Resolution_P_nTRDtracklets_MCPt", mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDntracklets());
1771 fHistograms->FillHistogram("Resolution_P_nTRDclusters_ESDPt",fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDncls());
1772 fHistograms->FillHistogram("Resolution_P_nTRDclusters_MCPt",mcPpt, fV0Reader->GetPositiveESDTrack()->GetTRDncls());
1773 fHistograms->FillHistogram("Resolution_P_TRDsignal_ESDPt", fV0Reader->GetPositiveTrackPt(), fV0Reader->GetPositiveESDTrack()->GetTRDsignal());
1774 }
1775
1776
9c1cb6f7 1777 Double_t resdR = 0.;
d7d7e825 1778 if(fV0Reader->GetNegativeMCParticle()->R() != 0){
9c1cb6f7 1779 resdR = ((fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R())/fV0Reader->GetNegativeMCParticle()->R())*100.;
d7d7e825 1780 }
9c1cb6f7 1781 Double_t resdRAbs = 0.;
037dc2db 1782 resdRAbs = (fV0Reader->GetXYRadius() - fV0Reader->GetNegativeMCParticle()->R());
1783
1784 fHistograms->FillHistogram("Resolution_dRAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdRAbs);
d7d7e825 1785 fHistograms->FillHistogram("Resolution_dR", fV0Reader->GetNegativeMCParticle()->R(), resdR);
1786 fHistograms->FillHistogram("Resolution_MC_R", fV0Reader->GetNegativeMCParticle()->R());
1787 fHistograms->FillHistogram("Resolution_ESD_R", fV0Reader->GetXYRadius());
ca6d4600 1788 fHistograms->FillHistogram("Resolution_R_dPt", fV0Reader->GetNegativeMCParticle()->R(), resdPt);
1789
9c1cb6f7 1790 Double_t resdPhiAbs=0.;
1791 resdPhiAbs=0.;
037dc2db 1792 resdPhiAbs= (fV0Reader->GetMotherCandidatePhi()-fV0Reader->GetNegativeMCParticle()->Phi());
1793 fHistograms->FillHistogram("Resolution_dPhiAbs_VS_R", fV0Reader->GetNegativeMCParticle()->R(), resdPhiAbs);
1794
d7d7e825 1795 }//if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22)
d7d7e825 1796 }//if(fDoMCTruth)
1797 }//while(fV0Reader->NextV0)
1798 fHistograms->FillHistogram("ESD_NumberOfSurvivingV0s", nSurvivingV0s);
1799 fHistograms->FillHistogram("ESD_NumberOfV0s", fV0Reader->GetNumberOfV0s());
b5832f95 1800 fHistograms->FillHistogram("ESD_NumberOfContributorsVtx", fV0Reader->GetNumberOfContributorsVtx());
cb90a330 1801
1802 fV0Reader->ResetV0IndexNumber();
d7d7e825 1803}
1804
1805void AliAnalysisTaskGammaConversion::FillAODWithConversionGammas(){
1806 // Fill AOD with reconstructed Gamma
a0b94e5c 1807
6c84d371 1808 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
1809 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
d7d7e825 1810 //Create AOD particle object from AliKFParticle
d7d7e825 1811 //You could add directly AliKFParticle objects to the AOD, avoiding dependences with PartCorr
1812 //but this means that I have to work a little bit more in my side.
1813 //AODPWG4Particle objects are simpler and lighter, I think
332f1f44 1814 /*
1815 AliKFParticle * gammakf = dynamic_cast<AliKFParticle*>(fKFReconstructedGammasTClone->At(gammaIndex));
d7d7e825 1816 AliAODPWG4Particle gamma = AliAODPWG4Particle(gammakf->Px(),gammakf->Py(),gammakf->Pz(), gammakf->E());
332f1f44 1817 //gamma.SetLabel(-1);//How to get the MC label of the reconstructed gamma?
1818 gamma.SetTrackLabel( fElectronv1[gammaIndex], fElectronv2[gammaIndex] ); //How to get the MC label of the 2 electrons that form the gamma?
a0b94e5c 1819 gamma.SetDetector("CTS"); //tag the gamma as reconstructed in the central barrel
332f1f44 1820 gamma.SetPdg(AliPID::kEleCon); //photon id
a0b94e5c 1821 gamma.SetTag(-1); //Here I usually put a flag saying that montecarlo says it is prompt, decay fragmentation photon, or hadrons or whatever
332f1f44 1822 gamma.SetChi2(gammakf->Chi2());
a0b94e5c 1823 Int_t i = fAODBranch->GetEntriesFast();
1824 new((*fAODBranch)[i]) AliAODPWG4Particle(gamma);
d7d7e825 1825 */
332f1f44 1826
6c84d371 1827 AliKFParticle * gammakf = (AliKFParticle *)fKFReconstructedGammasTClone->At(gammaIndex);
d7d7e825 1828 AliGammaConversionAODObject aodObject;
1829 aodObject.SetPx(gammakf->GetPx());
1830 aodObject.SetPy(gammakf->GetPy());
1831 aodObject.SetPz(gammakf->GetPz());
1832 aodObject.SetLabel1(fElectronv1[gammaIndex]);
1833 aodObject.SetLabel2(fElectronv2[gammaIndex]);
332f1f44 1834 aodObject.SetChi2(gammakf->Chi2());
04bf4381 1835 aodObject.SetE(gammakf->E());
1836 Int_t i = fAODGamma->GetEntriesFast();
1837 new((*fAODGamma)[i]) AliGammaConversionAODObject(aodObject);
a0b94e5c 1838 }
1839
d7d7e825 1840}
1841
6272370b 1842void AliAnalysisTaskGammaConversion::ProcessGammasForOmegaMesonAnalysis(){
1843 // omega meson analysis pi0+gamma decay
1844 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
48682642 1845 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
6272370b 1846 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
48682642 1847
6272370b 1848 AliKFParticle * omegaCandidateGammaDaughter = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
1849 if(fGammav1[firstPi0Index]==firstGammaIndex || fGammav2[firstPi0Index]==firstGammaIndex){
1850 continue;
1851 }
1852
48682642 1853 AliKFParticle omegaCandidate(*omegaCandidatePi0Daughter,*omegaCandidateGammaDaughter);
6272370b 1854 Double_t massOmegaCandidate = 0.;
1855 Double_t widthOmegaCandidate = 0.;
1856
48682642 1857 omegaCandidate.GetMass(massOmegaCandidate,widthOmegaCandidate);
6272370b 1858
e9aea39f 1859 if ( massOmegaCandidate > 733 && massOmegaCandidate < 833 ) {
c59360eb 1860 AddOmegaToAOD(&omegaCandidate, massOmegaCandidate, firstPi0Index, firstGammaIndex);
1861 }
1862
48682642 1863 fHistograms->FillHistogram("ESD_Omega_InvMass_vs_Pt",massOmegaCandidate ,omegaCandidate.GetPt());
6272370b 1864 fHistograms->FillHistogram("ESD_Omega_InvMass",massOmegaCandidate);
48682642 1865
1866 //delete omegaCandidate;
1867
1868 }// end of omega reconstruction in pi0+gamma channel
1869
1870 if(fDoJet == kTRUE){
1871 AliKFParticle* negPiKF=NULL;
1872 AliKFParticle* posPiKF=NULL;
1873
1874 // look at the pi+pi+pi0 channel
1875 for(Int_t iCh=0;iCh<fChargedParticles->GetEntriesFast();iCh++){
1876 AliESDtrack* posTrack = (AliESDtrack*)(fChargedParticles->At(iCh));
1877 if (posTrack->GetSign()<0) continue;
c0d9051e 1878 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(posTrack,AliPID::kPion))>2.) continue;
48682642 1879 if (posPiKF) delete posPiKF; posPiKF=NULL;
1880 posPiKF = new AliKFParticle( *(posTrack) ,211);
1881
1882 for(Int_t jCh=0;jCh<fChargedParticles->GetEntriesFast();jCh++){
1883 AliESDtrack* negTrack = (AliESDtrack*)(fChargedParticles->At(jCh));
1884 if( negTrack->GetSign()>0) continue;
c0d9051e 1885 if(TMath::Abs(fV0Reader->GetESDpid()->NumberOfSigmasTPC(negTrack,AliPID::kPion))>2.) continue;
48682642 1886 if (negPiKF) delete negPiKF; negPiKF=NULL;
1887 negPiKF = new AliKFParticle( *(negTrack) ,-211);
1888 AliKFParticle omegaCandidatePipPinPi0(*omegaCandidatePi0Daughter,*posPiKF,*negPiKF);
1889 Double_t massOmegaCandidatePipPinPi0 = 0.;
1890 Double_t widthOmegaCandidatePipPinPi0 = 0.;
1891
1892 omegaCandidatePipPinPi0.GetMass(massOmegaCandidatePipPinPi0,widthOmegaCandidatePipPinPi0);
ef2e2748 1893
1894 if ( massOmegaCandidatePipPinPi0 > 733 && massOmegaCandidatePipPinPi0 < 833 ) {
1895 AddOmegaToAOD(&omegaCandidatePipPinPi0, massOmegaCandidatePipPinPi0, -1, -1);
1896 }
48682642 1897
1898 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass_vs_Pt",massOmegaCandidatePipPinPi0 ,omegaCandidatePipPinPi0.GetPt());
1899 fHistograms->FillHistogram("ESD_OmegaPipPinPi0_InvMass",massOmegaCandidatePipPinPi0);
1900
1901 // delete omegaCandidatePipPinPi0;
1902 }
1903 }
1904 } // checking ig gammajet because in that case the chargedparticle list is created
1905
1906
6272370b 1907
6272370b 1908 }
48682642 1909
1910 if(fCalculateBackground){
5ce758b0 1911
1912 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
1913
1914 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
1915 Int_t mbin = 0;
1916 if(fUseTrackMultiplicityForBG == kTRUE){
1917 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
1918 }
1919 else{
1920 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
1921 }
1922
1923 AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
1924
48682642 1925 // Background calculation for the omega
1926 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
5ce758b0 1927 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);
1928
1929 if(fMoveParticleAccordingToVertex == kTRUE){
1930 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
1931 }
48682642 1932 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
1933 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
5ce758b0 1934
1935 if(fMoveParticleAccordingToVertex == kTRUE){
1936 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
1937 }
1938
48682642 1939 for(Int_t firstPi0Index=0;firstPi0Index<fKFReconstructedPi0sTClone->GetEntriesFast();firstPi0Index++){
1940 AliKFParticle * omegaCandidatePi0Daughter = (AliKFParticle *)fKFReconstructedPi0sTClone->At(firstPi0Index);
1941 AliKFParticle * omegaBckCandidate = new AliKFParticle(*omegaCandidatePi0Daughter,previousGoodV0);
1942 Double_t massOmegaBckCandidate = 0.;
1943 Double_t widthOmegaBckCandidate = 0.;
1944
1945 omegaBckCandidate->GetMass(massOmegaBckCandidate,widthOmegaBckCandidate);
ef2e2748 1946
1947
48682642 1948 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass_vs_Pt",massOmegaBckCandidate ,omegaBckCandidate->GetPt());
1949 fHistograms->FillHistogram("ESD_Omega_Bck_InvMass",massOmegaBckCandidate);
1950
1951 delete omegaBckCandidate;
1952 }
1953 }
1954 }
1955 } // end of checking if background calculation is available
6272370b 1956}
d7d7e825 1957
04bf4381 1958
1959void AliAnalysisTaskGammaConversion::AddOmegaToAOD(AliKFParticle * omegakf, Double_t mass, Int_t omegaDaughter, Int_t gammaDaughter) {
1960 //See header file for documentation
1961 AliGammaConversionAODObject omega;
1962 omega.SetPx(omegakf->Px());
1963 omega.SetPy(omegakf->Py());
1964 omega.SetPz(omegakf->Pz());
1965 omega.SetChi2(omegakf->GetChi2());
1966 omega.SetE(omegakf->E());
1967 omega.SetIMass(mass);
1968 omega.SetLabel1(omegaDaughter);
1969 //dynamic_cast<AliGammaConversionAODObject*>(fAODBranch->At(daughter1))->SetTagged(kTRUE);
1970 omega.SetLabel2(gammaDaughter);
1971 new((*fAODOmega)[fAODOmega->GetEntriesFast()]) AliGammaConversionAODObject(omega);
1972}
1973
1974
1975
d7d7e825 1976void AliAnalysisTaskGammaConversion::ProcessGammasForNeutralMesonAnalysis(){
1977 // see header file for documentation
1978
6c84d371 1979 // for(UInt_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammas.size();firstGammaIndex++){
1980 // for(UInt_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammas.size();secondGammaIndex++){
037dc2db 1981
5ce758b0 1982 fESDEvent = fV0Reader->GetESDEvent();
1983
037dc2db 1984 if(fKFReconstructedGammasTClone->GetEntriesFast()>fV0Reader->GetNumberOfV0s()){
1985 cout<<"Warning, number of entries in the tclone is bigger than number of v0s"<<endl;
1986 }
1987
6c84d371 1988 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
1989 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();secondGammaIndex++){
d7d7e825 1990
6c84d371 1991 // AliKFParticle * twoGammaDecayCandidateDaughter0 = &fKFReconstructedGammas[firstGammaIndex];
1992 // AliKFParticle * twoGammaDecayCandidateDaughter1 = &fKFReconstructedGammas[secondGammaIndex];
a0b94e5c 1993
6c84d371 1994 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
1995 AliKFParticle * twoGammaDecayCandidateDaughter1 = (AliKFParticle *)fKFReconstructedGammasTClone->At(secondGammaIndex);
a0b94e5c 1996
d7d7e825 1997 if(fElectronv1[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv1[firstGammaIndex]==fElectronv2[secondGammaIndex]){
1998 continue;
1999 }
2000 if(fElectronv2[firstGammaIndex]==fElectronv1[secondGammaIndex] || fElectronv2[firstGammaIndex]==fElectronv2[secondGammaIndex]){
2001 continue;
2002 }
a0b94e5c 2003
d7d7e825 2004 AliKFParticle *twoGammaCandidate = new AliKFParticle(*twoGammaDecayCandidateDaughter0,*twoGammaDecayCandidateDaughter1);
2005
2006 Double_t massTwoGammaCandidate = 0.;
2007 Double_t widthTwoGammaCandidate = 0.;
2008 Double_t chi2TwoGammaCandidate =10000.;
2009 twoGammaCandidate->GetMass(massTwoGammaCandidate,widthTwoGammaCandidate);
d707e3cf 2010 // if(twoGammaCandidate->GetNDF()>0){
2011 // chi2TwoGammaCandidate = twoGammaCandidate->GetChi2()/twoGammaCandidate->GetNDF();
2012 chi2TwoGammaCandidate = twoGammaCandidate->GetChi2();
d7d7e825 2013
70ef88b5 2014 fHistograms->FillHistogram("ESD_Mother_Chi2",chi2TwoGammaCandidate);
2e2da371 2015 if((chi2TwoGammaCandidate>0 && chi2TwoGammaCandidate<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
d7d7e825 2016
2017 TVector3 momentumVectorTwoGammaCandidate(twoGammaCandidate->GetPx(),twoGammaCandidate->GetPy(),twoGammaCandidate->GetPz());
2018 TVector3 spaceVectorTwoGammaCandidate(twoGammaCandidate->GetX(),twoGammaCandidate->GetY(),twoGammaCandidate->GetZ());
2019
2020 Double_t openingAngleTwoGammaCandidate = twoGammaDecayCandidateDaughter0->GetAngle(*twoGammaDecayCandidateDaughter1);
2021 Double_t rapidity;
96ade8ef 2022 if(twoGammaCandidate->GetE() - twoGammaCandidate->GetPz() <= 0 || twoGammaCandidate->GetE() + twoGammaCandidate->GetPz() <= 0){
2023 cout << "Error: |Pz| > E !!!! " << endl;
d7d7e825 2024 rapidity=0;
2025 }
2026 else{
2027 rapidity = 0.5*(TMath::Log((twoGammaCandidate->GetE() +twoGammaCandidate->GetPz()) / (twoGammaCandidate->GetE()-twoGammaCandidate->GetPz())));
2028 }
dc2883e4 2029
2030 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut()){
2031 delete twoGammaCandidate;
2032 continue; // rapidity cut
2033 }
2034
2035
48682642 2036 Double_t alfa=0.0;
2037 if( (twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()) != 0){
2038 alfa=TMath::Abs((twoGammaDecayCandidateDaughter0->GetE()-twoGammaDecayCandidateDaughter1->GetE())
2039 /(twoGammaDecayCandidateDaughter0->GetE()+twoGammaDecayCandidateDaughter1->GetE()));
2040 }
d7d7e825 2041
1e7846f4 2042 if(openingAngleTwoGammaCandidate < fMinOpeningAngleGhostCut){
2043 delete twoGammaCandidate;
2044 continue; // minimum opening angle to avoid using ghosttracks
2045 }
2046
d7d7e825 2047 fHistograms->FillHistogram("ESD_Mother_GammaDaughter_OpeningAngle", openingAngleTwoGammaCandidate);
2048 fHistograms->FillHistogram("ESD_Mother_Energy", twoGammaCandidate->GetE());
2049 fHistograms->FillHistogram("ESD_Mother_Pt", momentumVectorTwoGammaCandidate.Pt());
2050 fHistograms->FillHistogram("ESD_Mother_Eta", momentumVectorTwoGammaCandidate.Eta());
2051 fHistograms->FillHistogram("ESD_Mother_Rapid", rapidity);
2052 fHistograms->FillHistogram("ESD_Mother_Phi", spaceVectorTwoGammaCandidate.Phi());
2053 fHistograms->FillHistogram("ESD_Mother_Mass", massTwoGammaCandidate);
48682642 2054 fHistograms->FillHistogram("ESD_Mother_alfa", alfa);
d7d7e825 2055 fHistograms->FillHistogram("ESD_Mother_R", spaceVectorTwoGammaCandidate.Pt()); // Pt in Space == R!!!
2056 fHistograms->FillHistogram("ESD_Mother_ZR", twoGammaCandidate->GetZ(), spaceVectorTwoGammaCandidate.Pt());
2057 fHistograms->FillHistogram("ESD_Mother_XY", twoGammaCandidate->GetX(), twoGammaCandidate->GetY());
2058 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
67381a40 2059 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
9c1cb6f7 2060 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2061 }
3c45d101 2062 if(alfa<0.1){
2063 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_E_alpha",massTwoGammaCandidate ,twoGammaCandidate->GetE());
2064 }
d7d7e825 2065 fHistograms->FillHistogram("ESD_Mother_InvMass",massTwoGammaCandidate);
9640a3d1 2066
9c1cb6f7 2067 if(fCalculateBackground){
5ce758b0 2068 /* Kenneth, just for testing*/
2069 AliGammaConversionBGHandler * bgHandlerTest = fV0Reader->GetBGHandler();
2070
2071 Int_t zbin= bgHandlerTest->GetZBinIndex(fV0Reader->GetVertexZ());
2072 Int_t mbin=0;
2073 Int_t multKAA=0;
2074 if(fUseTrackMultiplicityForBG == kTRUE){
2075 multKAA=fV0Reader->CountESDTracks();
2076 mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2077 }
2078 else{// means we use #v0s for multiplicity
2079 multKAA=fV0Reader->GetNGoodV0s();
2080 mbin = bgHandlerTest->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2081 }
2082 // cout<<"Filling bin number "<<zbin<<" and "<<mbin<<endl;
2083 // cout<<"Corresponding to z = "<<fV0Reader->GetVertexZ()<<" and m = "<<multKAA<<endl;
2084 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass",zbin,mbin),massTwoGammaCandidate);
2085 fHistograms->FillHistogram(Form("%d%dESD_Mother_InvMass_vs_Pt",zbin,mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2086 /* end Kenneth, just for testing*/
2087 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2088 }
2089 /* if(fCalculateBackground){
9c1cb6f7 2090 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2091 Int_t mbin= bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2092 fHistograms->FillHistogram(Form("%dESD_Mother_InvMass_vs_Pt",mbin),massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
5ce758b0 2093 }*/
ebcfaa7e 2094 // if(fDoNeutralMesonV0MCCheck){
2095 if(fDoMCTruth){
037dc2db 2096 //Kenneth: Checking the eta of the gamma to check the difference between 0.9 and 1.2
2097 Int_t indexKF1 = fKFReconstructedGammasV0Index.at(firstGammaIndex);
2098 if(indexKF1<fV0Reader->GetNumberOfV0s()){
2099 fV0Reader->GetV0(indexKF1);//updates to the correct v0
2100 Double_t eta1 = fV0Reader->GetMotherCandidateEta();
2101 Bool_t isRealPi0=kFALSE;
10e3319b 2102 Bool_t isRealEta=kFALSE;
037dc2db 2103 Int_t gamma1MotherLabel=-1;
2104 if(fV0Reader->HasSameMCMother() == kTRUE){
2105 //cout<<"This v0 is a real v0!!!!"<<endl;
2106 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2107 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
2108 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
2109 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
2110 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2111 gamma1MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
2112 }
2113 }
2114 }
2115 }
2116 Int_t indexKF2 = fKFReconstructedGammasV0Index.at(secondGammaIndex);
2117 if(indexKF1 == indexKF2){
2118 cout<<"index of the two KF particles are the same.... should not happen"<<endl;
2119 }
2120 if(indexKF2<fV0Reader->GetNumberOfV0s()){
2121 fV0Reader->GetV0(indexKF2);
2122 Double_t eta2 = fV0Reader->GetMotherCandidateEta();
2123 Int_t gamma2MotherLabel=-1;
2124 if(fV0Reader->HasSameMCMother() == kTRUE){
2125 TParticle * negativeMC = (TParticle*)fV0Reader->GetNegativeMCParticle();
2126 TParticle * positiveMC = (TParticle*)fV0Reader->GetPositiveMCParticle();
2127 if(TMath::Abs(negativeMC->GetPdgCode())==11 && TMath::Abs(positiveMC->GetPdgCode())==11){
2128 if(negativeMC->GetUniqueID() == 5 && positiveMC->GetUniqueID() ==5){
2129 if(fV0Reader->GetMotherMCParticle()->GetPdgCode() == 22){
2130 gamma2MotherLabel=fV0Reader->GetMotherMCParticle()->GetFirstMother();
2131 }
2132 }
2133 }
2134 }
2135 if(gamma1MotherLabel>=0 && gamma1MotherLabel==gamma2MotherLabel){
2136 if(fV0Reader->CheckIfPi0IsMother(gamma1MotherLabel)){
2137 isRealPi0=kTRUE;
2138 }
10e3319b 2139 if(fV0Reader->CheckIfEtaIsMother(gamma1MotherLabel)){
2140 isRealEta=kTRUE;
2141 }
2142
037dc2db 2143 }
2144
2145 if(TMath::Abs(eta1)>0.9 && TMath::Abs(eta2)>0.9){
ebcfaa7e 2146 // fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
2147 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
10e3319b 2148 if(isRealPi0 || isRealEta){
037dc2db 2149 fHistograms->FillHistogram("ESD_TruePi0_InvMass_1212",massTwoGammaCandidate);
2150 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_1212",openingAngleTwoGammaCandidate);
2151 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 2152 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2153 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2154 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2155 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2156 }
037dc2db 2157 }
bd6d9fa3 2158 if(!isRealPi0 && !isRealEta){
2159 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
2160 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
2161 }else{
2162 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
2163 }
2164 }
2165
037dc2db 2166 }
2167 else if(TMath::Abs(eta1)>0.9 || TMath::Abs(eta2)>0.9){
ebcfaa7e 2168 // fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
2169 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
10e3319b 2170 if(isRealPi0 || isRealEta){
037dc2db 2171 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0912",massTwoGammaCandidate);
2172 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0912",openingAngleTwoGammaCandidate);
2173 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 2174 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2175 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2176 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2177 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2178 }
037dc2db 2179 }
bd6d9fa3 2180 if(!isRealPi0 && !isRealEta){
2181 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
2182 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
2183 }else{
2184 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
2185 }
2186 }
037dc2db 2187 }
2188 else{
ebcfaa7e 2189 // fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
2190 // fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
10e3319b 2191 if(isRealPi0 || isRealEta){
037dc2db 2192 fHistograms->FillHistogram("ESD_TruePi0_InvMass_0909",massTwoGammaCandidate);
2193 fHistograms->FillHistogram("ESD_TruePi0_OpeningAngle_0909",openingAngleTwoGammaCandidate);
2194 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
ebcfaa7e 2195 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2196 fHistograms->FillHistogram("ESD_TruePi0_InvMass",massTwoGammaCandidate);
d707e3cf 2197 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2198 fHistograms->FillHistogram("ESD_TruePi0_InvMass_vs_Pt_alpha",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2199 }
bd6d9fa3 2200 if(gamma1MotherLabel > fV0Reader->GetMCStack()->GetNprimary()){
2201 fHistograms->FillHistogram("ESD_TruePi0Sec_InvMass",massTwoGammaCandidate);
2202 }
2203 }
2204 if(!isRealPi0 && !isRealEta){
2205 if(gamma1MotherLabel>-1 && gamma2MotherLabel>-1){
2206 fHistograms->FillHistogram("ESD_TrueBckGG_InvMass",massTwoGammaCandidate);
2207 }else{
2208 fHistograms->FillHistogram("ESD_TrueBckCont_InvMass",massTwoGammaCandidate);
2209 }
037dc2db 2210 }
2211 }
2212 }
2213 }
2214 }
9640a3d1 2215 if ( TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())<0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())<0.9 ){
2216 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_Fiducial",massTwoGammaCandidate ,momentumVectorTwoGammaCandidate.Pt());
2217 fHistograms->FillHistogram("ESD_Mother_InvMass_Fiducial",massTwoGammaCandidate);
2218 }
ebcfaa7e 2219
2220 if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 && TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
2221 fHistograms->FillHistogram("ESD_Mother_InvMass_1212",massTwoGammaCandidate);
2222 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt1212",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2223 }
2224 else if(TMath::Abs(twoGammaDecayCandidateDaughter0->GetEta())>0.9 || TMath::Abs(twoGammaDecayCandidateDaughter1->GetEta())>0.9){
2225 fHistograms->FillHistogram("ESD_Mother_InvMass_0912",massTwoGammaCandidate);
2226 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0912",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2227 }
2228 else{
2229 fHistograms->FillHistogram("ESD_Mother_InvMass_0909",massTwoGammaCandidate);
2230 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt0909",massTwoGammaCandidate,momentumVectorTwoGammaCandidate.Pt());
2231 }
6272370b 2232 Double_t lowMassPi0=0.1;
2233 Double_t highMassPi0=0.15;
2234 if (massTwoGammaCandidate > lowMassPi0 && massTwoGammaCandidate < highMassPi0 ){
2235 new((*fKFReconstructedPi0sTClone)[fKFReconstructedPi0sTClone->GetEntriesFast()]) AliKFParticle(*twoGammaCandidate);
2236 fGammav1.push_back(firstGammaIndex);
2237 fGammav2.push_back(secondGammaIndex);
c59360eb 2238 AddPionToAOD(twoGammaCandidate, massTwoGammaCandidate, firstGammaIndex, secondGammaIndex);
6272370b 2239 }
2240
2e2da371 2241 }
d707e3cf 2242 //}
2243 delete twoGammaCandidate;
d7d7e825 2244 }
2245 }
2246}
2247
04bf4381 2248void AliAnalysisTaskGammaConversion::AddPionToAOD(AliKFParticle * pionkf, Double_t mass, Int_t daughter1, Int_t daughter2) {
2249 //See header file for documentation
2250 AliGammaConversionAODObject pion;
2251 pion.SetPx(pionkf->Px());
2252 pion.SetPy(pionkf->Py());
2253 pion.SetPz(pionkf->Pz());
2254 pion.SetChi2(pionkf->GetChi2());
2255 pion.SetE(pionkf->E());
2256 pion.SetIMass(mass);
2257 pion.SetLabel1(daughter1);
2258 //dynamic_cast<AliGammaConversionAODObject*>(fAODBranch->At(daughter1))->SetTagged(kTRUE);
2259 pion.SetLabel2(daughter2);
2260 new((*fAODPi0)[fAODPi0->GetEntriesFast()]) AliGammaConversionAODObject(pion);
2261
2262}
2263
2264
2265
6272370b 2266void AliAnalysisTaskGammaConversion::ProcessConvPHOSGammasForNeutralMesonAnalysis(){
48682642 2267/*
6272370b 2268 // see header file for documentation
2269 // Analyse Pi0 with one photon from Phos and 1 photon from conversions
2270
2271
2272
2273 Double_t vtx[3];
2274 vtx[0] = fV0Reader->GetPrimaryVertex()->GetX();
2275 vtx[1] = fV0Reader->GetPrimaryVertex()->GetY();
2276 vtx[2] = fV0Reader->GetPrimaryVertex()->GetZ();
2277
2278
2279 // Loop over all CaloClusters and consider only the PHOS ones:
2280 AliESDCaloCluster *clu;
2281 TLorentzVector pPHOS;
2282 TLorentzVector gammaPHOS;
2283 TLorentzVector gammaGammaConv;
2284 TLorentzVector pi0GammaConvPHOS;
2285 TLorentzVector gammaGammaConvBck;
2286 TLorentzVector pi0GammaConvPHOSBck;
2287
2288
2289 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
2290 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2291 if ( !clu->IsPHOS() || clu->E()<0.1 ) continue;
2292 clu ->GetMomentum(pPHOS ,vtx);
2293 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2294 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2295 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz(),0.);
2296 gammaPHOS.SetXYZM(pPHOS.Px(),pPHOS.Py(),pPHOS.Pz(),0.);
2297 pi0GammaConvPHOS=gammaGammaConv+gammaPHOS;
2298 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS",pi0GammaConvPHOS.M());
2299 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvPHOS",pi0GammaConvPHOS.M(),pi0GammaConvPHOS.Pt());
2300
2301 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),twoGammaDecayCandidateDaughter0->Py(),twoGammaDecayCandidateDaughter0->Pz());
2302 TVector3 v3D1(gammaPHOS.Px(),gammaPHOS.Py(),gammaPHOS.Pz());
2303 Double_t opanConvPHOS= v3D0.Angle(v3D1);
2304 if ( opanConvPHOS < 0.35){
2305 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanLow",pi0GammaConvPHOS.M());
2306 }else{
2307 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvPHOS_OpanHigh",pi0GammaConvPHOS.M());
2308 }
2309
2310 }
2311
2312 // Now the LorentVector pPHOS is obtained and can be paired with the converted proton
2313 }
2314 //==== End of the PHOS cluster selection ============
2315 TLorentzVector pEMCAL;
2316 TLorentzVector gammaEMCAL;
2317 TLorentzVector pi0GammaConvEMCAL;
2318 TLorentzVector pi0GammaConvEMCALBck;
2319
2320 for (Int_t i=0; i<fV0Reader->GetESDEvent()->GetNumberOfCaloClusters(); i++) {
2321 clu = fV0Reader->GetESDEvent()->GetCaloCluster(i);
2322 if ( !clu->IsEMCAL() || clu->E()<0.1 ) continue;
2323 if (clu->GetNCells() <= 1) continue;
2324 if ( clu->GetTOF()*1e9 < 550 || clu->GetTOF()*1e9 > 750) continue;
2325
2326 clu ->GetMomentum(pEMCAL ,vtx);
2327 for(Int_t firstGammaIndex=0;firstGammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();firstGammaIndex++){
2328 AliKFParticle * twoGammaDecayCandidateDaughter0 = (AliKFParticle *)fKFReconstructedGammasTClone->At(firstGammaIndex);
2329 gammaGammaConv.SetXYZM(twoGammaDecayCandidateDaughter0->Px(),
2330 twoGammaDecayCandidateDaughter0->Py(),
2331 twoGammaDecayCandidateDaughter0->Pz(),0.);
2332 gammaEMCAL.SetXYZM(pEMCAL.Px(),pEMCAL.Py(),pEMCAL.Pz(),0.);
2333 pi0GammaConvEMCAL=gammaGammaConv+gammaEMCAL;
2334 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL",pi0GammaConvEMCAL.M());
2335 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL",pi0GammaConvEMCAL.M(),pi0GammaConvEMCAL.Pt());
2336 TVector3 v3D0(twoGammaDecayCandidateDaughter0->Px(),
2337 twoGammaDecayCandidateDaughter0->Py(),
2338 twoGammaDecayCandidateDaughter0->Pz());
2339 TVector3 v3D1(gammaEMCAL.Px(),gammaEMCAL.Py(),gammaEMCAL.Pz());
2340
2341
2342 Double_t opanConvEMCAL= v3D0.Angle(v3D1);
2343 if ( opanConvEMCAL < 0.35){
2344 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanLow",pi0GammaConvEMCAL.M());
2345 }else{
2346 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_OpanHigh",pi0GammaConvEMCAL.M());
2347 }
2348
2349 }
48682642 2350 if(fCalculateBackground){
2351 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2352 AliGammaConversionKFVector * previousEventV0s = fV0Reader->GetBGGoodV0s(nEventsInBG);
2353 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2354 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2355 gammaGammaConvBck.SetXYZM(previousGoodV0.Px(),
2356 previousGoodV0.Py(),
2357 previousGoodV0.Pz(),0.);
2358 pi0GammaConvEMCALBck=gammaGammaConvBck+gammaEMCAL;
2359 fHistograms->FillHistogram("ESD_Mother_InvMass_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M());
2360 fHistograms->FillHistogram("ESD_Mother_InvMass_vs_Pt_GammaConvEMCAL_Bck",pi0GammaConvEMCALBck.M(),
2361 pi0GammaConvEMCALBck.Pt());
2362 }
6272370b 2363 }
48682642 2364
2365 // Now the LorentVector pEMCAL is obtained and can be paired with the converted proton
2366 } // end of checking if background photons are available
2367 }
6272370b 2368 //==== End of the PHOS cluster selection ============
48682642 2369*/
6272370b 2370}
2371
5ce758b0 2372void AliAnalysisTaskGammaConversion::MoveParticleAccordingToVertex(AliKFParticle *particle,AliGammaConversionBGHandler::GammaConversionVertex *vertex){
2373 //see header file for documentation
2374
2375 Double_t dx = vertex->fX - fESDEvent->GetPrimaryVertex()->GetX();
2376 Double_t dy = vertex->fY - fESDEvent->GetPrimaryVertex()->GetY();
2377 Double_t dz = vertex->fZ - fESDEvent->GetPrimaryVertex()->GetZ();
2378
2379 // cout<<"dx, dy, dz: ["<<dx<<","<<dy<<","<<dz<<"]"<<endl;
2380 particle->X() = particle->GetX() - dx;
2381 particle->Y() = particle->GetY() - dy;
2382 particle->Z() = particle->GetZ() - dz;
2383}
2384
111d75df 2385void AliAnalysisTaskGammaConversion::RotateKFParticle(AliKFParticle * kfParticle,Double_t angle){
2386
2387 Double_t c = cos(angle);
2388 Double_t s = sin(angle);
2389
2390 Double_t A[7][ 7];
2391 for( Int_t i=0; i<7; i++ ){
2392 for( Int_t j=0; j<7; j++){
2393 A[i][j] = 0;
2394 }
2395 }
2396 for( int i=0; i<7; i++ ){
2397 A[i][i] = 1;
2398 }
2399 A[0][0] = c; A[0][1] = s;
2400 A[1][0] = -s; A[1][1] = c;
2401 A[3][3] = c; A[3][4] = s;
2402 A[4][3] = -s; A[4][4] = c;
2403
2404 Double_t AC[7][7];
2405 Double_t Ap[7];
2406
2407 for( Int_t i=0; i<7; i++ ){
2408 Ap[i] = 0;
2409 for( Int_t k=0; k<7; k++){
2410 Ap[i]+=A[i][k] * kfParticle->GetParameter(k);
2411 }
2412 }
2413
2414 for( Int_t i=0; i<7; i++){
2415 kfParticle->Parameter(i) = Ap[i];
2416 }
2417
2418 for( Int_t i=0; i<7; i++ ){
2419 for( Int_t j=0; j<7; j++ ){
2420 AC[i][j] = 0;
2421 for( Int_t k=0; k<7; k++ ){
2422 AC[i][j]+= A[i][k] * kfParticle->GetCovariance(k,j);
2423 }
2424 }
2425 }
2426
2427 for( Int_t i=0; i<7; i++ ){
2428 for( Int_t j=0; j<=i; j++ ){
2429 Double_t xx = 0;
2430 for( Int_t k=0; k<7; k++){
2431 xx+= AC[i][k]*A[j][k];
2432 }
2433 kfParticle->Covariance(i,j) = xx;
2434 }
2435 }
2436}
2437
2438
d7d7e825 2439void AliAnalysisTaskGammaConversion::CalculateBackground(){
2440 // see header file for documentation
5e55d806 2441
2442
2443 TClonesArray * currentEventV0s = fV0Reader->GetCurrentEventGoodV0s();
2444
5ce758b0 2445 AliGammaConversionBGHandler * bgHandler = fV0Reader->GetBGHandler();
2446
2447 Int_t zbin= bgHandler->GetZBinIndex(fV0Reader->GetVertexZ());
2448 Int_t mbin = 0;
2449 if(fUseTrackMultiplicityForBG == kTRUE){
2450 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks());
2451 }
2452 else{
2453 mbin = bgHandler->GetMultiplicityBinIndex(fV0Reader->GetNGoodV0s());
2454 }
2455
111d75df 2456 if(fDoRotation == kTRUE){
3c45d101 2457 TRandom3 *random = new TRandom3();
111d75df 2458 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2459 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2460 for(Int_t iCurrent2=iCurrent+1;iCurrent2<currentEventV0s->GetEntriesFast();iCurrent2++){
2461 for(Int_t nRandom=0;nRandom<nRandomEventsForBG;nRandom++){
2462
2463 AliKFParticle currentEventGoodV02 = *(AliKFParticle *)(currentEventV0s->At(iCurrent2));
3bfbe89a 2464
2465 if(fCheckBGProbability == kTRUE){
2466 Double_t massBGprob =0.;
2467 Double_t widthBGprob = 0.;
2468 AliKFParticle *backgroundCandidateProb = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
2469 backgroundCandidateProb->GetMass(massBGprob,widthBGprob);
2470 if(massBGprob>0.1 && massBGprob<0.14){
2471 if(random->Rndm()>bgHandler->GetBGProb(zbin,mbin)){
2472 delete backgroundCandidateProb;
2473 continue;
2474 }
2475 }
2476 delete backgroundCandidateProb;
2477 }
111d75df 2478
2479 Double_t nRadiansPM = nDegreesPMBackground*TMath::Pi()/180;
2480
2481 Double_t rotationValue = random->Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
5ce758b0 2482
111d75df 2483 RotateKFParticle(&currentEventGoodV02,rotationValue);
5ce758b0 2484
111d75df 2485 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,currentEventGoodV02);
5e55d806 2486
5ce758b0 2487 Double_t massBG =0.;
2488 Double_t widthBG = 0.;
2489 Double_t chi2BG =10000.;
2490 backgroundCandidate->GetMass(massBG,widthBG);
111d75df 2491
111d75df 2492 // if(backgroundCandidate->GetNDF()>0){
2493 chi2BG = backgroundCandidate->GetChi2();
2494 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
2495
2496 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2497 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2498
2499 Double_t openingAngleBG = currentEventGoodV0.GetAngle(currentEventGoodV02);
2500
2501 Double_t rapidity;
2502 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
2503 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2504
dc2883e4 2505 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2506 delete backgroundCandidate;
2507 continue; // rapidity cut
2508 }
2509
111d75df 2510
2511 Double_t alfa=0.0;
2512 if( (currentEventGoodV0.GetE()+currentEventGoodV02.GetE()) != 0){
2513 alfa=TMath::Abs((currentEventGoodV0.GetE()-currentEventGoodV02.GetE())
2514 /(currentEventGoodV0.GetE()+currentEventGoodV02.GetE()));
2515 }
2516
2517
2518 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2519 delete backgroundCandidate;
2520 continue; // minimum opening angle to avoid using ghosttracks
2521 }
2522
2523 // original
2524 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2525 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2526 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2527 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2528 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2529 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2530 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2531 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2532 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2533 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2534 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2535 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
2536
2537 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2538 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
2539 }
3c45d101 2540 if(alfa<0.1){
2541 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
2542 }
2543
111d75df 2544 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
2545 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
2546 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
2547 }
2548
2549 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
2550 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
2551 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
2552 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
2553 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
2554 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
2555 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
2556 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2557 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2558 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
2559 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2560 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
2561
2562 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(currentEventGoodV02.GetEta())<0.9 ){
2563 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2564 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
2565 }
2e2da371 2566 }
111d75df 2567 //}
5ce758b0 2568 delete backgroundCandidate;
2569 }
2570 }
2571 }
111d75df 2572 delete random;
5ce758b0 2573 }
111d75df 2574 else{ // means no rotation
2575 AliGammaConversionBGHandler::GammaConversionVertex *bgEventVertex = NULL;
2576
2577 if(fUseTrackMultiplicityForBG){
2578 // cout<<"Using charged track multiplicity for background calculation"<<endl;
2579 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
5ce758b0 2580
111d75df 2581 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);//fV0Reader->GetBGGoodV0s(nEventsInBG);
2582
5ce758b0 2583 if(fMoveParticleAccordingToVertex == kTRUE){
2584 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
2585 }
2586
2587 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2588 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2589 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2590 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
111d75df 2591 AliKFParticle previousGoodV0test = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
5ce758b0 2592
111d75df 2593 //cout<<"Primary Vertex event: ["<<fESDEvent->GetPrimaryVertex()->GetX()<<","<<fESDEvent->GetPrimaryVertex()->GetY()<<","<<fESDEvent->GetPrimaryVertex()->GetZ()<<"]"<<endl;
2594 //cout<<"BG prim Vertex event: ["<<bgEventVertex->fX<<","<<bgEventVertex->fY<<","<<bgEventVertex->fZ<<"]"<<endl;
2595
2596 //cout<<"XYZ of particle before transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
5ce758b0 2597 if(fMoveParticleAccordingToVertex == kTRUE){
2598 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
2599 }
111d75df 2600 //cout<<"XYZ of particle after transport: ["<<previousGoodV0.X()<<","<<previousGoodV0.Y()<<","<<previousGoodV0.Z()<<"]"<<endl;
5ce758b0 2601
2602 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
111d75df 2603
5ce758b0 2604 Double_t massBG =0.;
2605 Double_t widthBG = 0.;
2606 Double_t chi2BG =10000.;
2607 backgroundCandidate->GetMass(massBG,widthBG);
111d75df 2608 // if(backgroundCandidate->GetNDF()>0){
2609 // chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
2610 chi2BG = backgroundCandidate->GetChi2();
2611 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
2612
2613 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2614 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2615
2616 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
2617
2618 Double_t rapidity;
2619
2620 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() <= 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() <= 0){
2621 cout << "Error: |Pz| > E !!!! " << endl;
2622 rapidity=0;
2623 } else {
2624 rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2625 }
dc2883e4 2626 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2627 delete backgroundCandidate;
2628 continue; // rapidity cut
2629 }
2630
2631
111d75df 2632 Double_t alfa=0.0;
2633 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
2634 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
2635 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
2636 }
2637
2638
2639 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2640 delete backgroundCandidate;
2641 continue; // minimum opening angle to avoid using ghosttracks
2642 }
2643
2644 // original
2645 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2646 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2647 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2648 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2649 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2650 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2651 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2652 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2653 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2654 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2655 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2656 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
2657
2658 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2659 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
2660 }
3c45d101 2661 if(alfa<0.1){
2662 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
2663 }
2664
111d75df 2665 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2666 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
2667 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
2668 }
2669
2670
2671 // test
2672 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
2673 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
2674 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
2675 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
2676 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
2677 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
2678 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
2679 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2680 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2681 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
2682 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2683 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
2684
2685 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2686 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2687 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
2688 }
2689 // }
2690 }
2691 delete backgroundCandidate;
2692 }
2693 }
2694 }
2695 }
2696 else{ // means using #V0s for multiplicity
2697
2698 // cout<<"Using the v0 multiplicity to calculate background"<<endl;
2699
2700 fHistograms->FillHistogram("ESD_Background_z_m",zbin,mbin);
2701 fHistograms->FillHistogram("ESD_Mother_multpilicityVSv0s",fV0Reader->CountESDTracks(),fV0Reader->GetNumberOfV0s());
2702
2703 for(Int_t nEventsInBG=0;nEventsInBG <fV0Reader->GetNBGEvents();nEventsInBG++){
2704 AliGammaConversionKFVector * previousEventV0s = bgHandler->GetBGGoodV0s(zbin,mbin,nEventsInBG);// fV0Reader->GetBGGoodV0s(nEventsInBG);
2705 if(previousEventV0s){
2706
2707 if(fMoveParticleAccordingToVertex == kTRUE){
2708 bgEventVertex = bgHandler->GetBGEventVertex(zbin,mbin,nEventsInBG);
2709 }
2710
2711 for(Int_t iCurrent=0;iCurrent<currentEventV0s->GetEntriesFast();iCurrent++){
2712 AliKFParticle currentEventGoodV0 = *(AliKFParticle *)(currentEventV0s->At(iCurrent));
2713 for(UInt_t iPrevious=0;iPrevious<previousEventV0s->size();iPrevious++){
2714 AliKFParticle previousGoodV0 = (AliKFParticle)(*(previousEventV0s->at(iPrevious)));
2715
2716 if(fMoveParticleAccordingToVertex == kTRUE){
2717 MoveParticleAccordingToVertex(&previousGoodV0,bgEventVertex);
5ce758b0 2718 }
111d75df 2719
2720 AliKFParticle *backgroundCandidate = new AliKFParticle(currentEventGoodV0,previousGoodV0);
2721 Double_t massBG =0.;
2722 Double_t widthBG = 0.;
2723 Double_t chi2BG =10000.;
2724 backgroundCandidate->GetMass(massBG,widthBG);
2725 /* if(backgroundCandidate->GetNDF()>0){
2726 chi2BG = backgroundCandidate->GetChi2()/backgroundCandidate->GetNDF();
2727 {//remember to remove
2728 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2729 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2730
2731 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
2732 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle_nochi2", openingAngleBG);
2733 }
2734 */
2e2da371 2735 chi2BG = backgroundCandidate->GetChi2();
2736 if((chi2BG>0 && chi2BG<fV0Reader->GetChi2CutMeson()) || fApplyChi2Cut == kFALSE){
5ce758b0 2737 TVector3 momentumVectorbackgroundCandidate(backgroundCandidate->GetPx(),backgroundCandidate->GetPy(),backgroundCandidate->GetPz());
2738 TVector3 spaceVectorbackgroundCandidate(backgroundCandidate->GetX(),backgroundCandidate->GetY(),backgroundCandidate->GetZ());
2739
2740 Double_t openingAngleBG = currentEventGoodV0.GetAngle(previousGoodV0);
2741
2742 Double_t rapidity;
2743 if(backgroundCandidate->GetE() - backgroundCandidate->GetPz() == 0 || backgroundCandidate->GetE() + backgroundCandidate->GetPz() == 0) rapidity=0;
2744 else rapidity = 0.5*(TMath::Log((backgroundCandidate->GetE() +backgroundCandidate->GetPz()) / (backgroundCandidate->GetE()-backgroundCandidate->GetPz())));
2745
dc2883e4 2746 if(TMath::Abs(rapidity) > fV0Reader->GetRapidityMesonCut() ){
2747 delete backgroundCandidate;
2748 continue; // rapidity cut
2749 }
2750
2751
5ce758b0 2752 Double_t alfa=0.0;
2753 if( (currentEventGoodV0.GetE()+previousGoodV0.GetE()) != 0){
2754 alfa=TMath::Abs((currentEventGoodV0.GetE()-previousGoodV0.GetE())
2755 /(currentEventGoodV0.GetE()+previousGoodV0.GetE()));
2756 }
2757
2758
2759 if(openingAngleBG < fMinOpeningAngleGhostCut ){
2760 delete backgroundCandidate;
2761 continue; // minimum opening angle to avoid using ghosttracks
2762 }
2763
2764 fHistograms->FillHistogram("ESD_Background_GammaDaughter_OpeningAngle", openingAngleBG);
2765 fHistograms->FillHistogram("ESD_Background_Energy", backgroundCandidate->GetE());
2766 fHistograms->FillHistogram("ESD_Background_Pt", momentumVectorbackgroundCandidate.Pt());
2767 fHistograms->FillHistogram("ESD_Background_Eta", momentumVectorbackgroundCandidate.Eta());
2768 fHistograms->FillHistogram("ESD_Background_Rapidity", rapidity);
2769 fHistograms->FillHistogram("ESD_Background_Phi", spaceVectorbackgroundCandidate.Phi());
2770 fHistograms->FillHistogram("ESD_Background_Mass", massBG);
2771 fHistograms->FillHistogram("ESD_Background_R", spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2772 fHistograms->FillHistogram("ESD_Background_ZR", backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2773 fHistograms->FillHistogram("ESD_Background_XY", backgroundCandidate->GetX(), backgroundCandidate->GetY());
2774 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt",massBG,momentumVectorbackgroundCandidate.Pt());
2775 fHistograms->FillHistogram("ESD_Background_InvMass",massBG);
2776
2777 if(alfa>fV0Reader->GetAlphaMinCutMeson() && alfa<fV0Reader->GetAlphaCutMeson()){
2778 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_alpha",massBG,momentumVectorbackgroundCandidate.Pt());
2779 }
3c45d101 2780 if(alfa<0.1){
2781 fHistograms->FillHistogram("ESD_Background_InvMass_vs_E_alpha",massBG ,backgroundCandidate->GetE());
2782 }
2783
5ce758b0 2784 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2785 fHistograms->FillHistogram("ESD_Background_InvMass_vs_Pt_Fiducial",massBG,momentumVectorbackgroundCandidate.Pt());
2786 fHistograms->FillHistogram("ESD_Background_InvMass_Fiducial",massBG);
2787 }
2788
2789 if(massBG>0.5 && massBG<0.6){
2790 fHistograms->FillHistogram("ESD_Background_alfa_pt0506",momentumVectorbackgroundCandidate.Pt(),alfa);
2791 }
2792 if(massBG>0.3 && massBG<0.4){
2793 fHistograms->FillHistogram("ESD_Background_alfa_pt0304",momentumVectorbackgroundCandidate.Pt(),alfa);
2794 }
037dc2db 2795
5ce758b0 2796 // test
2797 fHistograms->FillHistogram(Form("%d%dESD_Background_GammaDaughter_OpeningAngle",zbin,mbin), openingAngleBG);
2798 fHistograms->FillHistogram(Form("%d%dESD_Background_Energy",zbin,mbin), backgroundCandidate->GetE());
2799 fHistograms->FillHistogram(Form("%d%dESD_Background_Pt",zbin,mbin), momentumVectorbackgroundCandidate.Pt());
2800 fHistograms->FillHistogram(Form("%d%dESD_Background_Eta",zbin,mbin), momentumVectorbackgroundCandidate.Eta());
2801 fHistograms->FillHistogram(Form("%d%dESD_Background_Rapidity",zbin,mbin), rapidity);
2802 fHistograms->FillHistogram(Form("%d%dESD_Background_Phi",zbin,mbin), spaceVectorbackgroundCandidate.Phi());
2803 fHistograms->FillHistogram(Form("%d%dESD_Background_Mass",zbin,mbin), massBG);
2804 fHistograms->FillHistogram(Form("%d%dESD_Background_R",zbin,mbin), spaceVectorbackgroundCandidate.Pt()); // Pt in Space == R!!!!
2805 fHistograms->FillHistogram(Form("%d%dESD_Background_ZR",zbin,mbin), backgroundCandidate->GetZ(), spaceVectorbackgroundCandidate.Pt());
2806 fHistograms->FillHistogram(Form("%d%dESD_Background_XY",zbin,mbin), backgroundCandidate->GetX(), backgroundCandidate->GetY());
2807 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2808 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass",zbin,mbin),massBG);
2809
2810 if ( TMath::Abs(currentEventGoodV0.GetEta())<0.9 && TMath::Abs(previousGoodV0.GetEta())<0.9 ){
2811 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_vs_Pt_Fiducial",zbin,mbin),massBG,momentumVectorbackgroundCandidate.Pt());
2812 fHistograms->FillHistogram(Form("%d%dESD_Background_InvMass_Fiducial",zbin,mbin),massBG);
2813 }
2e2da371 2814 // }
5ce758b0 2815 }
111d75df 2816 delete backgroundCandidate;
2817 }
5ce758b0 2818 }
2819 }
d7d7e825 2820 }
111d75df 2821 } // end else (means use #v0s as multiplicity)
2822 } // end no rotation
d7d7e825 2823}
2824
2825
d7d7e825 2826void AliAnalysisTaskGammaConversion::ProcessGammasForGammaJetAnalysis(){
2827 //ProcessGammasForGammaJetAnalysis
a0b94e5c 2828
d7d7e825 2829 Double_t distIsoMin;
a0b94e5c 2830
d7d7e825 2831 CreateListOfChargedParticles();
a0b94e5c 2832
2833
6c84d371 2834 // for(UInt_t gammaIndex=0;gammaIndex<fKFReconstructedGammas.size();gammaIndex++){
2835 for(Int_t gammaIndex=0;gammaIndex<fKFReconstructedGammasTClone->GetEntriesFast();gammaIndex++){
2836 AliKFParticle * currentGamma = (AliKFParticle*)fKFReconstructedGammasTClone->At(gammaIndex);
01b7fdcc 2837 TVector3 momentumVectorCurrentGamma(currentGamma->GetPx(),currentGamma->GetPy(),currentGamma->GetPz());
a0b94e5c 2838
01b7fdcc 2839 if( momentumVectorCurrentGamma.Pt()> fMinPtForGammaJet){
d7d7e825 2840 distIsoMin=GetMinimumDistanceToCharge(gammaIndex);
a0b94e5c 2841
d7d7e825 2842 if (distIsoMin > fMinIsoConeSize && fLeadingChargedIndex>=0){
01b7fdcc 2843 CalculateJetCone(gammaIndex);
d7d7e825 2844 }
2845 }
2846 }
2847}
2848
6272370b 2849//____________________________________________________________________
2850Bool_t AliAnalysisTaskGammaConversion::IsGoodImpPar(AliESDtrack *const track)
2851{
2852//
2853// check whether particle has good DCAr(Pt) impact
2854// parameter. Only for TPC+ITS tracks (7*sigma cut)
2855// Origin: Andrea Dainese
2856//
2857
2858Float_t d0z0[2],covd0z0[3];
2859track->GetImpactParameters(d0z0,covd0z0);
2860Float_t sigma= 0.0050+0.0060/TMath::Power(track->Pt(),0.9);
2861Float_t d0max = 7.*sigma;
2862if(TMath::Abs(d0z0[0]) < d0max) return kTRUE;
2863
2864return kFALSE;
2865}
2866
2867
d7d7e825 2868void AliAnalysisTaskGammaConversion::CreateListOfChargedParticles(){
2869 // CreateListOfChargedParticles
a0b94e5c 2870
d7d7e825 2871 fESDEvent = fV0Reader->GetESDEvent();
037dc2db 2872 Int_t numberOfESDTracks=0;
d7d7e825 2873 for(Int_t iTracks = 0; iTracks < fESDEvent->GetNumberOfTracks(); iTracks++){
2874 AliESDtrack* curTrack = fESDEvent->GetTrack(iTracks);
a0b94e5c 2875
d7d7e825 2876 if(!curTrack){
2877 continue;
2878 }
d707e3cf 2879 // Not needed if Standard function used.
2880// if(!IsGoodImpPar(curTrack)){
2881// continue;
2882// }
a0b94e5c 2883
d7d7e825 2884 if(fEsdTrackCuts->AcceptTrack(curTrack) ){
6c84d371 2885 new((*fChargedParticles)[fChargedParticles->GetEntriesFast()]) AliESDtrack(*curTrack);
2886 // fChargedParticles.push_back(curTrack);
d7d7e825 2887 fChargedParticlesId.push_back(iTracks);
037dc2db 2888 numberOfESDTracks++;
d7d7e825 2889 }
2890 }
037dc2db 2891 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracks",numberOfESDTracks);
e40fd7e2 2892
2893 if (fV0Reader->GetNumberOfContributorsVtx()>=1){
2894 fHistograms->FillHistogram("ESD_NumberOfGoodESDTracksVtx",numberOfESDTracks);
2895 }
d7d7e825 2896}
9c1cb6f7 2897void AliAnalysisTaskGammaConversion::RecalculateV0ForGamma(){
2898
2899 Double_t massE=0.00051099892;
2900 TLorentzVector curElecPos;
2901 TLorentzVector curElecNeg;
2902 TLorentzVector curGamma;
2903
2904 TLorentzVector curGammaAt;
2905 TLorentzVector curElecPosAt;
2906 TLorentzVector curElecNegAt;
2907 AliKFVertex primVtxGamma(*(fESDEvent->GetPrimaryVertex()));
2908 AliKFVertex primVtxImprovedGamma = primVtxGamma;
2909
2910 const AliESDVertex *vtxT3D=fESDEvent->GetPrimaryVertex();
2911
2912 Double_t xPrimaryVertex=vtxT3D->GetXv();
2913 Double_t yPrimaryVertex=vtxT3D->GetYv();
2914 Double_t zPrimaryVertex=vtxT3D->GetZv();
70ef88b5 2915 // Float_t primvertex[3]={xPrimaryVertex,yPrimaryVertex,zPrimaryVertex};
9c1cb6f7