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