]>
Commit | Line | Data |
---|---|---|
7f3c7cc6 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: Pedro González, Pedro Ladrón de Guevara, Ernesto López Torres, * | |
5 | * Eulogio Serradilla * | |
6 | * Version 2 * | |
7 | * * | |
8 | * Permission to use, copy, modify and distribute this software and its * | |
9 | * documentation strictly for non-commercial purposes is hereby granted * | |
10 | * without fee, provided that the above copyright notice appears in all * | |
11 | * copies and that both the copyright notice and this permission notice * | |
12 | * appear in the supporting documentation. The authors make no claims * | |
13 | * about the suitability of this software for any purpose. It is * | |
14 | * provided "as is" without express or implied warranty. * | |
15 | **************************************************************************/ | |
16 | ||
17 | // Analysis task for pi0->e+e-gamma (Dalitz decay) | |
18 | ||
19 | #include <vector> | |
20 | ||
21 | #include "TParticle.h" | |
22 | #include "TPDGCode.h" | |
23 | #include "TMCProcess.h" | |
24 | #include "TDatabasePDG.h" | |
25 | #include "TList.h" | |
26 | #include "TChain.h" | |
27 | #include "TDirectory.h" | |
28 | ||
29 | #include "AliStack.h" | |
30 | #include "AliAnalysisManager.h" | |
31 | #include "AliESDInputHandler.h" | |
32 | #include "AliESDtrack.h" | |
33 | #include "AliMCEvent.h" | |
34 | #include "AliStack.h" | |
35 | #include "AliMCEventHandler.h" | |
36 | #include "AliPID.h" | |
37 | #include "AliLog.h" | |
38 | #include "AliESDtrackCuts.h" | |
39 | #include "AliESDpidCuts.h" | |
40 | #include "AliMCEvent.h" | |
41 | #include "AliESDv0.h" | |
42 | #include "AliESDEvent.h" | |
43 | #include "AliESDpid.h" | |
44 | #include "AliKFParticle.h" | |
45 | #include "AliMCEventHandler.h" | |
46 | #include "AliGammaConversionHistograms.h" | |
47 | #include "AliV0Reader.h" | |
48 | #include "AliKFVertex.h" | |
49 | ||
50 | #include "AliAnalysisTaskGammaConvDalitz.h" | |
51 | #include "TH1.h" | |
52 | ||
53 | ClassImp( AliAnalysisTaskGammaConvDalitz ) | |
54 | ||
55 | //----------------------------------------------------------------------------------------------- | |
56 | AliAnalysisTaskGammaConvDalitz::AliAnalysisTaskGammaConvDalitz(): | |
57 | AliAnalysisTaskSE(), | |
58 | fStack(0), | |
7811d09e | 59 | fGCMCEvent(0), |
7f3c7cc6 | 60 | fESDEvent(0), |
61 | fEposCandidateIndex(), | |
62 | fEnegCandidateIndex(), | |
63 | fGammaCandidatePosIndex(), | |
64 | fGammaCandidateNegIndex(), | |
65 | fGammaCandidates(0), | |
66 | fGammaPool(0), | |
67 | fPoolMaxSize(10), | |
68 | fGamPoolPos(0), | |
69 | fBGEventHandler(0), | |
70 | fOutputContainer(0), | |
71 | fMCTruth(0), | |
72 | fV0Reader(0), | |
73 | fESDpid(0), | |
74 | fESDtrackCuts(0), | |
75 | fITSsaTrackCuts(0), | |
76 | fESDpidCuts(0), | |
77 | fHistograms(0), | |
78 | fStandalone(kFALSE), | |
79 | fDoMC(kFALSE), | |
80 | fComputeBkg(kTRUE), | |
81 | fUseBayesPID(kFALSE), | |
82 | fUseTrackIndexCut(kTRUE), | |
83 | fUsePsiPairCut(kTRUE), | |
84 | fUseMassCut(kFALSE), | |
85 | fUseGammaCut(kFALSE), | |
86 | fReadMagFieldSign(kTRUE), | |
87 | fUseAliKF(kFALSE), | |
88 | fMagFieldSign(1), | |
b9295e41 | 89 | fkElectronMass(0.00051099891), |
7f3c7cc6 | 90 | fPsiPairCut(0.45), |
91 | fDeltaPhiCutMin(0.), | |
92 | fDeltaPhiCutMax(0.12), | |
93 | fMassCutMin(0.), | |
94 | fMassCutMax(0.1), | |
95 | fNSigmaBelowElecTPCbethe(-2.), | |
96 | fNSigmaAboveElecTPCbethe(3.), | |
97 | fNSigmaAbovePionTPCbethe(3.), | |
98 | fNSigmaAboveKaonTPCbethe(3.), | |
99 | fNSigmaAboveProtonTPCbethe(3.), | |
100 | fTrkSelectionCriteria(kGlobalTrack) | |
101 | { | |
7811d09e | 102 | // |
103 | // Default constructor | |
104 | // | |
7f3c7cc6 | 105 | AdoptITSsaTrackCuts(); |
106 | AdoptESDtrackCuts(); | |
107 | AdoptESDpidCuts(); | |
108 | ||
109 | fGammaPool = new TClonesArray("AliKFParticle", fPoolMaxSize); | |
110 | fGammaPool->SetOwner(kTRUE); | |
111 | } | |
112 | ||
113 | //----------------------------------------------------------------------------------------------- | |
114 | AliAnalysisTaskGammaConvDalitz::AliAnalysisTaskGammaConvDalitz( const char* name ): | |
115 | AliAnalysisTaskSE( name ), | |
116 | fStack(0), | |
7811d09e | 117 | fGCMCEvent(0), |
7f3c7cc6 | 118 | fESDEvent(0), |
119 | fEposCandidateIndex(), | |
120 | fEnegCandidateIndex(), | |
121 | fGammaCandidatePosIndex(), | |
122 | fGammaCandidateNegIndex(), | |
123 | fGammaCandidates(0), | |
124 | fGammaPool(0), | |
125 | fPoolMaxSize(10), | |
126 | fGamPoolPos(0), | |
127 | fBGEventHandler(0), | |
128 | fOutputContainer(0), | |
129 | fMCTruth(0), | |
130 | fV0Reader(0), | |
131 | fESDpid(0), | |
132 | fESDtrackCuts(0), | |
133 | fITSsaTrackCuts(0), | |
134 | fESDpidCuts(0), | |
135 | fHistograms(0), | |
136 | fStandalone(kFALSE), | |
137 | fDoMC(kFALSE), | |
138 | fComputeBkg(kTRUE), | |
139 | fUseBayesPID(kFALSE), | |
140 | fUseTrackIndexCut(kTRUE), | |
141 | fUsePsiPairCut(kTRUE), | |
142 | fUseMassCut(kFALSE), | |
143 | fUseGammaCut(kFALSE), | |
144 | fReadMagFieldSign(kTRUE), | |
145 | fUseAliKF(kFALSE), | |
146 | fMagFieldSign(1), | |
b9295e41 | 147 | fkElectronMass(0.00051099891), |
7f3c7cc6 | 148 | fPsiPairCut(0.45), |
149 | fDeltaPhiCutMin(0.), | |
150 | fDeltaPhiCutMax(0.12), | |
151 | fMassCutMin(0.), | |
152 | fMassCutMax(0.1), | |
153 | fNSigmaBelowElecTPCbethe(-2.), | |
154 | fNSigmaAboveElecTPCbethe(3.), | |
155 | fNSigmaAbovePionTPCbethe(3.), | |
156 | fNSigmaAboveKaonTPCbethe(3.), | |
157 | fNSigmaAboveProtonTPCbethe(3.), | |
158 | fTrkSelectionCriteria(kGlobalTrack) | |
159 | { | |
160 | // Common I/O in slot 0 | |
161 | DefineInput (0, TChain::Class()); | |
162 | ||
163 | // Your private output | |
164 | DefineOutput(1, TList::Class()); | |
165 | // DefineOutput(2, AliCFContainer::Class()); // for CF | |
166 | ||
167 | AdoptITSsaTrackCuts(); | |
168 | AdoptESDtrackCuts(); | |
169 | AdoptESDpidCuts(); | |
b9295e41 | 170 | // fkElectronMass = TDatabasePDG::Instance()->GetParticle( ::kElectron )->Mass(); // |
7f3c7cc6 | 171 | |
172 | fGammaPool = new TClonesArray("AliKFParticle", fPoolMaxSize); | |
173 | fGammaPool->SetOwner(kTRUE); | |
174 | } | |
175 | ||
176 | //----------------------------------------------------------------------------------------------- | |
177 | AliAnalysisTaskGammaConvDalitz::~AliAnalysisTaskGammaConvDalitz() | |
178 | { | |
7811d09e | 179 | // |
180 | // virtual destructor | |
181 | // | |
182 | ||
7f3c7cc6 | 183 | if( fOutputContainer ) delete fOutputContainer; |
184 | if( fHistograms ) delete fHistograms; | |
185 | if( fStandalone && fV0Reader ) delete fV0Reader; | |
186 | if( fITSsaTrackCuts ) delete fITSsaTrackCuts; | |
187 | if( fESDtrackCuts ) delete fESDtrackCuts; | |
188 | if( fESDpidCuts ) delete fESDpidCuts; | |
189 | if( fGammaCandidates) delete fGammaCandidates; | |
190 | if( fGammaPool ) delete fGammaPool; | |
191 | } | |
192 | ||
193 | //----------------------------------------------------------------------------------------------- | |
194 | void AliAnalysisTaskGammaConvDalitz::ConnectInputData(Option_t *option) | |
195 | { | |
7811d09e | 196 | // |
197 | // Connect Input Data | |
198 | // | |
7f3c7cc6 | 199 | if( fDebug ) AliInfo("=> ConnectInputData"); |
200 | ||
201 | AliAnalysisTaskSE::ConnectInputData(option); | |
202 | ||
203 | if( fV0Reader == 0 ) | |
204 | { | |
205 | AliFatal("There is not pointer to AliV0Reader object!!!"); | |
206 | } | |
207 | } | |
208 | ||
209 | //----------------------------------------------------------------------------------------------- | |
210 | void AliAnalysisTaskGammaConvDalitz::UserCreateOutputObjects() | |
211 | { | |
7811d09e | 212 | // |
213 | // Create ouput objects | |
214 | // | |
7f3c7cc6 | 215 | if( fDebug ) AliInfo("=> UserCreateOutputObjects"); |
216 | ||
217 | // Create the output container | |
218 | if( fOutputContainer != 0 ) | |
219 | { | |
220 | delete fOutputContainer; | |
221 | } | |
222 | ||
223 | fOutputContainer = new TList(); | |
224 | ||
225 | // Add the histograms to the output container | |
226 | fHistograms->GetOutputContainer( fOutputContainer ); | |
227 | fOutputContainer->SetOwner(kTRUE); | |
228 | ||
229 | PostData( 1, fOutputContainer ); | |
230 | } | |
231 | ||
232 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 233 | void AliAnalysisTaskGammaConvDalitz::UserExec(Option_t */*option*/) |
234 | { | |
7811d09e | 235 | // |
236 | // Execute analysis for current event | |
237 | // | |
7f3c7cc6 | 238 | if( fDebug ) AliInfo("=> UserExec"); |
239 | ||
240 | if( fV0Reader == 0 ) | |
241 | { | |
242 | AliFatal("no pointer to AliV0Reader"); | |
243 | return; | |
244 | } | |
245 | ||
7f3c7cc6 | 246 | // Create list of gamma candidates in standalone mode |
247 | // otherwise use the created ones by AliAnalysisTaskGammaConversion | |
248 | if( fStandalone ) | |
249 | { | |
7811d09e | 250 | |
251 | AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager(); | |
252 | AliESDInputHandler *esdHandler=0; | |
253 | if ( (esdHandler=dynamic_cast<AliESDInputHandler*>(man->GetInputEventHandler())) && esdHandler->GetESDpid() ){ | |
254 | AliV0Reader::SetESDpid(esdHandler->GetESDpid()); | |
255 | } else { | |
256 | //load esd pid bethe bloch parameters depending on the existance of the MC handler | |
257 | // yes: MC parameters | |
258 | // no: data parameters | |
259 | if (!AliV0Reader::GetESDpid()){ | |
9cb96e64 | 260 | if (MCEvent() ) { |
7811d09e | 261 | AliV0Reader::InitESDpid(); |
262 | } else { | |
263 | AliV0Reader::InitESDpid(1); | |
264 | } | |
265 | } | |
266 | } | |
267 | ||
268 | ||
d1758e4f | 269 | if (MCEvent() ) { |
7811d09e | 270 | |
271 | // To avoid crashes due to unzip errors. Sometimes the trees are not there. | |
272 | AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()); | |
273 | ||
274 | if (!mcHandler){ | |
275 | AliError("Could not retrive MC event handler!"); | |
276 | return; | |
277 | } | |
278 | ||
279 | if (!mcHandler->InitOk() ){ | |
280 | ||
281 | return; | |
282 | } | |
283 | if (!mcHandler->TreeK() ){ | |
284 | return; | |
285 | } | |
286 | if (!mcHandler->TreeTR() ) { | |
287 | return; | |
288 | } | |
289 | } | |
290 | ||
291 | ||
292 | ||
293 | ||
294 | fV0Reader->SetInputAndMCEvent( InputEvent(), MCEvent() ); | |
7f3c7cc6 | 295 | fV0Reader->Initialize(); |
296 | } | |
297 | ||
298 | if( fV0Reader->CheckForPrimaryVertex() == kFALSE ) | |
299 | { | |
300 | if( fDebug ) AliInfo("no contributors to primary vertex"); | |
301 | return; | |
302 | } | |
303 | ||
7811d09e | 304 | if( fV0Reader->CheckForPrimaryVertexZ() == kFALSE ) |
305 | { | |
306 | ||
307 | if( fDebug ) AliInfo("z vertex out of range"); | |
308 | return; | |
309 | } | |
310 | ||
7f3c7cc6 | 311 | // Get Pointers |
312 | fBGEventHandler = fV0Reader->GetBGHandler(); | |
313 | fESDpid = fV0Reader->GetESDpid(); | |
314 | fESDEvent = fV0Reader->GetESDEvent(); | |
18922343 | 315 | if(fDoMC && MCEvent()) |
7f3c7cc6 | 316 | { |
7811d09e | 317 | fStack= MCEvent()->Stack(); |
318 | fGCMCEvent=MCEvent(); | |
7f3c7cc6 | 319 | } |
320 | ||
321 | // Read the magnetic field sign from ESD | |
322 | if ( fReadMagFieldSign == kTRUE ) | |
323 | { | |
324 | fMagFieldSign = (fESDEvent->GetMagneticField() < 0) ? 1 : -1; | |
325 | } | |
326 | ||
327 | // Process MC information | |
328 | if(fDoMC) | |
329 | { | |
330 | ProcessMCData(); | |
331 | } | |
332 | ||
333 | if( fStandalone ){ | |
334 | while(fV0Reader->NextV0()){}; //SelectGammas | |
335 | fV0Reader->ResetV0IndexNumber(); | |
336 | } | |
7811d09e | 337 | |
7f3c7cc6 | 338 | CreateListOfDalitzPairCandidates(); |
339 | ProcessGammaElectronsForDalitzAnalysis(); | |
7811d09e | 340 | |
341 | if ( fStandalone ){ | |
342 | ||
343 | fV0Reader->UpdateEventByEventData(); | |
344 | ||
345 | } | |
7f3c7cc6 | 346 | |
347 | PostData( 1, fOutputContainer ); | |
348 | } | |
349 | ||
350 | ||
351 | void AliAnalysisTaskGammaConvDalitz::Terminate(Option_t */*option*/) | |
352 | { | |
7811d09e | 353 | // |
7f3c7cc6 | 354 | if( fDebug ) AliInfo("Not to do anything in Terminate"); |
355 | } | |
356 | ||
357 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 358 | void AliAnalysisTaskGammaConvDalitz::AdoptITSsaTrackCuts( AliESDtrackCuts* esdCuts ) |
359 | { | |
7811d09e | 360 | // |
361 | // set user ITSsa track cuts | |
362 | // | |
7f3c7cc6 | 363 | if( fITSsaTrackCuts ) delete fITSsaTrackCuts; |
364 | ||
365 | if( esdCuts ) | |
366 | { | |
367 | fITSsaTrackCuts = esdCuts; | |
368 | } | |
369 | else | |
370 | { | |
7811d09e | 371 | // default cuts |
7f3c7cc6 | 372 | fITSsaTrackCuts = new AliESDtrackCuts("Default ITSsa track cuts for Pi0 Dalitz decay"); |
373 | ||
374 | fITSsaTrackCuts->SetEtaRange( -0.9, 0.9 ); | |
375 | fITSsaTrackCuts->SetAcceptKinkDaughters(kFALSE); | |
376 | ||
377 | fITSsaTrackCuts->SetMinNClustersITS(2); | |
378 | fITSsaTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kFirst); | |
379 | fITSsaTrackCuts->SetRequireITSRefit(kTRUE); | |
380 | ||
381 | fITSsaTrackCuts->SetRequireSigmaToVertex(kTRUE); | |
382 | fITSsaTrackCuts->SetMaxNsigmaToVertex(3); | |
383 | ||
384 | fITSsaTrackCuts->SetRequireITSStandAlone(kTRUE); | |
385 | } | |
386 | } | |
387 | ||
388 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 389 | void AliAnalysisTaskGammaConvDalitz::AdoptESDtrackCuts( AliESDtrackCuts* esdCuts ) |
390 | { | |
7811d09e | 391 | // |
392 | // set user global track cuts | |
393 | // | |
7f3c7cc6 | 394 | if( fESDtrackCuts ) delete fESDtrackCuts; |
395 | ||
396 | if( esdCuts ) | |
397 | { | |
398 | fESDtrackCuts = esdCuts; | |
399 | } | |
400 | else | |
401 | { | |
7811d09e | 402 | //default cuts |
7f3c7cc6 | 403 | fESDtrackCuts = new AliESDtrackCuts("Default global track cuts for Pi0 Dalitz decay"); |
404 | ||
405 | fESDtrackCuts->SetEtaRange( -0.9, 0.9 ); | |
406 | fESDtrackCuts->SetAcceptKinkDaughters(kFALSE); | |
407 | ||
408 | fESDtrackCuts->SetMinNClustersITS(2); | |
409 | fESDtrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kFirst); | |
410 | fESDtrackCuts->SetRequireITSRefit(kTRUE); | |
411 | ||
412 | fESDtrackCuts->SetRequireSigmaToVertex(kTRUE); | |
413 | fESDtrackCuts->SetMaxNsigmaToVertex(3); | |
414 | ||
415 | fESDtrackCuts->SetMinNClustersTPC(80); | |
416 | fESDtrackCuts->SetMaxChi2PerClusterTPC(4.); | |
417 | fESDtrackCuts->SetMaxCovDiagonalElements(2,2,0.5,0.5,2); | |
418 | fESDtrackCuts->SetRequireTPCRefit(kTRUE); | |
419 | } | |
420 | } | |
421 | ||
422 | //----------------------------------------------------------------------------------------------- | |
423 | void AliAnalysisTaskGammaConvDalitz::AdoptESDpidCuts( AliESDpidCuts* esdPIDCuts ) | |
424 | { | |
7811d09e | 425 | // |
426 | // set user pid cuts | |
427 | // | |
7f3c7cc6 | 428 | if( fESDpidCuts ) delete fESDpidCuts; |
429 | if( esdPIDCuts ) | |
430 | { | |
431 | fESDpidCuts = esdPIDCuts; | |
432 | } | |
7811d09e | 433 | else // default cuts |
7f3c7cc6 | 434 | { |
435 | fESDpidCuts = new AliESDpidCuts("Electrons", "Electron PID cuts"); | |
436 | // fESDpidCuts->SetTPCnSigmaCut(AliPID::kElectron, 3.); | |
437 | fESDpidCuts->SetTPCnSigmaCut(AliPID::kElectron, -4., 6.); | |
438 | } | |
439 | } | |
440 | ||
441 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 442 | void AliAnalysisTaskGammaConvDalitz::ProcessMCData() |
443 | { | |
7811d09e | 444 | // |
445 | // Process generation | |
446 | // | |
7f3c7cc6 | 447 | if( fDebug ) AliInfo("=> ProcessMCData"); |
448 | ||
449 | fHistograms->FillTable("Table_Generation", 0); //number of events | |
450 | ||
451 | for ( Int_t i = 0; i < fStack->GetNtrack(); i++ ) | |
452 | { | |
453 | TParticle* iParticle = fStack->Particle( i ); | |
454 | if( !iParticle ) continue; | |
455 | ||
456 | if ( i >= fStack->GetNprimary() ) continue; // is primary? | |
457 | if ( iParticle->GetPdgCode() != ::kPi0 ) continue; // is Pi0? | |
458 | ||
459 | if( iParticle->GetNDaughters() == 2 && | |
460 | fStack->Particle(iParticle->GetFirstDaughter())->GetPdgCode() == ::kGamma && | |
461 | fStack->Particle(iParticle->GetLastDaughter())->GetPdgCode() == ::kGamma ) | |
462 | { | |
463 | fHistograms->FillTable("Table_Generation", 1); // pi0 -> gg | |
464 | } | |
465 | ||
466 | if ( iParticle->GetNDaughters() != 3 ) continue; // Num == 3 (e+,e-,gamma) | |
467 | ||
468 | // Check for Pi0 Dalitz decay | |
469 | TParticle* eposPi0 = 0; | |
470 | TParticle* enegPi0 = 0; | |
471 | TParticle* gammaPi0 = 0; | |
472 | ||
473 | for( Int_t idxPi0 = iParticle->GetFirstDaughter(); idxPi0 <= iParticle->GetLastDaughter(); idxPi0++ ) | |
474 | { | |
475 | switch(fStack->Particle(idxPi0)->GetPdgCode()) | |
476 | { | |
477 | case ::kPositron: | |
478 | eposPi0 = fStack->Particle(idxPi0); | |
479 | break; | |
480 | case ::kElectron: | |
481 | enegPi0 = fStack->Particle(idxPi0); | |
482 | break; | |
483 | case ::kGamma: | |
484 | gammaPi0 = fStack->Particle(idxPi0); | |
485 | break; | |
486 | } | |
487 | } | |
488 | ||
489 | if (eposPi0==0 || enegPi0==0 || gammaPi0==0) continue; | |
490 | ||
491 | // found a Pi0 Dalitz decay | |
492 | ||
493 | fHistograms->FillTable("Table_Generation", 2); | |
494 | fHistograms->FillHistogram("MC_Pi0Dalitz_P", iParticle->P()); | |
495 | fHistograms->FillHistogram("MC_Pi0Dalitz_Pt", iParticle->Pt()); | |
496 | fHistograms->FillHistogram("MC_Pi0Dalitz_Eta", iParticle->Eta()); | |
497 | fHistograms->FillHistogram("MC_Pi0Dalitz_Pt_vs_Y", Rapidity(iParticle), iParticle->Pt()); | |
498 | fHistograms->FillHistogram("MC_EposDalitz_Pt", eposPi0->Pt()); | |
499 | fHistograms->FillHistogram("MC_EposDalitz_Eta", eposPi0->Eta()); | |
500 | fHistograms->FillHistogram("MC_EnegDalitz_Pt", enegPi0->Pt()); | |
501 | fHistograms->FillHistogram("MC_EnegDalitz_Eta", enegPi0->Eta()); | |
502 | fHistograms->FillHistogram("MC_GammaPi0Dalitz_Pt", gammaPi0->Pt()); | |
503 | fHistograms->FillHistogram("MC_GammaPi0Dalitz_Eta", gammaPi0->Eta()); | |
504 | ||
505 | // Angle between the gamma and the plane e+e- | |
506 | TVector3 ePosMom( eposPi0->Px(), eposPi0->Py(), eposPi0->Pz() ); | |
507 | TVector3 eNegMom( enegPi0->Px(), enegPi0->Py(), enegPi0->Pz() ); | |
508 | TVector3 gamMom( gammaPi0->Px(), gammaPi0->Py() , gammaPi0->Pz() ); | |
509 | TVector3 planeEposEneg = eNegMom.Cross( ePosMom ); | |
510 | Double_t anglePlaneGamma = planeEposEneg.Angle(gamMom); | |
511 | ||
512 | fHistograms->FillHistogram("MC_EposEnegDalitz_Angle", ePosMom.Angle(eNegMom) ); | |
513 | ||
514 | fHistograms->FillHistogram("MC_EposEnegDalitz_GammaPi0_Angle", anglePlaneGamma); | |
515 | fHistograms->FillHistogram("MC_EposEnegDalitz_GammaPi0_Angle_vs_P", anglePlaneGamma, gammaPi0->P()); | |
516 | fHistograms->FillHistogram("MC_EposEnegDalitz_GammaPi0_Angle_vs_Pt", anglePlaneGamma, gammaPi0->Pt()); | |
517 | ||
518 | ||
519 | // check for gamma conversion | |
520 | Bool_t daugGammaElectron = kFALSE; | |
521 | Bool_t daugGammaPositron = kFALSE; // acceptance | |
522 | Bool_t daugGammaElectronAll = kFALSE; | |
523 | Bool_t daugGammaPositronAll = kFALSE; | |
524 | ||
525 | // is the gamma converted? -> has 2 daughter e+e- | |
526 | // are e+ e- from gamma in the acceptance for the V0s | |
527 | if( gammaPi0->GetNDaughters() >= 2 ) | |
528 | { | |
529 | for( Int_t tIndex=gammaPi0->GetFirstDaughter(); tIndex<=gammaPi0->GetLastDaughter(); ++tIndex ) | |
530 | { | |
531 | TParticle* tmpDaughter = fStack->Particle(tIndex); | |
7811d09e | 532 | |
7f3c7cc6 | 533 | if( tmpDaughter->GetUniqueID() != kPPair ) continue; // check if the daughters come from a conversion |
534 | if( tmpDaughter->GetPdgCode() == ::kElectron ) | |
535 | { // e+ | |
536 | daugGammaElectronAll = kTRUE; | |
7811d09e | 537 | |
7f3c7cc6 | 538 | if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() && |
539 | ((TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue()) < tmpDaughter->R() && | |
540 | (tmpDaughter->R()< fV0Reader->GetMaxRCut() ) ) | |
541 | { | |
542 | daugGammaElectron = kTRUE; | |
543 | } | |
544 | } | |
545 | else if( tmpDaughter->GetPdgCode() == ::kPositron ) | |
546 | { | |
547 | daugGammaPositronAll = kTRUE; | |
548 | if( TMath::Abs(tmpDaughter->Eta()) <= fV0Reader->GetEtaCut() && | |
549 | ((TMath::Abs(tmpDaughter->Vz()) * fV0Reader->GetLineCutZRSlope()) - fV0Reader->GetLineCutZValue()) < tmpDaughter->R() && | |
550 | (tmpDaughter->R() < fV0Reader->GetMaxRCut() ) ) | |
551 | { | |
552 | daugGammaPositron = kTRUE; | |
553 | } | |
554 | } | |
555 | } | |
556 | } | |
557 | ||
558 | ||
559 | if( daugGammaElectronAll && daugGammaPositronAll ) | |
560 | { | |
561 | TParticle* tmpDaughter = fStack->Particle( gammaPi0->GetFirstDaughter() ); | |
562 | fHistograms->FillHistogram("MC_GC_GammaPi0Dalitz_All_Z_vs_R",tmpDaughter->Vz(),tmpDaughter->R() ); | |
563 | } | |
564 | ||
565 | Float_t etaMin, etaMax; | |
566 | fESDtrackCuts->GetEtaRange( etaMin, etaMax ); | |
567 | ||
568 | // e+e- pair within acceptance | |
569 | if ( TMath::Abs( eposPi0->Eta() ) < etaMax && TMath::Abs( enegPi0->Eta() ) < etaMax ) | |
570 | { | |
571 | fHistograms->FillHistogram("MC_Acceptance_EposDalitz_Pt", eposPi0->Pt()); | |
572 | fHistograms->FillHistogram("MC_Acceptance_EposDalitz_Eta", eposPi0->Eta()); | |
573 | fHistograms->FillHistogram("MC_Acceptance_EnegDalitz_Pt", enegPi0->Pt()); | |
574 | fHistograms->FillHistogram("MC_Acceptance_EnegDalitz_Eta", enegPi0->Eta()); | |
575 | fHistograms->FillHistogram("MC_Acceptance_DalitzPair_EposPt_vs_EnegPt", eposPi0->Pt(), enegPi0->Pt()); | |
576 | } | |
577 | ||
578 | // Pi0 (e+e-gamma) within acceptance | |
579 | //cout<<"Gamma Eta Cut"<<fV0Reader->GetEtaCut()<<endl; | |
580 | ||
581 | if ( ( TMath::Abs( gammaPi0->Eta() ) < fV0Reader->GetEtaCut() && gammaPi0->R() < fV0Reader->GetMaxRCut() ) && | |
582 | TMath::Abs( eposPi0->Eta() ) < etaMax && TMath::Abs( enegPi0->Eta() ) < etaMax ) | |
583 | { | |
584 | fHistograms->FillTable("Table_Generation",3); // | |
585 | fHistograms->FillHistogram("MC_Acceptance_Pi0Dalitz_Pt",iParticle->Pt()); | |
586 | fHistograms->FillHistogram("MC_Acceptance_Pi0Dalitz_Eta",iParticle->Eta()); | |
587 | fHistograms->FillHistogram("MC_Acceptance_Pi0Dalitz_Radius",iParticle->R()); | |
588 | fHistograms->FillHistogram("MC_Acceptance_GammaPi0Dalitz_Pt",gammaPi0->Pt()); | |
589 | fHistograms->FillHistogram("MC_Acceptance_GammaPi0Dalitz_Eta",gammaPi0->Eta()); | |
7811d09e | 590 | fHistograms->FillHistogram("MC_Acceptance_EposPi0Dalitz_Pt",eposPi0->Pt()); |
591 | fHistograms->FillHistogram("MC_Acceptance_EposPi0Dalitz_Eta",eposPi0->Eta()); | |
592 | fHistograms->FillHistogram("MC_Acceptance_EnegPi0Dalitz_Pt",enegPi0->Pt()); | |
593 | fHistograms->FillHistogram("MC_Acceptance_EnegPi0Dalitz_Eta",enegPi0->Eta()); | |
7f3c7cc6 | 594 | fHistograms->FillHistogram("MC_Acceptance_DalitzPair_OpeningAngle", ePosMom.Angle(eNegMom) ); |
595 | fHistograms->FillHistogram("MC_Acceptance_Pi0Dalitz_Pt_vs_Y", Rapidity(iParticle), iParticle->Pt()); | |
596 | ||
597 | // Pi0 within acceptance with gamma converted | |
598 | ||
599 | if ( daugGammaElectron && daugGammaPositron ) | |
600 | { | |
601 | fHistograms->FillTable("Table_Generation",4); // | |
602 | ||
603 | fHistograms->FillHistogram("MC_Acceptance_GC_Pi0Dalitz_Pt",iParticle->Pt()); | |
604 | fHistograms->FillHistogram("MC_Acceptance_GC_Pi0Dalitz_Eta",iParticle->Eta()); | |
7811d09e | 605 | fHistograms->FillHistogram("MC_Acceptance_GC_EposPi0Dalitz_Pt",eposPi0->Pt()); |
606 | fHistograms->FillHistogram("MC_Acceptance_GC_EposPi0Dalitz_Eta",eposPi0->Eta()); | |
607 | fHistograms->FillHistogram("MC_Acceptance_GC_EnegPi0Dalitz_Pt",enegPi0->Pt()); | |
608 | fHistograms->FillHistogram("MC_Acceptance_GC_EnegPi0Dalitz_Eta",enegPi0->Eta()); | |
7f3c7cc6 | 609 | fHistograms->FillHistogram("MC_Acceptance_GC_GammaPi0Dalitz_Pt",gammaPi0->Pt()); |
610 | fHistograms->FillHistogram("MC_Acceptance_GC_GammaPi0Dalitz_Eta",gammaPi0->Eta()); | |
611 | //fHistograms->FillHistogram("MC_Acceptance_GC_Gamma_Angle",anglePlaneGamma); | |
612 | //fHistograms->FillHistogram("MC_Acceptance_GC_Gamma_Angle_vs_Pt",anglePlaneGamma,gammaPi0->Pt()); | |
613 | TParticle* tmpDaughter = fStack->Particle( gammaPi0->GetFirstDaughter() ); | |
614 | fHistograms->FillHistogram("MC_Acceptance_GC_GammaPi0Dalitz_Z_vs_R",tmpDaughter->Vz(),tmpDaughter->R() ); | |
615 | fHistograms->FillHistogram("MC_Acceptance_GC_Pi0Dalitz_Pt_vs_Y", Rapidity(iParticle), iParticle->Pt()); | |
616 | } | |
617 | } | |
618 | } | |
619 | } | |
620 | ||
621 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 622 | void AliAnalysisTaskGammaConvDalitz::CreateListOfDalitzPairCandidates() |
623 | { | |
7811d09e | 624 | // |
625 | // Dalitz pair candidates | |
626 | // | |
7f3c7cc6 | 627 | if( fDebug ) AliInfo("=> CreateListOfDalitzPairCandidates"); |
628 | ||
629 | fEposCandidateIndex.clear(); | |
630 | fEnegCandidateIndex.clear(); | |
631 | ||
632 | fHistograms->FillTable("Table_Cuts", 0); | |
633 | ||
634 | for( Int_t i = 0; i < fESDEvent->GetNumberOfTracks(); ++i ) | |
635 | { | |
636 | AliESDtrack* iTrack = fESDEvent->GetTrack(i); | |
637 | if ( !iTrack ) continue; | |
638 | ||
639 | ||
640 | Double_t p[3]; | |
641 | ||
642 | if ( !iTrack->GetConstrainedPxPyPz(p) ) continue; | |
643 | ||
644 | TVector3 iMom(p[0],p[1],p[2]); | |
645 | ||
646 | // | |
647 | // Check track cuts and find track type | |
648 | // | |
649 | ||
650 | Bool_t isTrackAccepted = 0; | |
651 | Int_t trackType = -1; | |
652 | switch(fTrkSelectionCriteria) | |
653 | { | |
654 | case kITSsaTrack: | |
655 | isTrackAccepted = fITSsaTrackCuts->AcceptTrack( iTrack ); | |
656 | trackType = kITSsaTrack; | |
657 | break; | |
658 | ||
659 | case kGlobalTrack: | |
660 | isTrackAccepted = fESDtrackCuts->AcceptTrack( iTrack ); | |
661 | trackType = kGlobalTrack; | |
662 | break; | |
663 | ||
664 | case kITSsaGlobalTrack: | |
665 | if(fITSsaTrackCuts->AcceptTrack( iTrack ) || fESDtrackCuts->AcceptTrack( iTrack )) | |
666 | { | |
667 | isTrackAccepted = kTRUE; | |
668 | if(fITSsaTrackCuts->AcceptTrack( iTrack )) trackType = kITSsaTrack; | |
669 | else trackType = kGlobalTrack; | |
670 | } | |
671 | break; | |
672 | } | |
673 | ||
674 | if(!isTrackAccepted) continue; | |
675 | ||
676 | // | |
677 | // PID | |
678 | // | |
679 | ||
680 | Int_t pid=-1; | |
681 | Int_t pidMC=-1; | |
682 | ||
683 | if(fUseBayesPID) | |
684 | { | |
685 | pid = GetBayesPid(iTrack,trackType); | |
686 | } | |
687 | else | |
688 | { | |
689 | pid = GetNSigmaPid(iTrack,trackType); | |
690 | } | |
691 | ||
692 | if( fDoMC ) | |
693 | { | |
694 | pidMC = GetMonteCarloPid(iTrack); | |
695 | // pid table | |
696 | Int_t iLabel = TMath::Abs(iTrack->GetLabel()); | |
697 | TParticle* iParticle = fStack->Particle(iLabel); | |
698 | FillPidTable(iParticle, pid); | |
699 | } | |
7811d09e | 700 | |
7f3c7cc6 | 701 | // ITS standalone tracks |
702 | if( trackType == kITSsaTrack) | |
703 | { | |
704 | Double_t mom = iTrack->GetP(); | |
705 | Double_t signal = iTrack->GetITSsignal(); | |
706 | ||
707 | fHistograms->FillHistogram( "ESD_ITSsa_dEdx_vs_P", mom, signal ); | |
708 | ||
709 | if( pid == AliPID::kElectron ) | |
710 | { | |
711 | ||
712 | fHistograms->FillHistogram( "ESD_ITSsa_PidCut_dEdx_vs_P", mom, signal ); | |
713 | if(fDoMC && pid == pidMC) | |
714 | { | |
715 | fHistograms->FillHistogram( "MC_ESD_ITSsa_PidCut_dEdx_vs_P", mom, signal ); | |
716 | } | |
717 | } | |
718 | ||
719 | if( fDoMC && pidMC == AliPID::kElectron) | |
720 | { | |
721 | fHistograms->FillHistogram( "MC_ESD_ITSsa_Electron_dEdx_vs_P", mom, signal ); | |
722 | } | |
723 | } | |
724 | ||
725 | else // global tracks | |
726 | { | |
727 | const AliExternalTrackParam *in = iTrack->GetInnerParam(); | |
728 | Double_t mom = in->GetP(); | |
729 | Double_t signal = iTrack->GetTPCsignal(); | |
730 | ||
731 | fHistograms->FillHistogram( "ESD_TPC_dEdx_vs_P", mom, signal ); | |
732 | ||
733 | if( fDoMC && pidMC == AliPID::kElectron ) | |
734 | { | |
735 | fHistograms->FillHistogram( "MC_ESD_TPC_Electron_dEdx_vs_P", mom, signal ); | |
736 | } | |
737 | ||
738 | if( pid == AliPID::kElectron ) | |
739 | { | |
740 | fHistograms->FillHistogram( "ESD_TPC_PidCut_dEdx_vs_P", mom, signal ); | |
741 | if(fDoMC && pid == pidMC) | |
742 | { | |
743 | fHistograms->FillHistogram( "MC_ESD_TPC_PidCut_dEdx_vs_P", mom, signal ); | |
744 | } | |
745 | } | |
746 | } | |
747 | ||
748 | if( AliPID::kElectron != pid) continue; | |
749 | ||
750 | // electron track candidates from here | |
751 | ||
752 | if( iTrack->GetSign() > 0 ) | |
753 | { | |
754 | fEposCandidateIndex.push_back(i); | |
755 | } | |
756 | else | |
757 | { | |
758 | fEnegCandidateIndex.push_back(i); | |
759 | } | |
760 | } | |
761 | ||
762 | // gamma candidates | |
763 | GetGammaCandidates(fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
764 | ||
765 | if(fDoMC) | |
766 | { | |
767 | TClonesArray* pi0Dalitz = FindPi0Dalitz(fEposCandidateIndex, fEnegCandidateIndex, fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
768 | ((TH1*)fHistograms->GetValue("Table_Cuts"))->Fill(1.,(Double_t)pi0Dalitz->GetEntriesFast()); | |
769 | delete pi0Dalitz; | |
770 | } | |
771 | ||
772 | if(fUseTrackIndexCut) // remove repeated tracks | |
773 | { | |
774 | ESDtrackIndexCut(fEposCandidateIndex,fEnegCandidateIndex, fGammaCandidates); | |
775 | ||
776 | if(fDoMC) | |
777 | { | |
778 | TClonesArray* pi0Dalitz = FindPi0Dalitz(fEposCandidateIndex, fEnegCandidateIndex, fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
779 | ((TH1*)fHistograms->GetValue("Table_Cuts"))->Fill(2.,(Double_t)pi0Dalitz->GetEntriesFast()); | |
780 | delete pi0Dalitz; | |
781 | } | |
782 | } | |
783 | ||
784 | if(fUsePsiPairCut) // remove electrons from gamma conversions | |
785 | { | |
786 | PsiPairCut(fEposCandidateIndex,fEnegCandidateIndex); | |
787 | ||
788 | if(fDoMC) | |
789 | { | |
790 | TClonesArray* pi0Dalitz = FindPi0Dalitz(fEposCandidateIndex, fEnegCandidateIndex, fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
791 | ((TH1*)fHistograms->GetValue("Table_Cuts"))->Fill(3.,(Double_t)pi0Dalitz->GetEntriesFast()); | |
792 | delete pi0Dalitz; | |
793 | } | |
794 | } | |
795 | ||
796 | if( fUseMassCut ) | |
797 | { | |
798 | MassCut(fEposCandidateIndex, fEnegCandidateIndex); | |
799 | ||
800 | if(fDoMC) | |
801 | { | |
802 | TClonesArray* pi0Dalitz = FindPi0Dalitz(fEposCandidateIndex, fEnegCandidateIndex, fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
803 | ((TH1*)fHistograms->GetValue("Table_Cuts"))->Fill(4.,(Double_t)pi0Dalitz->GetEntriesFast()); | |
804 | delete pi0Dalitz; | |
805 | } | |
806 | } | |
807 | ||
808 | if(fUseGammaCut) | |
809 | { | |
810 | AngleEposEnegGammaCut(fEposCandidateIndex,fEnegCandidateIndex,fV0Reader->GetCurrentEventGoodV0s(), fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
811 | ||
812 | if(fDoMC) | |
813 | { | |
814 | TClonesArray* pi0Dalitz = FindPi0Dalitz(fEposCandidateIndex, fEnegCandidateIndex, fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
815 | ((TH1*)fHistograms->GetValue("Table_Cuts"))->Fill(5.,(Double_t)pi0Dalitz->GetEntriesFast()); | |
816 | delete pi0Dalitz; | |
817 | } | |
818 | } | |
819 | } | |
820 | ||
821 | //----------------------------------------------------------------------------------------------- | |
822 | void AliAnalysisTaskGammaConvDalitz::ProcessGammaElectronsForDalitzAnalysis() | |
823 | { | |
7811d09e | 824 | // |
825 | // Process gamma and electrons for pi0 Dalitz decay | |
826 | // | |
7f3c7cc6 | 827 | if( fDebug ) AliInfo("=> ProcessGammaElectronsForDalitzAnalysis"); |
828 | ||
829 | fHistograms->FillTable( "Table_Reconstruction", 0); // number of events | |
830 | ||
831 | TClonesArray* ePosCandidates = IndexToAliKFParticle(fEposCandidateIndex, ::kPositron); | |
832 | ||
833 | for(Int_t i=0; i < ePosCandidates->GetEntriesFast(); ++i) | |
834 | { | |
835 | AliKFParticle* epos = (AliKFParticle*) ePosCandidates->At(i); | |
836 | fHistograms->FillHistogram("ESD_EposCandidates_Pt", epos->GetPt()); | |
837 | fHistograms->FillHistogram("ESD_EposCandidates_Eta", epos->GetEta()); | |
838 | fHistograms->FillTable( "Table_Reconstruction", 1); | |
839 | } | |
840 | ||
841 | TClonesArray* eNegCandidates = IndexToAliKFParticle(fEnegCandidateIndex, ::kElectron); | |
842 | ||
843 | for(Int_t i=0; i < eNegCandidates->GetEntriesFast(); ++i) | |
844 | { | |
845 | AliKFParticle* eneg = (AliKFParticle*) eNegCandidates->At(i); | |
846 | fHistograms->FillHistogram("ESD_EnegCandidates_Pt", eneg->GetPt()); | |
847 | fHistograms->FillHistogram("ESD_EnegCandidates_Eta", eneg->GetEta()); | |
848 | fHistograms->FillTable( "Table_Reconstruction", 2); | |
849 | } | |
850 | ||
851 | TClonesArray* dalitzPairCandidates = FindDalitzPair(ePosCandidates, eNegCandidates); | |
852 | for(Int_t i=0; i < dalitzPairCandidates->GetEntriesFast(); ++i) | |
853 | { | |
854 | TLorentzVector* dalitz = (TLorentzVector*)dalitzPairCandidates->At(i); | |
855 | ||
856 | fHistograms->FillHistogram("ESD_DalitzPairCandidates_Pt", dalitz->Pt()); | |
857 | fHistograms->FillHistogram("ESD_DalitzPairCandidates_Mass", dalitz->M()); | |
858 | } | |
859 | ||
860 | // gamma candidates | |
861 | for(Int_t i=0; i < fGammaCandidates->GetEntriesFast(); ++i) | |
862 | { | |
863 | AliKFParticle* gamma = (AliKFParticle*) fGammaCandidates->At(i); | |
864 | fHistograms->FillHistogram("ESD_GammaCandidates_Pt", gamma->GetPt()); | |
865 | fHistograms->FillHistogram("ESD_GammaCandidates_Eta", gamma->GetEta()); | |
866 | } | |
867 | ||
868 | // psi pair for all candidates | |
869 | //if(fUsePsiPairCut) | |
870 | FillPsiPair(ePosCandidates,eNegCandidates,"ESD_EposEneg_PsiPair_vs_DPhi"); | |
871 | ||
872 | // Angle epos,eneg gamma | |
873 | FillAngle(ePosCandidates, fGammaCandidates, "ESD_EposEneg_GammaCandidates_Angle"); | |
874 | FillAngle(eNegCandidates, fGammaCandidates, "ESD_EposEneg_GammaCandidates_Angle"); | |
875 | ||
876 | TClonesArray* pi0Candidates = FindPi0Dalitz(ePosCandidates, eNegCandidates, fGammaCandidates); | |
877 | for(Int_t i=0; i < pi0Candidates->GetEntriesFast(); ++i) | |
878 | { | |
879 | TLorentzVector* pi0 = (TLorentzVector*)pi0Candidates->At(i); | |
880 | ||
881 | fHistograms->FillHistogram("ESD_Pi0_P", pi0->P()); | |
882 | fHistograms->FillHistogram("ESD_Pi0_Pt", pi0->Pt()); | |
883 | fHistograms->FillHistogram("ESD_Pi0_Eta", pi0->Eta()); | |
884 | fHistograms->FillHistogram("ESD_Pi0_Y", pi0->Rapidity()); | |
885 | fHistograms->FillHistogram("ESD_Pi0_Phi", pi0->Phi()); | |
886 | fHistograms->FillHistogram("ESD_Pi0_Pt_vs_Y", pi0->Rapidity(), pi0->Pt()); | |
887 | fHistograms->FillHistogram("ESD_Pi0_Mass", pi0->M()); | |
888 | fHistograms->FillHistogram("ESD_Pi0_Mass_vs_Pt", pi0->Pt(), pi0->M()); | |
889 | fHistograms->FillHistogram("ESD_Pi0_Mass_vs_Y", pi0->Rapidity(),pi0->M()); | |
890 | fHistograms->FillHistogram("ESD_Pi0_Mass_vs_Eta", pi0->Eta(), pi0->M()); | |
891 | } | |
892 | ||
893 | delete dalitzPairCandidates; | |
894 | delete pi0Candidates; | |
895 | ||
896 | if(fComputeBkg) | |
897 | { | |
898 | // 1) e+e- with with gammas used in the signal | |
899 | TClonesArray* pi0Bkg = FindPi0Dalitz(ePosCandidates, eNegCandidates, fGammaPool); | |
900 | ||
901 | for(Int_t i=0; i < pi0Bkg->GetEntriesFast(); ++i) | |
902 | { | |
903 | TLorentzVector* pi0 = (TLorentzVector*)pi0Bkg->At(i); | |
904 | fHistograms->FillHistogram("ESD_BKG_PrevGamma", pi0->M()); | |
905 | fHistograms->FillHistogram("ESD_BKG_PrevGamma_vs_Pt", pi0->Pt(), pi0->M()); | |
906 | } | |
907 | ||
908 | if(ePosCandidates->GetEntriesFast() > 0 && | |
909 | eNegCandidates->GetEntriesFast() > 0 && | |
910 | fGammaCandidates->GetEntriesFast() > 0) | |
911 | { | |
912 | UpdateGammaPool(fGammaCandidates); | |
913 | } | |
914 | ||
915 | delete pi0Bkg; | |
916 | ||
917 | // 2) e+e- with gammas from a pool of events | |
7811d09e | 918 | TClonesArray* gammaBGHandler = GammasFromBGHandler(); |
7f3c7cc6 | 919 | pi0Bkg = FindPi0Dalitz(ePosCandidates, eNegCandidates, gammaBGHandler); |
920 | ||
921 | for(Int_t i=0; i < pi0Bkg->GetEntriesFast(); ++i) | |
922 | { | |
923 | TLorentzVector* pi0 = (TLorentzVector*)pi0Bkg->At(i); | |
924 | fHistograms->FillHistogram("ESD_BKG_BGHandler", pi0->M()); | |
925 | fHistograms->FillHistogram("ESD_BKG_BGHandler_vs_Pt", pi0->Pt(), pi0->M()); | |
926 | } | |
927 | ||
928 | delete pi0Bkg; | |
929 | ||
930 | // 3) e+ with e-, gamma from a pool of events | |
931 | TClonesArray* elecBGHandler = ElectronFromBGHandler(); | |
932 | pi0Bkg = FindPi0Dalitz(ePosCandidates, elecBGHandler, gammaBGHandler); | |
933 | ||
934 | for(Int_t i=0; i < pi0Bkg->GetEntriesFast(); ++i) | |
935 | { | |
936 | TLorentzVector* pi0 = (TLorentzVector*)pi0Bkg->At(i); | |
937 | fHistograms->FillHistogram("ESD_BKG_Electron", pi0->M()); | |
938 | fHistograms->FillHistogram("ESD_BKG_Electron_vs_Pt", pi0->Pt(), pi0->M()); | |
939 | } | |
940 | ||
941 | if(eNegCandidates->GetEntriesFast() > 0) | |
942 | { | |
943 | UpdateElectronPool(eNegCandidates); | |
944 | } | |
945 | ||
946 | delete gammaBGHandler; | |
947 | delete elecBGHandler; | |
948 | delete pi0Bkg; | |
7811d09e | 949 | |
7f3c7cc6 | 950 | } |
951 | ||
952 | delete ePosCandidates; | |
953 | delete eNegCandidates; | |
954 | ||
955 | if(fDoMC) | |
956 | { | |
957 | TClonesArray* ePosPi0Dalitz = FindElectronFromPi0Dalitz(fEposCandidateIndex, ::kPositron); | |
958 | for(Int_t i=0; i < ePosPi0Dalitz->GetEntriesFast(); ++i) | |
959 | { | |
960 | AliKFParticle* epos = (AliKFParticle*) ePosPi0Dalitz->At(i); | |
961 | fHistograms->FillHistogram("MC_ESD_EposDalitz_Pt", epos->GetPt()); | |
962 | fHistograms->FillHistogram("MC_ESD_EposDalitz_Eta", epos->GetEta()); | |
963 | fHistograms->FillTable( "Table_Reconstruction", 3); | |
964 | } | |
965 | ||
966 | TClonesArray* eNegPi0Dalitz = FindElectronFromPi0Dalitz(fEnegCandidateIndex, ::kElectron); | |
967 | for(Int_t i=0; i < eNegPi0Dalitz->GetEntriesFast(); ++i) | |
968 | { | |
969 | AliKFParticle* eneg = (AliKFParticle*) eNegPi0Dalitz->At(i); | |
970 | fHistograms->FillHistogram("MC_ESD_EnegDalitz_Pt", eneg->GetPt()); | |
971 | fHistograms->FillHistogram("MC_ESD_EnegDalitz_Eta", eneg->GetEta()); | |
972 | fHistograms->FillTable( "Table_Reconstruction", 4); | |
973 | } | |
974 | ||
975 | TClonesArray* dalitzPair = FindDalitzPair(fEposCandidateIndex, fEnegCandidateIndex); | |
976 | for(Int_t i=0; i < dalitzPair->GetEntriesFast(); ++i) | |
977 | { | |
978 | TLorentzVector* dalitz = (TLorentzVector*) dalitzPair->At(i); | |
979 | fHistograms->FillHistogram("MC_ESD_DalitzPair_Pt", dalitz->Pt()); | |
980 | fHistograms->FillHistogram("MC_ESD_DalitzPair_Mass", dalitz->M()); | |
981 | fHistograms->FillHistogram( "Table_Reconstruction", 5 ); | |
982 | } | |
983 | ||
984 | // psi pair for dalitz pairs | |
985 | //if(fUsePsiPairCut) | |
986 | FillPsiPair(ePosPi0Dalitz,eNegPi0Dalitz,"MC_ESD_DalitzPair_PsiPair_vs_DPhi"); | |
987 | ||
988 | delete ePosPi0Dalitz; | |
989 | delete eNegPi0Dalitz; | |
990 | delete dalitzPair; | |
991 | ||
992 | // all gammas | |
993 | TClonesArray* gamma = FindGamma(fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
994 | for(Int_t i=0; i < gamma->GetEntriesFast(); ++i) | |
995 | { | |
996 | AliKFParticle* iGamma = (AliKFParticle*) gamma->At(i); | |
997 | fHistograms->FillHistogram("MC_ESD_Gamma_Pt", iGamma->GetPt()); | |
998 | fHistograms->FillHistogram("MC_ESD_Gamma_Eta", iGamma->GetEta()); | |
999 | } | |
1000 | ||
1001 | delete gamma; | |
1002 | ||
1003 | // gamma from pi0 dalitz | |
1004 | TClonesArray* gammaPi0Dalitz = FindGammaFromPi0Dalitz(fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
1005 | for(Int_t i=0; i < gammaPi0Dalitz->GetEntriesFast(); ++i) | |
1006 | { | |
1007 | AliKFParticle* iGamma = (AliKFParticle*) gammaPi0Dalitz->At(i); | |
1008 | fHistograms->FillHistogram("MC_ESD_GammaPi0Dalitz_Pt", iGamma->GetPt()); | |
1009 | fHistograms->FillHistogram("MC_ESD_GammaPi0Dalitz_Eta", iGamma->GetEta()); | |
1010 | fHistograms->FillTable( "Table_Reconstruction", 6); | |
1011 | } | |
1012 | ||
1013 | delete gammaPi0Dalitz; | |
1014 | ||
1015 | TClonesArray* pi0Dalitz = FindPi0Dalitz(fEposCandidateIndex, fEnegCandidateIndex, fGammaCandidates, fGammaCandidatePosIndex, fGammaCandidateNegIndex); | |
1016 | ||
1017 | for(Int_t i=0; i < pi0Dalitz->GetEntriesFast(); ++i) | |
1018 | { | |
1019 | TLorentzVector* pi0 = (TLorentzVector*) pi0Dalitz->At(i); | |
1020 | ||
1021 | fHistograms->FillHistogram("MC_ESD_Pi0_P", pi0->P()); | |
1022 | fHistograms->FillHistogram("MC_ESD_Pi0_Pt", pi0->Pt()); | |
1023 | fHistograms->FillHistogram("MC_ESD_Pi0_Eta", pi0->Eta()); | |
1024 | fHistograms->FillHistogram("MC_ESD_Pi0_Y", pi0->Rapidity()); | |
1025 | fHistograms->FillHistogram("MC_ESD_Pi0_Phi", pi0->Phi()); | |
1026 | fHistograms->FillHistogram("MC_ESD_Pi0_Pt_vs_Y", pi0->Rapidity(), pi0->Pt()); | |
1027 | fHistograms->FillHistogram("MC_ESD_Pi0_Mass", pi0->M()); | |
1028 | fHistograms->FillHistogram("MC_ESD_Pi0_Mass_vs_Pt", pi0->Pt(), pi0->M()); | |
1029 | fHistograms->FillHistogram("MC_ESD_Pi0_Mass_vs_Y", pi0->Rapidity(),pi0->M()); | |
1030 | fHistograms->FillHistogram("MC_ESD_Pi0_Mass_vs_Eta", pi0->Eta(), pi0->M()); | |
1031 | fHistograms->FillHistogram( "Table_Reconstruction", 7); | |
1032 | } | |
1033 | delete pi0Dalitz; | |
1034 | ||
1035 | // psi pair for electrons from gamma conversions assuming they came from main vertex | |
1036 | // if(fUsePsiPairCut) | |
1037 | for(UInt_t i=0; i < fEposCandidateIndex.size(); ++i) | |
1038 | { | |
1039 | AliESDtrack* posTrack = fESDEvent->GetTrack(fEposCandidateIndex[i]); | |
1040 | Int_t posLabel = TMath::Abs(posTrack->GetLabel()); | |
1041 | ||
1042 | for(UInt_t j=0; j < fEnegCandidateIndex.size(); ++j) | |
1043 | { | |
1044 | AliESDtrack* negTrack = fESDEvent->GetTrack(fEnegCandidateIndex[j]); | |
1045 | Int_t negLabel = TMath::Abs(negTrack->GetLabel()); | |
1046 | ||
1047 | if(!IsFromGammaConversion(posLabel,negLabel)) continue; | |
1048 | ||
1049 | Double_t psiPair = GetPsiPair(posTrack, negTrack); | |
1050 | Double_t deltaPhi = fMagFieldSign * TVector2::Phi_mpi_pi( negTrack->GetConstrainedParam()->Phi()-posTrack->GetConstrainedParam()->Phi()); | |
1051 | ||
1052 | fHistograms->FillHistogram("MC_ESD_EposEnegGamma_PsiPair_vs_DPhi", deltaPhi, psiPair); | |
1053 | } | |
1054 | } | |
1055 | // FIXME: eta -> e+e-gamma | |
1056 | } | |
1057 | } | |
1058 | ||
1059 | //-------------------------------------------------------------------------- | |
1060 | Double_t AliAnalysisTaskGammaConvDalitz::Rapidity(const TParticle* p) const | |
1061 | { | |
7811d09e | 1062 | // |
1063 | // Get rapidity | |
1064 | // | |
1065 | const double kEPSILON=1.e-16; | |
7f3c7cc6 | 1066 | |
7811d09e | 1067 | if(p->Energy() - TMath::Abs(p->Pz()) < kEPSILON ) |
7f3c7cc6 | 1068 | { |
1069 | return 1.e10; | |
1070 | } | |
1071 | return 0.5*TMath::Log( (p->Energy()+p->Pz()) / (p->Energy()-p->Pz()) ); | |
1072 | } | |
1073 | ||
1074 | //-------------------------------------------------------------------------- | |
1075 | void AliAnalysisTaskGammaConvDalitz::FillPsiPair(const TClonesArray* pos, const TClonesArray* neg, const TString& hName) | |
1076 | { | |
7811d09e | 1077 | // |
1078 | // Fill histogram with psipair(pos,neg) | |
1079 | // | |
7f3c7cc6 | 1080 | for(Int_t i=0; i < pos->GetEntriesFast(); ++i ) |
1081 | { | |
1082 | AliKFParticle* posKF = (AliKFParticle*) pos->At(i); | |
1083 | for( Int_t j=0; j < neg->GetEntriesFast(); ++j ) | |
1084 | { | |
1085 | AliKFParticle* negKF = (AliKFParticle*) neg->At(j); | |
1086 | Double_t psiPair = GetPsiPair(posKF, negKF); | |
1087 | Double_t deltaPhi = fMagFieldSign * TVector2::Phi_mpi_pi( negKF->GetPhi() - posKF->GetPhi()); | |
1088 | fHistograms->FillHistogram(hName, deltaPhi, psiPair); | |
1089 | } | |
1090 | } | |
1091 | } | |
1092 | ||
1093 | //-------------------------------------------------------------------------- | |
1094 | void AliAnalysisTaskGammaConvDalitz::FillAngle(const TClonesArray* x, const TClonesArray* y, const TString& hName) | |
1095 | { | |
7811d09e | 1096 | // |
1097 | // Fill histogram with angle(x,y) | |
1098 | // | |
7f3c7cc6 | 1099 | for(Int_t i=0; i < x->GetEntriesFast(); ++i ) |
1100 | { | |
1101 | AliKFParticle* xKF = (AliKFParticle*) x->At(i); | |
1102 | TVector3 xMom(xKF->Px(),xKF->Py(),xKF->Pz()); | |
1103 | for( Int_t j=0; j < y->GetEntriesFast(); ++j ) | |
1104 | { | |
1105 | AliKFParticle* yKF = (AliKFParticle*) y->At(j); | |
1106 | TVector3 yMom(yKF->Px(),yKF->Py(),yKF->Pz()); | |
1107 | fHistograms->FillHistogram(hName, xMom.Angle(yMom)); | |
1108 | } | |
1109 | } | |
1110 | } | |
1111 | ||
1112 | //-------------------------------------------------------------------------- | |
7f3c7cc6 | 1113 | void AliAnalysisTaskGammaConvDalitz::FillPidTable(const TParticle* p, Int_t pid) |
1114 | { | |
7811d09e | 1115 | // |
1116 | // Fill table with pid info | |
1117 | // | |
7f3c7cc6 | 1118 | Int_t iGen=-1; |
1119 | switch(TMath::Abs(p->GetPdgCode())) | |
1120 | { | |
1121 | case ::kElectron: iGen=0; break; | |
1122 | case ::kMuonMinus: iGen=1; break; | |
1123 | case ::kPiPlus: iGen=2; break; | |
1124 | case ::kKPlus: iGen=3; break; | |
1125 | case ::kProton: iGen=4; break; | |
1126 | default: iGen=-1; | |
1127 | } | |
1128 | ||
1129 | int jRec=-1; | |
1130 | if(pid > -1 && pid < 5) jRec = pid; | |
1131 | ||
1132 | if ((iGen > -1) && (jRec > -1)) | |
1133 | { | |
1134 | fHistograms->FillTable("Table_PID", iGen, jRec); | |
1135 | // sum | |
1136 | fHistograms->FillTable("Table_PID", iGen, 5); | |
1137 | fHistograms->FillTable("Table_PID", 5, jRec); | |
1138 | } | |
1139 | } | |
1140 | ||
1141 | //-------------------------------------------------------------------------- | |
1142 | void AliAnalysisTaskGammaConvDalitz::GetGammaCandidates(TClonesArray*& gamma, vector<Int_t>& posIndex, vector<Int_t>& negIndex) | |
1143 | { | |
7811d09e | 1144 | // |
1145 | // Make a copy of gamma candidates from V0reader | |
1146 | // | |
7f3c7cc6 | 1147 | posIndex.clear(); |
1148 | negIndex.clear(); | |
1149 | ||
1150 | if(gamma) delete gamma; | |
1151 | ||
1152 | TClonesArray* gammaV0 = fV0Reader->GetCurrentEventGoodV0s(); | |
1153 | ||
1154 | gamma = new TClonesArray("AliKFParticle", gammaV0->GetEntriesFast()); | |
1155 | gamma->SetOwner(kTRUE); | |
1156 | ||
1157 | // make a copy | |
1158 | for(Int_t i=0; i < gammaV0->GetEntriesFast(); ++i) | |
1159 | { | |
1160 | AliKFParticle* gamKF = (AliKFParticle*)gammaV0->At(i); | |
1161 | new ((*gamma)[i]) AliKFParticle(*gamKF); | |
1162 | posIndex.push_back(fV0Reader->GetPindex(i)); | |
1163 | negIndex.push_back(fV0Reader->GetNindex(i)); | |
1164 | } | |
1165 | } | |
1166 | ||
1167 | //-------------------------------------------------------------------------- | |
1168 | TClonesArray* AliAnalysisTaskGammaConvDalitz::IndexToAliKFParticle(const vector<Int_t>& index, Int_t PDG) | |
1169 | { | |
7811d09e | 1170 | // |
1171 | // Convert track index vector to AliKFParticle array | |
1172 | // | |
7f3c7cc6 | 1173 | TClonesArray* indexKF = new TClonesArray("AliKFParticle",index.size()); |
1174 | indexKF->SetOwner(kTRUE); | |
1175 | ||
1176 | for(UInt_t i = 0; i < index.size(); ++i) | |
1177 | { | |
1178 | AliESDtrack* t = fESDEvent->GetTrack(index[i]); | |
1179 | new((*indexKF)[i]) AliKFParticle(*t->GetConstrainedParam(), PDG); | |
1180 | } | |
1181 | ||
1182 | return indexKF; | |
1183 | } | |
1184 | ||
1185 | //-------------------------------------------------------------------------- | |
1186 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindElectronFromPi0Dalitz(const vector<Int_t>& candidates, Int_t PDG) | |
1187 | { | |
7811d09e | 1188 | // |
1189 | // Find true electrons from pi0 Dalitz decay candidates with MC | |
1190 | // | |
7f3c7cc6 | 1191 | TClonesArray* elec = new TClonesArray("AliKFParticle"); |
1192 | elec->SetOwner(kTRUE); | |
1193 | ||
1194 | for(UInt_t i=0, j=0; i < candidates.size(); ++i) | |
1195 | { | |
1196 | AliESDtrack* track = fESDEvent->GetTrack(candidates[i]); | |
1197 | Int_t trackLabel = TMath::Abs(track->GetLabel()); | |
1198 | ||
1199 | if( fStack->Particle(trackLabel)->GetPdgCode() != PDG ) continue; | |
1200 | if( !IsPi0DalitzDaughter(trackLabel) ) continue; | |
1201 | ||
1202 | new ((*elec)[j++]) AliKFParticle(*track->GetConstrainedParam(), PDG); | |
1203 | } | |
1204 | ||
1205 | return elec; | |
1206 | } | |
1207 | ||
1208 | //-------------------------------------------------------------------------- | |
1209 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindGammaFromPi0Dalitz(const TClonesArray* gamma, const vector<Int_t>& posIdx, const vector<Int_t>& negIdx) | |
1210 | { | |
7811d09e | 1211 | // |
1212 | // Find true gammas from pi0 Dalitz decay candidates with MC | |
1213 | // | |
7f3c7cc6 | 1214 | TClonesArray* gammaPi0 = new TClonesArray("AliKFParticle"); |
1215 | gammaPi0->SetOwner(kTRUE); | |
1216 | ||
1217 | for(Int_t i=0, j=0; i < gamma->GetEntriesFast(); ++i) | |
1218 | { | |
1219 | AliKFParticle* gamKF = (AliKFParticle*)gamma->At(i); | |
1220 | ||
1221 | Int_t labelv1 = TMath::Abs((fESDEvent->GetTrack(posIdx[i]))->GetLabel()); | |
1222 | Int_t labelv2 = TMath::Abs((fESDEvent->GetTrack(negIdx[i]))->GetLabel()); | |
1223 | ||
1224 | if( !HaveSameMother(labelv1,labelv2) ) continue; | |
1225 | ||
1226 | Int_t labelGamma = TMath::Abs(fStack->Particle(labelv1)->GetMother(0)); | |
1227 | ||
1228 | if( fStack->Particle(labelGamma)->GetPdgCode() != ::kGamma ) continue; | |
1229 | ||
1230 | if( !IsPi0DalitzDaughter( labelGamma) ) continue; | |
1231 | ||
1232 | new ((*gammaPi0)[j++]) AliKFParticle(*gamKF); | |
1233 | } | |
1234 | ||
1235 | return gammaPi0; | |
1236 | } | |
1237 | ||
1238 | //-------------------------------------------------------------------------- | |
1239 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindGamma(const TClonesArray* gamma, const vector<Int_t>& posIdx, const vector<Int_t>& negIdx) | |
1240 | { | |
7811d09e | 1241 | // |
1242 | // Find true gammas from gamma candidates with MC | |
1243 | // | |
7f3c7cc6 | 1244 | TClonesArray* gammaConv = new TClonesArray("AliKFParticle"); |
1245 | gammaConv->SetOwner(kTRUE); | |
1246 | ||
1247 | for(Int_t i=0, j=0; i < gamma->GetEntriesFast(); ++i) | |
1248 | { | |
1249 | AliKFParticle* gamKF = (AliKFParticle*)gamma->At(i); | |
1250 | ||
1251 | Int_t labelv1 = TMath::Abs((fESDEvent->GetTrack(posIdx[i]))->GetLabel()); | |
1252 | Int_t labelv2 = TMath::Abs((fESDEvent->GetTrack(negIdx[i]))->GetLabel()); | |
1253 | ||
1254 | if( !HaveSameMother(labelv1,labelv2) ) continue; | |
1255 | ||
1256 | Int_t labelGamma = TMath::Abs(fStack->Particle(labelv1)->GetMother(0)); | |
1257 | ||
1258 | if( fStack->Particle(labelGamma)->GetPdgCode() != ::kGamma ) continue; | |
1259 | ||
1260 | new ((*gammaConv)[j++]) AliKFParticle(*gamKF); | |
1261 | } | |
1262 | ||
1263 | return gammaConv; | |
1264 | } | |
1265 | ||
1266 | //-------------------------------------------------------------- | |
7f3c7cc6 | 1267 | void AliAnalysisTaskGammaConvDalitz::ESDtrackIndexCut(vector<Int_t>& pos, vector<Int_t>& neg, const TClonesArray* gamma) |
1268 | { | |
7811d09e | 1269 | // |
1270 | // Remove repeated electron candidate tracks | |
1271 | // according to the gamma candidate array | |
1272 | // | |
7f3c7cc6 | 1273 | vector<Bool_t> posTag(pos.size(),kTRUE); |
1274 | vector<Bool_t> negTag(neg.size(),kTRUE); | |
1275 | ||
1276 | for(Int_t i=0; i < gamma->GetEntriesFast(); ++i) | |
1277 | { | |
1278 | Int_t gamPosIndex = fGammaCandidatePosIndex[i]; | |
1279 | Int_t gamNegIndex = fGammaCandidateNegIndex[i]; | |
1280 | ||
1281 | for( UInt_t j=0; j < pos.size(); ++j ) | |
1282 | { | |
1283 | if(pos[j] == gamPosIndex || pos[j] == gamNegIndex) posTag[j] = kFALSE; | |
1284 | } | |
1285 | ||
1286 | for( UInt_t j=0; j < neg.size(); ++j ) | |
1287 | { | |
1288 | if(neg[j] == gamPosIndex || neg[j] == gamNegIndex) negTag[j] = kFALSE; | |
1289 | } | |
1290 | } | |
1291 | ||
1292 | CleanArray(pos, posTag); | |
1293 | CleanArray(neg, negTag); | |
1294 | } | |
1295 | ||
1296 | //-------------------------------------------------------------------------- | |
7f3c7cc6 | 1297 | void AliAnalysisTaskGammaConvDalitz::PsiPairCut(vector<Int_t>& pos, vector<Int_t>& neg) |
1298 | { | |
7811d09e | 1299 | // |
1300 | // Remove electron candidates from gamma conversions | |
1301 | // according to the Psi pair angle | |
1302 | // | |
7f3c7cc6 | 1303 | vector<Bool_t> posTag(pos.size(), kTRUE); |
1304 | vector<Bool_t> negTag(neg.size(), kTRUE); | |
1305 | ||
1306 | for( UInt_t i=0; i < pos.size(); ++i ) | |
1307 | { | |
1308 | AliESDtrack* posTrack = fESDEvent->GetTrack(pos[i]); | |
1309 | ||
1310 | for( UInt_t j=0; j < neg.size(); ++j ) | |
1311 | { | |
1312 | AliESDtrack* negTrack = fESDEvent->GetTrack(neg[j]); | |
1313 | ||
1314 | Double_t psiPair = GetPsiPair(posTrack, negTrack); | |
1315 | Double_t deltaPhi = fMagFieldSign * TVector2::Phi_mpi_pi( negTrack->GetConstrainedParam()->Phi()-posTrack->GetConstrainedParam()->Phi()); | |
1316 | ||
1317 | if(IsFromGammaConversion( psiPair, deltaPhi )) | |
1318 | { | |
1319 | posTag[i] = kFALSE; | |
1320 | negTag[j] = kFALSE; | |
1321 | } | |
1322 | } | |
1323 | } | |
1324 | ||
1325 | CleanArray(pos, posTag); | |
1326 | CleanArray(neg, negTag); | |
1327 | } | |
1328 | ||
1329 | //----------------------------------------------------------------------------------- | |
7f3c7cc6 | 1330 | void AliAnalysisTaskGammaConvDalitz::MassCut(vector<Int_t>& pos, vector<Int_t>& neg) |
1331 | { | |
7811d09e | 1332 | // |
1333 | // Remove electron candidates pairs | |
1334 | // which have mass not in the range (fMassCutMin,fMassCutMax) | |
1335 | // | |
7f3c7cc6 | 1336 | vector<Bool_t> posTag(pos.size(), kTRUE); |
1337 | vector<Bool_t> negTag(neg.size(), kTRUE); | |
1338 | ||
1339 | for( UInt_t i=0; i < pos.size(); ++i ) | |
1340 | { | |
1341 | AliESDtrack* posTrack = fESDEvent->GetTrack(pos[i]); | |
1342 | ||
1343 | Double_t posMom[3]; posTrack->GetConstrainedPxPyPz(posMom); | |
1344 | TLorentzVector posLV; | |
b9295e41 | 1345 | posLV.SetXYZM(posMom[0],posMom[1],posMom[2],fkElectronMass); |
7f3c7cc6 | 1346 | |
1347 | for( UInt_t j=0; j < neg.size(); ++j ) | |
1348 | { | |
1349 | AliESDtrack* negTrack = fESDEvent->GetTrack(neg[j]); | |
1350 | ||
1351 | Double_t negMom[3]; negTrack->GetConstrainedPxPyPz(negMom); | |
1352 | TLorentzVector negLV; | |
b9295e41 | 1353 | negLV.SetXYZM(negMom[0],negMom[1],negMom[2],fkElectronMass); |
7f3c7cc6 | 1354 | |
1355 | TLorentzVector posnegLV = posLV + negLV; | |
1356 | ||
1357 | if( (posnegLV.M() < fMassCutMin) || (posnegLV.M() > fMassCutMax) ) | |
1358 | { | |
1359 | posTag[i] = kFALSE; | |
1360 | negTag[j] = kFALSE; | |
1361 | } | |
1362 | } | |
1363 | } | |
1364 | ||
1365 | CleanArray(pos, posTag); | |
1366 | CleanArray(neg, negTag); | |
1367 | } | |
1368 | ||
1369 | //----------------------------------------------------------------------------------------------- | |
1370 | void AliAnalysisTaskGammaConvDalitz::CleanArray(vector<Int_t>& x, const vector<Bool_t>& tag) | |
1371 | { | |
7811d09e | 1372 | // |
1373 | // Clean the x array according to the tag parameter | |
1374 | // | |
7f3c7cc6 | 1375 | vector<Int_t> tmp; |
1376 | ||
1377 | for(UInt_t i=0; i< x.size(); ++i) | |
1378 | { | |
1379 | if(tag[i]) tmp.push_back(x[i]); | |
1380 | } | |
1381 | ||
1382 | x = tmp; | |
1383 | } | |
1384 | ||
1385 | //-------------------------------------------------------------------------- | |
7f3c7cc6 | 1386 | void AliAnalysisTaskGammaConvDalitz::AngleEposEnegGammaCut( const vector<Int_t>& posIdx, const vector<Int_t>& negIdx, const TClonesArray* candidates, TClonesArray*& gamma, vector<Int_t>& posGamIdx, vector<Int_t>& negGamIdx) |
1387 | { | |
7811d09e | 1388 | // |
1389 | // Remove gamma candidates according to | |
1390 | // the angle between the plane e+,e- and the gamma | |
1391 | // | |
7f3c7cc6 | 1392 | vector<Bool_t> gammaTag(candidates->GetEntriesFast(), kTRUE); |
1393 | ||
1394 | for( UInt_t iPos=0; iPos < posIdx.size(); ++iPos ) | |
1395 | { | |
1396 | AliESDtrack* posTrack = fESDEvent->GetTrack(posIdx[iPos]); | |
1397 | Double_t posMom[3]; posTrack->GetConstrainedPxPyPz(posMom); | |
1398 | TVector3 xMom(posMom[0],posMom[1],posMom[2]); | |
1399 | ||
1400 | for( UInt_t iNeg=0; iNeg < negIdx.size(); ++iNeg ) | |
1401 | { | |
1402 | AliESDtrack* negTrack = fESDEvent->GetTrack(negIdx[iNeg]); | |
1403 | Double_t negMom[3]; negTrack->GetConstrainedPxPyPz(negMom); | |
1404 | TVector3 yMom(negMom[0],negMom[1],negMom[2]); | |
1405 | ||
1406 | // normal vector to x+y- plane | |
1407 | TVector3 planePosNeg = xMom.Cross(yMom); | |
1408 | for(Int_t i=0; i < candidates->GetEntriesFast(); ++i) | |
1409 | { | |
1410 | AliKFParticle* gamKF = (AliKFParticle*)candidates->At(i); | |
1411 | TVector3 gamMom(gamKF->Px(),gamKF->Py(),gamKF->Pz()); | |
1412 | if (planePosNeg.Angle(gamMom) < 1. || planePosNeg.Angle(gamMom) > 2.) | |
1413 | { | |
1414 | gammaTag[i] = kFALSE; | |
1415 | } | |
1416 | } | |
1417 | } | |
1418 | } | |
1419 | ||
1420 | // Rebuild gamma candidates array | |
1421 | ||
1422 | if(gamma) delete gamma; | |
1423 | gamma = new TClonesArray("AliKFParticle"); | |
1424 | gamma->SetOwner(kTRUE); | |
1425 | ||
1426 | posGamIdx.clear(); | |
1427 | negGamIdx.clear(); | |
1428 | ||
1429 | for(Int_t i=0, j=0; i < candidates->GetEntriesFast(); ++i) | |
1430 | { | |
1431 | AliKFParticle* iGamma = (AliKFParticle*)candidates->At(i); | |
1432 | if(gammaTag[i]) | |
1433 | { | |
1434 | new ((*gamma)[j++]) AliKFParticle(*iGamma); | |
1435 | posGamIdx.push_back(fV0Reader->GetPindex(i)); | |
1436 | negGamIdx.push_back(fV0Reader->GetNindex(i)); | |
1437 | } | |
1438 | } | |
1439 | } | |
1440 | ||
1441 | //-------------------------------------------------------------------------- | |
7f3c7cc6 | 1442 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindDalitzPair(const TClonesArray* pos, const TClonesArray* neg) |
1443 | { | |
7811d09e | 1444 | // |
1445 | // Find Dalitz pair candidates | |
1446 | // | |
7f3c7cc6 | 1447 | TClonesArray* dalitz = new TClonesArray("TLorentzVector"); |
1448 | dalitz->SetOwner(kTRUE); | |
1449 | ||
1450 | for( Int_t iPos=0, j=0; iPos < pos->GetEntriesFast(); ++iPos ) | |
1451 | { | |
1452 | AliKFParticle* posKF = (AliKFParticle*)pos->At(iPos); | |
1453 | ||
1454 | TLorentzVector posLV; | |
b9295e41 | 1455 | posLV.SetXYZM(posKF->Px(),posKF->Py(),posKF->Pz(),fkElectronMass); |
7f3c7cc6 | 1456 | |
1457 | for( Int_t iNeg=0; iNeg < neg->GetEntriesFast(); ++iNeg ) | |
1458 | { | |
1459 | AliKFParticle* negKF = (AliKFParticle*)neg->At(iNeg); | |
1460 | ||
1461 | TLorentzVector negLV; | |
b9295e41 | 1462 | negLV.SetXYZM(negKF->Px(),negKF->Py(),negKF->Pz(),fkElectronMass); |
7f3c7cc6 | 1463 | |
1464 | if(fUseAliKF) | |
1465 | { | |
1466 | AliKFParticle posNeg( *posKF, *negKF); | |
1467 | ||
1468 | TLorentzVector posNegLV; | |
1469 | posNegLV.SetXYZM(posNeg.Px(), posNeg.Py(), posNeg.Pz(), posNeg.GetMass()); | |
1470 | new ((*dalitz)[j++]) TLorentzVector(posNegLV); | |
1471 | } | |
1472 | else | |
1473 | { | |
1474 | new ((*dalitz)[j++]) TLorentzVector(posLV + negLV); | |
1475 | } | |
1476 | } | |
1477 | } | |
1478 | ||
1479 | return dalitz; | |
1480 | } | |
1481 | ||
1482 | //-------------------------------------------------------------------------- | |
1483 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindPi0Dalitz(const TClonesArray* pos, const TClonesArray* neg, const TClonesArray* gamma) | |
1484 | { | |
7811d09e | 1485 | // |
1486 | // Find pi0 Dalitz decay candidates | |
1487 | // | |
7f3c7cc6 | 1488 | TClonesArray* pi0 = new TClonesArray("TLorentzVector"); |
1489 | pi0->SetOwner(kTRUE); | |
1490 | ||
1491 | for( Int_t iPos=0, j=0; iPos < pos->GetEntriesFast(); ++iPos ) | |
1492 | { | |
1493 | AliKFParticle* posKF = (AliKFParticle*)pos->At(iPos); | |
1494 | ||
1495 | TLorentzVector posLV; | |
b9295e41 | 1496 | posLV.SetXYZM(posKF->Px(),posKF->Py(),posKF->Pz(),fkElectronMass); |
7f3c7cc6 | 1497 | |
1498 | for( Int_t iNeg=0; iNeg < neg->GetEntriesFast(); ++iNeg ) | |
1499 | { | |
1500 | AliKFParticle* negKF = (AliKFParticle*)neg->At(iNeg); | |
1501 | ||
1502 | TLorentzVector negLV; | |
b9295e41 | 1503 | negLV.SetXYZM(negKF->Px(),negKF->Py(),negKF->Pz(),fkElectronMass); |
7f3c7cc6 | 1504 | |
1505 | for(Int_t iGam=0; iGam < gamma->GetEntriesFast(); ++iGam) | |
1506 | { | |
1507 | AliKFParticle* gamKF = (AliKFParticle*)gamma->At(iGam); | |
1508 | ||
1509 | if(fUseAliKF) | |
1510 | { | |
1511 | AliKFParticle posNegGam( *posKF, *negKF, *gamKF ); | |
1512 | TLorentzVector posNegGamLV; | |
1513 | posNegGamLV.SetXYZM(posNegGam.Px(),posNegGam.Py(),posNegGam.Pz(),posNegGam.GetMass()); | |
1514 | new ((*pi0)[j++]) TLorentzVector(posNegGamLV); | |
1515 | } | |
1516 | else | |
1517 | { | |
1518 | TLorentzVector gamLV; | |
1519 | gamLV.SetXYZM(gamKF->Px(),gamKF->Py(),gamKF->Pz(),0); | |
1520 | ||
1521 | new ((*pi0)[j++]) TLorentzVector(posLV + negLV + gamLV); | |
1522 | } | |
1523 | } | |
1524 | } | |
1525 | } | |
1526 | ||
1527 | return pi0; | |
1528 | } | |
1529 | ||
1530 | //-------------------------------------------------------------------------- | |
7f3c7cc6 | 1531 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindDalitzPair(const vector<Int_t>& posIdx, const vector<Int_t>& negIdx) |
1532 | { | |
7811d09e | 1533 | // |
1534 | // Find true Dalitz pairs from Dalitz pair candidats with MC | |
1535 | // | |
7f3c7cc6 | 1536 | TClonesArray* dalitz = new TClonesArray("TLorentzVector"); |
1537 | dalitz->SetOwner(kTRUE); | |
1538 | ||
1539 | for( UInt_t iPos=0, j=0; iPos < posIdx.size(); ++iPos ) | |
1540 | { | |
1541 | AliESDtrack* posTrack = fESDEvent->GetTrack(posIdx[iPos]); | |
1542 | Double_t posMom[3]; posTrack->GetConstrainedPxPyPz(posMom); | |
1543 | Int_t posLabel = TMath::Abs(posTrack->GetLabel()); | |
1544 | ||
1545 | TLorentzVector posLV; | |
b9295e41 | 1546 | posLV.SetXYZM(posMom[0],posMom[1],posMom[2],fkElectronMass); |
7f3c7cc6 | 1547 | |
1548 | AliKFParticle posKF( *posTrack->GetConstrainedParam(), ::kPositron ); | |
1549 | ||
1550 | for( UInt_t iNeg=0; iNeg < negIdx.size(); ++iNeg ) | |
1551 | { | |
1552 | AliESDtrack* negTrack = fESDEvent->GetTrack(negIdx[iNeg]); | |
1553 | Int_t negLabel = TMath::Abs(negTrack->GetLabel()); | |
1554 | ||
1555 | if(!IsDalitzPair(posLabel,negLabel)) continue; | |
1556 | ||
1557 | if(fUseAliKF) | |
1558 | { | |
1559 | AliKFParticle negKF( *negTrack->GetConstrainedParam(), ::kElectron ); | |
1560 | AliKFParticle posNeg( posKF, negKF); | |
1561 | ||
1562 | TLorentzVector posNegLV; | |
1563 | posNegLV.SetXYZM(posNeg.Px(),posNeg.Py(),posNeg.Pz(),posNeg.GetMass()); | |
1564 | ||
1565 | new ((*dalitz)[j++]) TLorentzVector(posNegLV); | |
1566 | } | |
1567 | else // TLorentzVector | |
1568 | { | |
1569 | Double_t negMom[3]; negTrack->GetConstrainedPxPyPz(negMom); | |
1570 | ||
1571 | TLorentzVector negLV; | |
b9295e41 | 1572 | negLV.SetXYZM(negMom[0],negMom[1],negMom[2],fkElectronMass); |
7f3c7cc6 | 1573 | |
1574 | new ((*dalitz)[j++]) TLorentzVector(posLV + negLV); | |
1575 | } | |
1576 | } | |
1577 | } | |
1578 | ||
1579 | return dalitz; | |
1580 | } | |
1581 | ||
1582 | //-------------------------------------------------------------------------- | |
7f3c7cc6 | 1583 | TClonesArray* AliAnalysisTaskGammaConvDalitz::FindPi0Dalitz(const vector<Int_t>& posIdx, const vector<Int_t>& negIdx, const TClonesArray* gamma, const vector<Int_t>& posGam, const vector<Int_t>& negGam) |
1584 | { | |
7811d09e | 1585 | // |
1586 | // Find true pi0 Dalitz decay from pi0 candidates with MC | |
1587 | // | |
7f3c7cc6 | 1588 | TClonesArray* pi0 = new TClonesArray("TLorentzVector"); |
1589 | pi0->SetOwner(kTRUE); | |
1590 | ||
1591 | for( UInt_t iPos=0, j=0; iPos < posIdx.size(); ++iPos ) | |
1592 | { | |
1593 | AliESDtrack* posTrack = fESDEvent->GetTrack(posIdx[iPos]); | |
1594 | Double_t posMom[3]; posTrack->GetConstrainedPxPyPz(posMom); | |
1595 | Int_t posLabel = TMath::Abs(posTrack->GetLabel()); | |
1596 | ||
1597 | TLorentzVector posLV; | |
b9295e41 | 1598 | posLV.SetXYZM(posMom[0],posMom[1],posMom[2],fkElectronMass); |
7f3c7cc6 | 1599 | |
1600 | AliKFParticle posKF( *posTrack->GetConstrainedParam(), ::kPositron ); | |
1601 | ||
1602 | for( UInt_t iNeg=0; iNeg < negIdx.size(); ++iNeg ) | |
1603 | { | |
1604 | AliESDtrack* negTrack = fESDEvent->GetTrack(negIdx[iNeg]); | |
1605 | Int_t negLabel = TMath::Abs(negTrack->GetLabel()); | |
1606 | ||
1607 | if(!IsDalitzPair(posLabel,negLabel)) continue; | |
1608 | ||
1609 | Double_t negMom[3]; negTrack->GetConstrainedPxPyPz(negMom); | |
1610 | ||
1611 | TLorentzVector negLV; | |
b9295e41 | 1612 | negLV.SetXYZM(negMom[0],negMom[1],negMom[2],fkElectronMass); |
7f3c7cc6 | 1613 | |
1614 | AliKFParticle negKF( *negTrack->GetConstrainedParam(), ::kElectron ); | |
1615 | ||
1616 | for(Int_t iGam=0; iGam < gamma->GetEntriesFast(); ++iGam) | |
1617 | { | |
1618 | AliKFParticle* gamKF = (AliKFParticle*)gamma->At(iGam); | |
1619 | ||
1620 | Int_t labelv1 = TMath::Abs((fESDEvent->GetTrack(posGam[iGam]))->GetLabel()); | |
1621 | Int_t labelv2 = TMath::Abs((fESDEvent->GetTrack(negGam[iGam]))->GetLabel()); | |
1622 | ||
1623 | if( !HaveSameMother(labelv1,labelv2) ) continue; | |
1624 | ||
1625 | Int_t labelGamma = TMath::Abs(fStack->Particle(labelv1)->GetMother(0)); | |
1626 | ||
1627 | if( fStack->Particle(labelGamma)->GetPdgCode() != ::kGamma ) continue; | |
1628 | if( !HaveSameMother(labelGamma, posLabel) ) continue; | |
1629 | ||
1630 | if(fUseAliKF) | |
1631 | { | |
1632 | AliKFParticle posNegGam( posKF, negKF, *gamKF ); | |
1633 | TLorentzVector posNegGamLV; | |
1634 | posNegGamLV.SetXYZM(posNegGam.Px(),posNegGam.Py(),posNegGam.Pz(),posNegGam.GetMass()); | |
1635 | new ((*pi0)[j++]) TLorentzVector(posNegGamLV); | |
1636 | } | |
1637 | else // TLorentzVector | |
1638 | { | |
1639 | TLorentzVector gamLV; | |
1640 | gamLV.SetXYZM(gamKF->Px(),gamKF->Py(),gamKF->Pz(),0); | |
1641 | ||
1642 | new ((*pi0)[j++]) TLorentzVector(posLV + negLV + gamLV); | |
1643 | } | |
1644 | } | |
1645 | } | |
1646 | } | |
1647 | ||
1648 | return pi0; | |
1649 | } | |
1650 | ||
1651 | //----------------------------------------------------------------------------------------------- | |
1652 | void AliAnalysisTaskGammaConvDalitz::UpdateGammaPool(const TClonesArray* gamma) | |
1653 | { | |
7811d09e | 1654 | // |
1655 | // Update gamma event pool for background computation | |
1656 | // | |
7f3c7cc6 | 1657 | if( fDebug ) AliInfo("=> UpdateGammaPool"); |
1658 | ||
1659 | // cycle | |
1660 | for(Int_t j=0; j< gamma->GetEntriesFast(); ++j) | |
1661 | { | |
1662 | if((AliKFParticle*)fGammaPool->At(fGamPoolPos)) delete (AliKFParticle*)fGammaPool->RemoveAt(fGamPoolPos); | |
1663 | new ((*fGammaPool)[fGamPoolPos]) AliKFParticle( *((AliKFParticle*)gamma->At(j))); | |
1664 | ++fGamPoolPos; | |
1665 | if(fGamPoolPos == fPoolMaxSize) | |
1666 | { | |
1667 | fGamPoolPos = 0; | |
1668 | } | |
1669 | } | |
1670 | } | |
1671 | ||
1672 | void AliAnalysisTaskGammaConvDalitz::UpdateElectronPool(TClonesArray* elec) // FIXME: const | |
1673 | { | |
7811d09e | 1674 | // |
1675 | // Update electron event pool for background computation | |
1676 | // | |
1677 | Int_t multiplicity = fV0Reader->CountESDTracks(); | |
1678 | ||
1679 | ||
1680 | fBGEventHandler->AddElectronEvent(elec,fESDEvent->GetPrimaryVertex()->GetZ(),multiplicity); | |
7f3c7cc6 | 1681 | } |
1682 | ||
1683 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 1684 | TClonesArray* AliAnalysisTaskGammaConvDalitz::GammasFromBGHandler() const |
1685 | { | |
7811d09e | 1686 | // |
1687 | // Gamma copy from events with same multiplicity and Z | |
1688 | // | |
7f3c7cc6 | 1689 | if( fDebug ) AliInfo("=> GammasFromBGHandler"); |
7811d09e | 1690 | |
1691 | Int_t zbin = fBGEventHandler->GetZBinIndex(fV0Reader->GetVertexZ()); | |
1692 | Int_t mbin = fBGEventHandler->GetMultiplicityBinIndex(fV0Reader->CountESDTracks()); | |
1693 | ||
1694 | ||
7f3c7cc6 | 1695 | |
1696 | TClonesArray* gammaPool = new TClonesArray("AliKFParticle"); | |
1697 | gammaPool->SetOwner(kTRUE); | |
1698 | ||
1699 | for( Int_t iEventBG=0; iEventBG < fV0Reader->GetNBGEvents(); ++iEventBG ) | |
1700 | { | |
7811d09e | 1701 | AliGammaConversionKFVector* gammaV0s = fBGEventHandler->GetBGGoodV0s(zbin,mbin,iEventBG); |
7f3c7cc6 | 1702 | for( UInt_t i = 0; i < gammaV0s->size(); ++i) |
1703 | { | |
1704 | new ((*gammaPool)[i]) AliKFParticle( *((AliKFParticle*)gammaV0s->at(i)) ); | |
1705 | } | |
1706 | } | |
1707 | ||
1708 | return gammaPool; | |
1709 | } | |
1710 | ||
1711 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 1712 | TClonesArray* AliAnalysisTaskGammaConvDalitz::ElectronFromBGHandler() const |
1713 | { | |
7811d09e | 1714 | // |
1715 | // Electron copy from events with same multiplicity and Z | |
1716 | // | |
7f3c7cc6 | 1717 | if( fDebug ) AliInfo("=> ElectronFromBGHandler"); |
1718 | ||
1719 | TClonesArray* electronPool = new TClonesArray("AliKFParticle"); | |
1720 | electronPool->SetOwner(kTRUE); | |
7811d09e | 1721 | |
1722 | Int_t multiplicity = fV0Reader->CountESDTracks(); | |
1723 | ||
1724 | ||
1725 | ||
7f3c7cc6 | 1726 | |
1727 | for( Int_t iEventBG=0; iEventBG < fV0Reader->GetNBGEvents(); ++iEventBG ) | |
1728 | { | |
7811d09e | 1729 | AliGammaConversionKFVector* electronNeg = fBGEventHandler->GetBGGoodENeg(iEventBG,fESDEvent->GetPrimaryVertex()->GetZ(),multiplicity); |
7f3c7cc6 | 1730 | for (UInt_t i = 0; i < electronNeg->size(); ++i ) |
1731 | { | |
1732 | new ((*electronPool)[i]) AliKFParticle( *((AliKFParticle*)electronNeg->at(i)) ); | |
1733 | } | |
1734 | } | |
1735 | ||
1736 | return electronPool; | |
1737 | } | |
1738 | ||
1739 | //----------------------------------------------------------------------------------------------- | |
1740 | Int_t AliAnalysisTaskGammaConvDalitz::GetMonteCarloPid(const AliESDtrack* t) const | |
1741 | { | |
7811d09e | 1742 | // |
1743 | // Get track pid according to MC | |
1744 | // | |
7f3c7cc6 | 1745 | Int_t label = TMath::Abs(t->GetLabel()); |
1746 | Int_t pdgCode = TMath::Abs(fStack->Particle(label)->GetPdgCode()); | |
1747 | ||
1748 | switch(pdgCode) | |
1749 | { | |
1750 | case ::kElectron: return AliPID::kElectron; | |
1751 | case ::kMuonMinus: return AliPID::kMuon; | |
1752 | case ::kPiPlus: return AliPID::kPion; | |
1753 | case ::kKPlus: return AliPID::kKaon; | |
1754 | case ::kProton: return AliPID::kProton; | |
1755 | } | |
1756 | ||
1757 | return -1; | |
1758 | } | |
1759 | ||
1760 | //----------------------------------------------------------------------------------------------- | |
7811d09e | 1761 | //FIXME PID ITS |
1762 | // NOTE prior should be estimated from data | |
1763 | // NOTE: move to config | |
1764 | ||
7f3c7cc6 | 1765 | Int_t AliAnalysisTaskGammaConvDalitz::GetBayesPid(const AliESDtrack* t, Int_t trackType ) const |
1766 | { | |
7811d09e | 1767 | // |
1768 | // Get track pid according to Bayes' formula | |
1769 | // | |
7f3c7cc6 | 1770 | double priors[AliPID::kSPECIES] = {0.009, 0.01, 0.82, 0.10, 0.05}; |
1771 | Double_t detectoProb[AliPID::kSPECIES]; | |
1772 | ||
1773 | if( trackType == kITSsaTrack ) // ITS standalone pid | |
1774 | { | |
1775 | t->GetITSpid( detectoProb ); | |
1776 | } | |
1777 | else // global track | |
1778 | { | |
1779 | t->GetESDpid( detectoProb ); | |
1780 | } | |
1781 | ||
1782 | AliPID bayesPID( detectoProb ); | |
1783 | return bayesPID.GetMostProbable( priors ); | |
1784 | } | |
1785 | ||
1786 | //----------------------------------------------------------------------------------------------- | |
1787 | Int_t AliAnalysisTaskGammaConvDalitz::GetNSigmaPid(const AliESDtrack* t, Int_t trackType ) const | |
1788 | { | |
7811d09e | 1789 | // |
1790 | // Get track pid according to a n-sigma cut around ITS and/or TPC signals | |
1791 | // | |
7f3c7cc6 | 1792 | if( trackType == kITSsaTrack) // ITS standalone tracks |
1793 | { | |
1794 | Double_t mom = t->GetP(); | |
1795 | ||
1796 | ||
1797 | // ITS Number of sigmas (FIXME: add new fESDpidCuts) | |
1798 | // NOTE there is not AliESDpidCuts::SetITSnSigmaCut yet | |
1799 | Double_t nElecSigma = fESDpid->NumberOfSigmasITS(t, AliPID::kElectron ); | |
1800 | Double_t nPionSigma = fESDpid->NumberOfSigmasITS(t, AliPID::kPion ); | |
1801 | ||
1802 | if( nElecSigma < 4. && nElecSigma > -3. && mom < .2 && nPionSigma < -1.5 ) | |
1803 | { | |
1804 | return AliPID::kElectron; | |
1805 | } | |
1806 | } | |
1807 | else // global track | |
1808 | { | |
1809 | Double_t nElecSigma = fESDpid->NumberOfSigmasTPC(t, AliPID::kElectron ); | |
1810 | Double_t nPionSigma = fESDpid->NumberOfSigmasTPC(t, AliPID::kPion ); | |
1811 | Double_t nKaonSigma = TMath::Abs(fESDpid->NumberOfSigmasTPC(t, AliPID::kKaon )); | |
1812 | Double_t nProtonSigma = TMath::Abs(fESDpid->NumberOfSigmasTPC(t, AliPID::kProton)); | |
1813 | if( nElecSigma > fNSigmaBelowElecTPCbethe && nElecSigma < fNSigmaAboveElecTPCbethe && | |
1814 | nPionSigma > fNSigmaAbovePionTPCbethe && //NOTE mom > 0.5 | |
1815 | nKaonSigma > fNSigmaAboveKaonTPCbethe && | |
1816 | nProtonSigma > fNSigmaAboveProtonTPCbethe ) | |
1817 | { | |
1818 | return AliPID::kElectron; | |
1819 | } | |
1820 | // NOTE: add other particle types | |
1821 | } | |
1822 | ||
1823 | return -1; | |
1824 | } | |
1825 | ||
1826 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 1827 | Bool_t AliAnalysisTaskGammaConvDalitz::IsDalitzPair( Int_t posLabel, Int_t negLabel ) const |
1828 | { | |
7811d09e | 1829 | // |
1830 | // Returns true if the two particles is a Dalitz pair | |
1831 | // | |
7f3c7cc6 | 1832 | if(!HaveSameMother(posLabel, negLabel)) return kFALSE; |
1833 | ||
1834 | TParticle* pos = fStack->Particle( posLabel ); | |
1835 | TParticle* neg = fStack->Particle( negLabel ); | |
1836 | ||
1837 | if( pos->GetPdgCode() != ::kPositron ) return kFALSE; | |
1838 | if( neg->GetPdgCode() != ::kElectron ) return kFALSE; | |
1839 | ||
1840 | //if( pos->GetUniqueID() != ::kPDecay ) return kFALSE; | |
1841 | //if( neg->GetUniqueID() != ::kPDecay ) return kFALSE; | |
1842 | ||
1843 | Int_t motherLabel = pos->GetMother(0); | |
1844 | if( motherLabel < 0 ) return kFALSE; | |
1845 | ||
1846 | TParticle* mother = fStack->Particle( motherLabel ); | |
1847 | ||
1848 | if( mother->GetPdgCode() != ::kPi0 ) return kFALSE; | |
1849 | if( mother->GetNDaughters() != 3) return kFALSE; | |
1850 | // NOTE: one of them must be a photon | |
1851 | ||
1852 | return kTRUE; | |
1853 | } | |
1854 | ||
1855 | //----------------------------------------------------------------------------------------------- | |
7f3c7cc6 | 1856 | Bool_t AliAnalysisTaskGammaConvDalitz::IsPi0DalitzDaughter( Int_t label ) const |
1857 | { | |
7811d09e | 1858 | // |
1859 | // Returns true if the particle comes from Pi0 -> e+ e- gamma | |
1860 | // | |
7f3c7cc6 | 1861 | Bool_t ePlusFlag = kFALSE; |
1862 | Bool_t eMinusFlag = kFALSE; | |
1863 | Bool_t gammaFlag = kFALSE; | |
1864 | ||
1865 | Int_t motherLabel = fStack->Particle( label )->GetMother(0); | |
1866 | ||
1867 | if( motherLabel < 0 ) return kFALSE; | |
1868 | ||
1869 | TParticle* mother = fStack->Particle( motherLabel ); | |
1870 | ||
1871 | if ( mother->GetPdgCode() != ::kPi0 ) return kFALSE; | |
1872 | ||
1873 | if ( mother->GetNDaughters() != 3 ) return kFALSE; | |
1874 | ||
1875 | for( Int_t idx = mother->GetFirstDaughter(); idx <= mother->GetLastDaughter(); ++idx ) | |
1876 | { | |
1877 | switch( fStack->Particle(idx)->GetPdgCode()) | |
1878 | { | |
1879 | case ::kPositron: | |
1880 | ePlusFlag = kTRUE; | |
1881 | break; | |
1882 | case ::kElectron: | |
1883 | eMinusFlag = kTRUE; | |
1884 | break; | |
1885 | case ::kGamma: | |
1886 | gammaFlag = kTRUE; | |
1887 | break; | |
1888 | } | |
1889 | } | |
1890 | ||
1891 | return ( ePlusFlag && eMinusFlag && gammaFlag ); | |
1892 | } | |
1893 | ||
1894 | //-------------------------------------------------------------------------- | |
1895 | Bool_t AliAnalysisTaskGammaConvDalitz::IsFromGammaConversion( Double_t psiPair, Double_t deltaPhi ) const | |
1896 | { | |
7811d09e | 1897 | // |
1898 | // Returns true if it is a gamma conversion according to psi pair value | |
1899 | // | |
7f3c7cc6 | 1900 | return ( (deltaPhi > fDeltaPhiCutMin && deltaPhi < fDeltaPhiCutMax) && |
1901 | TMath::Abs(psiPair) < ( fPsiPairCut - fPsiPairCut/fDeltaPhiCutMax * deltaPhi ) ); | |
1902 | } | |
1903 | ||
1904 | //-------------------------------------------------------------------------- | |
1905 | Bool_t AliAnalysisTaskGammaConvDalitz::IsFromGammaConversion( Int_t posLabel, Int_t negLabel ) const | |
1906 | { | |
7811d09e | 1907 | // |
1908 | // Returns true if it is a gamma conversion according to MC | |
1909 | // | |
7f3c7cc6 | 1910 | if( !HaveSameMother(posLabel,negLabel) ) return kFALSE; |
1911 | ||
1912 | TParticle* pos = fStack->Particle( posLabel ); | |
1913 | TParticle* neg = fStack->Particle( negLabel ); | |
1914 | ||
1915 | if( pos->GetPdgCode() != ::kPositron ) return kFALSE; | |
1916 | if( neg->GetPdgCode() != ::kElectron ) return kFALSE; | |
1917 | ||
1918 | if( pos->GetUniqueID() != kPPair ) return kFALSE; | |
1919 | ||
1920 | Int_t motherLabel = pos->GetMother(0); | |
1921 | if( motherLabel < 0 ) return kFALSE; | |
1922 | ||
1923 | TParticle* mother = fStack->Particle( motherLabel ); | |
1924 | ||
1925 | if( mother->GetPdgCode() != ::kGamma ) return kFALSE; | |
1926 | ||
1927 | return kTRUE; | |
1928 | } | |
1929 | ||
1930 | //----------------------------------------------------------------------------------------------- | |
1931 | Bool_t AliAnalysisTaskGammaConvDalitz::HaveSameMother( Int_t label1, Int_t label2 ) const | |
1932 | { | |
7811d09e | 1933 | // |
1934 | // Returns true if the two particle have the same mother | |
1935 | // | |
7f3c7cc6 | 1936 | if(fStack->Particle( label1 )->GetMother(0) < 0 ) return kFALSE; |
1937 | return (fStack->Particle( label1 )->GetMother(0) == fStack->Particle( label2 )->GetMother(0)); | |
1938 | } | |
1939 | ||
1940 | //----------------------------------------------------------------------------------------------- | |
1941 | Double_t AliAnalysisTaskGammaConvDalitz::GetPsiPair( const AliESDtrack* trackPos, const AliESDtrack* trackNeg ) const | |
1942 | { | |
7811d09e | 1943 | // |
1944 | // This angle is a measure for the contribution of the opening in polar | |
1945 | // direction Δ0 to the opening angle ξ Pair | |
1946 | // | |
1947 | // Ref. Measurement of photons via conversion pairs with the PHENIX experiment at RHIC | |
1948 | // Master Thesis. Thorsten Dahms. 2005 | |
1949 | // https://twiki.cern.ch/twiki/pub/ALICE/GammaPhysicsPublications/tdahms_thesis.pdf | |
1950 | // | |
b9295e41 | 1951 | Double_t momPos[3]; |
1952 | Double_t momNeg[3]; | |
1953 | if( trackPos->GetConstrainedPxPyPz(momPos) == 0 ) trackPos->GetPxPyPz( momPos ); | |
1954 | if( trackNeg->GetConstrainedPxPyPz(momNeg) == 0 ) trackNeg->GetPxPyPz( momNeg ); | |
7f3c7cc6 | 1955 | |
b9295e41 | 1956 | TVector3 posDaughter; |
1957 | TVector3 negDaughter; | |
7f3c7cc6 | 1958 | |
b9295e41 | 1959 | posDaughter.SetXYZ( momPos[0], momPos[1], momPos[2] ); |
1960 | negDaughter.SetXYZ( momNeg[0], momNeg[1], momNeg[2] ); | |
7f3c7cc6 | 1961 | |
b9295e41 | 1962 | Double_t deltaTheta = negDaughter.Theta() - posDaughter.Theta(); |
1963 | Double_t openingAngle = posDaughter.Angle( negDaughter ); //TMath::ACos( posDaughter.Dot(negDaughter)/(negDaughter.Mag()*posDaughter.Mag()) ); | |
7f3c7cc6 | 1964 | if( openingAngle < 1e-20 ) return 0.; |
1965 | Double_t psiAngle = TMath::ASin( deltaTheta/openingAngle ); | |
1966 | ||
1967 | return psiAngle; | |
1968 | } | |
1969 | ||
1970 | //----------------------------------------------------------------------------------------------- | |
7811d09e | 1971 | Double_t AliAnalysisTaskGammaConvDalitz::GetPsiPair(const AliKFParticle* xPos, const AliKFParticle* yNeg ) const |
7f3c7cc6 | 1972 | { |
7811d09e | 1973 | // |
1974 | // Get psi pair value | |
1975 | // | |
b9295e41 | 1976 | TVector3 pos(xPos->GetPx(), xPos->GetPy(), xPos->GetPz()); |
1977 | TVector3 neg(yNeg->GetPx(), yNeg->GetPy(), yNeg->GetPz()); | |
7f3c7cc6 | 1978 | |
1979 | Double_t deltaTheta = neg.Theta() - pos.Theta(); | |
1980 | Double_t openingAngle = pos.Angle( neg ); | |
1981 | ||
1982 | if( openingAngle < 1e-20 ) return 0.; | |
1983 | ||
1984 | return TMath::ASin( deltaTheta/openingAngle ); | |
1985 | } | |
1986 | ||
1987 | //----------------------------------------------------------------------------------------------- | |
1988 | Double_t AliAnalysisTaskGammaConvDalitz::GetPsiPair(const TLorentzVector* xPos, const TLorentzVector* yNeg ) const | |
1989 | { | |
7811d09e | 1990 | // |
1991 | // Get psi pair value | |
1992 | // | |
7f3c7cc6 | 1993 | Double_t deltaTheta = yNeg->Theta() - xPos->Theta(); |
1994 | Double_t openingAngle = xPos->Angle( yNeg->Vect() ); | |
1995 | ||
1996 | if( openingAngle < 1e-20 ) return 0.; | |
1997 | ||
1998 | return TMath::ASin( deltaTheta/openingAngle );; | |
1999 | } | |
2000 | ||
2001 | // tmp NOTE: Should go to AliV0Reader | |
2002 | //----------------------------------------------------------------------------------------------- | |
2003 | /* | |
2004 | Double_t AliAnalysisTaskGammaConvDalitz::GetPsiAngleV0( | |
2005 | Double_t radiusVO, // radius at XY of VO vertex | |
2006 | const AliExternalTrackParam* trackPos, // pos. track parm. at V0 vertex | |
2007 | const AliExternalTrackParam* trackNeg // neg. track parm. at V0 vertex | |
2008 | ) | |
2009 | { | |
2010 | // This angle is a measure for the contribution of the opening in polar | |
2011 | // direction Δ0 to the opening angle ξ Pair | |
2012 | ||
2013 | // Ref. Measurement of photons via conversion pairs with the PHENIX experiment at RHIC | |
2014 | // Master Thesis. Thorsten Dahms. 2005 | |
2015 | // https://twiki.cern.ch/twiki/pub/ALICE/GammaPhysicsPublications/tdahms_thesis.pdf | |
2016 | ||
2017 | const Double_t kForceDeltaPhi = 1.2; | |
2018 | static const Double_t kB2C = 0.299792458e-3; // Taken from AliVParticle | |
2019 | ||
2020 | Double_t psiAngle = 0.; // Default value | |
2021 | ||
2022 | static Double_t MagnField = fESDEvent->GetMagneticField(); | |
2023 | static Double_t MagnFieldG = fESDEvent->GetMagneticField()*kB2C; | |
2024 | ||
2025 | if( TMath::Abs( MagnField ) < 1e-20 ) return psiAngle; // Do nothing if 0 mag field | |
2026 | ||
2027 | // compute propagation radius for a fixed angle | |
2028 | Double_t Rpos = radiusVO + trackPos->Pt()*TMath::Sin( kForceDeltaPhi ) / MagnFieldG; | |
2029 | Double_t Rneg = radiusVO + trackNeg->Pt()*TMath::Sin( kForceDeltaPhi ) / MagnFieldG; | |
2030 | ||
7811d09e | 2031 | Double_t MomPos[3]; |
2032 | Double_t MomNeg[3]; | |
2033 | if( trackPos->GetPxPyPzAt( Rpos, MagnField, MomPos ) == 0 ) trackPos->GetPxPyPz( MomPos ); | |
2034 | if( trackNeg->GetPxPyPzAt( Rneg, MagnField, MomNeg ) == 0 ) trackNeg->GetPxPyPz( MomNeg ); | |
7f3c7cc6 | 2035 | |
7811d09e | 2036 | TVector3 PosDaughter; |
2037 | TVector3 NegDaughter; | |
2038 | PosDaughter.SetXYZ( MomPos[0], MomPos[1], MomPos[2] ); | |
2039 | NegDaughter.SetXYZ( MomNeg[0], MomNeg[1], MomNeg[2] ); | |
7f3c7cc6 | 2040 | |
7811d09e | 2041 | Double_t deltaTheta = NegDaughter.Theta() - PosDaughter.Theta(); |
2042 | Double_t chiPar = PosDaughter.Angle( NegDaughter ); //TMath::ACos( PosDaughter.Dot(NegDaughter) / (NegDaughter.Mag()*PosDaughter.Mag()) ); | |
7f3c7cc6 | 2043 | if( chiPar < 1e-20 ) return psiAngle; |
2044 | ||
2045 | psiAngle = TMath::ASin( deltaTheta / chiPar ); | |
2046 | ||
2047 | return psiAngle; | |
2048 | } | |
2049 | */ |