1 /*************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 // Utility class for V0 PID
17 // Identifies Electrons, Pions and Protons using gamma conversions and
18 // the decays of K0s and Lambdas
19 // Containers with samples of Electrons, Pions and Protons can be accessed
20 // via GetListOfElectrons() etc.
23 // Matus Kalisky <matus.kalisky@cern.ch>
24 // Markus Heide <mheide@uni-muenster.de>
25 // Markus Fasel <M.Fasel@gsi.de>
27 #include <TObjArray.h>
29 #include "AliAODEvent.h"
31 #include "AliESDEvent.h"
32 #include "AliESDtrack.h"
34 #include "AliKFVertex.h"
35 #include "AliVEvent.h"
36 #include "AliVTrack.h"
37 #include "AliMCParticle.h"
39 #include "AliMCEvent.h"
41 #include "AliHFEV0cuts.h"
42 #include "AliHFEV0info.h"
43 #include "AliHFEcollection.h"
45 #include "AliHFEV0pid.h"
48 AliHFEV0pid::AliHFEV0pid():
54 , fPrimaryVertex(NULL)
72 // Default constructor
75 //____________________________________________________________
76 AliHFEV0pid::AliHFEV0pid(const char *name):
82 , fPrimaryVertex(NULL)
99 // Standard constructor
101 fElectrons = new TObjArray();
102 fPionsK0 = new TObjArray();
103 fPionsL = new TObjArray();
104 fKaons = new TObjArray();
105 fProtons = new TObjArray();
107 fElectrons->SetOwner();
108 fPionsK0->SetOwner();
109 fProtons->SetOwner();
113 fGammas = new TObjArray();
114 fK0s = new TObjArray();
115 fLambdas = new TObjArray();
116 fAntiLambdas = new TObjArray();
118 fIndices = new AliHFEV0pidTrackIndex();
121 //____________________________________________________________
123 AliHFEV0pid::~AliHFEV0pid(){
128 if(fElectrons) delete fElectrons;
129 if(fPionsK0) delete fPionsK0;
130 if(fPionsL) delete fPionsL;
131 if(fKaons) delete fKaons;
132 if(fProtons) delete fProtons;
134 if(fGammas) delete fGammas;
135 if(fK0s) delete fK0s;
136 if(fLambdas) delete fLambdas;
137 if(fAntiLambdas) delete fAntiLambdas;
139 if(fIndices) delete fIndices;
140 if(fV0cuts) delete fV0cuts;
142 if(TESTBIT(fDestBits, 1)){
144 if(fOutput) delete fOutput;
148 //____________________________________________________________
149 void AliHFEV0pid::InitQA(){
151 // Initialize QA histograms
154 memset(&fDestBits, 0, sizeof(fDestBits));
155 SETBIT(fDestBits, 1);
157 fOutput = new TList();
160 fV0cuts = new AliHFEV0cuts();
161 fV0cuts->Init("V0cuts");
164 fQA = new AliHFEcollection("v0pidQA", "QA histograms for V0 selection");
166 fQA->CreateTH1F("h_nV0s", "No. of found and accepted V0s", 5, -0.5, 4.5);
168 // QA histograms for invariant mass
169 fQA->CreateTH1F("h_InvMassGamma", "Gamma invariant mass; inv mass [GeV/c^{2}]; counts", 100, 0, 0.25);
170 fQA->CreateTH1F("h_InvMassK0s", "K0s invariant mass; inv mass [GeV/c^{2}]; counts", 200, 0.4, 0.65);
171 fQA->CreateTH1F("h_InvMassL", "Lambda invariant mass; inv mass [GeV/c^{2}]; counts", 100, 1.05, 1.15);
172 fQA->CreateTH1F("h_InvMassAL", "Lambda invariant mass; inv mass [GeV/c^{2}]; counts", 100, 1.05, 1.15);
174 // QA histograms for p distribution (of the daughters)
175 fQA->CreateTH1F("h_P_electron", "P distribution of the gamma electrons; p (GeV/c); counts", 100, 0.1, 10);
176 fQA->CreateTH1F("h_P_K0pion", "P distribution of the K0 pions; p (GeV/c); counts", 100, 0.1, 10);
177 fQA->CreateTH1F("h_P_Lpion", "P distribution of the Lambda pions; p (GeV/c); counts", 100, 0.1, 10);
178 fQA->CreateTH1F("h_P_Lproton", "P distribution of the Lambda protons; p (GeV/c); counts", 100, 0.1, 10);
181 fQA->CreateTH1F("h_Pt_Gamma", "Pt of the gamma conversion; p_{T} (GeV/c); counts", 100, 0, 10);
182 fQA->CreateTH1F("h_Pt_K0", "Pt of the K0; p_{T} (GeV/c); counts", 100, 0, 10);
183 fQA->CreateTH1F("h_Pt_L", "Pt of the Lambda; p_{T} (GeV/c); counts", 100, 0, 10);
184 fQA->CreateTH1F("h_Pt_AL", "Pt of the Lambda; p_{T} (GeV/c); counts", 100, 0, 10);
186 // Armenteros plot V0 preselection
187 fQA->CreateTH2F("h_AP_all_V0s", "armenteros plot for all V0 candidates; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
188 fQA->CreateTH2F("h_AP_selected_V0s", "armenteros plot for all V0 candidates; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
193 fQA->CreateTH2F("h_AP_MC_all_V0s", "armenteros plot for all MC tagged V0s; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
194 fQA->CreateTH2F("h_AP_MC_Gamma", "armenteros plot for MC tagged Gammas; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
195 fQA->CreateTH2F("h_AP_MC_K0", "armenteros plot for MC tagged K0s; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
196 fQA->CreateTH2F("h_AP_MC_Lambda", "armenteros plot for MC tagged Lambdas; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
197 fQA->CreateTH2F("h_AP_MC_ALambda", "armenteros plot for MC tagged A-Lambdass; #alpha; Q_{T}", 200, -1, 1, 200, 0, 0.25);
201 // armenteros plots for different V0 momenta - MC signal only
202 fQA->CreateTH2Fvector1(12, "h_AP_MC_Gamma_p", "armenteros plot for MC tagged Gammas; #alpha; Q_{T}", 200, -1., 1., 200, 0., 0.25);
203 fQA->CreateTH2Fvector1(12, "h_AP_MC_K0_p", "armenteros plot for MC tagged K0s; #alpha; Q_{T}", 200, -1., 1., 200, 0., 0.25);
204 fQA->CreateTH2Fvector1(12, "h_AP_MC_Lambda_p", "armenteros plot for MC tagged Lambdas; #alpha; Q_{T}", 200, -1., 1., 200, 0., 0.25);
211 //____________________________________________________________
212 void AliHFEV0pid::Process(AliVEvent * const inputEvent){
215 // Find protons, pions and electrons using V0 decays and
216 // store the pointers in the TObjArray
221 Int_t nGamma = 0, nK0s = 0, nLambda = 0, nPhi = 0;
222 fInputEvent = inputEvent;
223 fNtracks = fInputEvent->GetNumberOfTracks();
224 fIndices->Init(fInputEvent->GetNumberOfV0s() * 2);
225 fPrimaryVertex = new AliKFVertex(*(fInputEvent->GetPrimaryVertex()));
226 if(!fPrimaryVertex) return;
227 fV0cuts->SetInputEvent(fInputEvent);
228 fV0cuts->SetPrimaryVertex(fPrimaryVertex);
229 if(fMCEvent) fV0cuts->SetMCEvent(fMCEvent);
230 Int_t check[fNtracks];
231 memset(check, 0, sizeof(Int_t)*fNtracks);
234 //BenchmarkV0finder();
236 for(Int_t iv0 = 0; iv0 < fInputEvent->GetNumberOfV0s(); iv0++){
237 if(!TString(fInputEvent->IsA()->GetName()).CompareTo("AliESDEvent")){
240 AliESDv0 *esdV0 = (static_cast<AliESDEvent *>(fInputEvent))->GetV0(iv0);
242 if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
243 v0status = ProcessV0(esdV0);
247 AliAODv0 *aodV0 = (static_cast<AliAODEvent *>(fInputEvent))->GetV0(iv0);
248 if(aodV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
249 v0status = ProcessV0(aodV0);
250 if(AliHFEV0cuts::kUndef != v0status){
254 case AliHFEV0cuts::kRecoGamma: nGamma++; break;
255 case AliHFEV0cuts::kRecoK0: nK0s++; break;
256 case AliHFEV0cuts::kRecoPhi: nPhi++; break;
257 case AliHFEV0cuts::kRecoLambda: nLambda++; break;
261 AliDebug(1, Form("Number of gammas : %d", nGamma));
262 AliDebug(1, Form("Number of K0s : %d", nK0s));
263 AliDebug(1, Form("Number of Phis : %d", nPhi));
264 AliDebug(1, Form("Number of Lambdas : %d", nLambda));
266 AliDebug(1, "Number of stored tracks:");
267 AliDebug(1, Form("Number of electrons : %d", fElectrons->GetEntries()));
268 AliDebug(1, Form("Number of K0 pions : %d", fPionsK0->GetEntries()));
269 AliDebug(1, Form("Number of Lambda pions : %d", fPionsL->GetEntries()));
270 AliDebug(1, Form("Number of Phi kaons : %d", fKaons->GetEntries()));
271 AliDebug(1, Form("Number of protons : %d", fProtons->GetEntries()));
273 delete fPrimaryVertex;
277 //____________________________________________________________
278 Int_t AliHFEV0pid::ProcessV0(TObject *v0){
281 // Apply general cut and special cuts for gamma, K0s, Lambda
283 if(!v0) return AliHFEV0cuts::kUndef;
285 AliESDv0* esdV0 = dynamic_cast<AliESDv0 *>(v0);
287 AliError("Unexpected v0 type.");
288 return AliHFEV0cuts::kUndef;
290 AliVTrack* daughter[2];
291 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(esdV0->GetPindex()));
292 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(esdV0->GetNindex()));
293 if(!daughter[0] || !daughter[1]) return AliHFEV0cuts::kUndef;
295 if(fMCEvent != NULL) fMCon = kTRUE;
296 //printf("-D: fMCEvent %x, fMCon: %i\n", fMCEvent, fMCon);
299 Int_t dMC[2] = {-1, -1};
300 Int_t idMC = AliHFEV0cuts::kUndef;
304 idMC = IdentifyV0(v0, dMC);
305 //printf("--D: V0 pid: %i, P: %i, N: %i\n", id, d[0], d[1]);
306 fV0cuts->SetCurrentV0id(idMC);
307 fV0cuts->SetDaughtersID(dMC);
309 // check common single track cuts
310 for(Int_t i=0; i<2; ++i){
311 if(!fV0cuts->TrackCutsCommon(static_cast<AliESDtrack*>(daughter[i]))) return AliHFEV0cuts::kUndef;
313 // check commom V0 cuts
315 if(!fV0cuts->V0CutsCommon(esdV0)) return AliHFEV0cuts::kUndef;
318 // preselect the V0 candidates based on the Armenteros plot
319 Int_t id = PreselectV0(esdV0, idMC);
321 if(AliHFEV0cuts::kRecoGamma == id && IsGammaConv(v0)){
322 fQA->Fill("h_nV0s", AliHFEV0cuts::kRecoGamma);
323 return AliHFEV0cuts::kRecoGamma;
325 else if(AliHFEV0cuts::kRecoK0 == id && IsK0s(v0)){
326 fQA->Fill("h_nV0s", AliHFEV0cuts::kRecoK0);
327 return AliHFEV0cuts::kRecoK0;
329 else if(AliHFEV0cuts::kRecoLambda == TMath::Abs(id) && IsLambda(v0)){
330 fQA->Fill("h_nV0s", AliHFEV0cuts::kRecoLambda);
331 return AliHFEV0cuts::kRecoLambda;
333 else return AliHFEV0cuts::kUndef;
338 //____________________________________________________________
339 void AliHFEV0pid::Flush(){
343 AliDebug(1, "Flushing containers");
347 fElectrons->Delete();
350 //____________________________________________________________
351 Int_t AliHFEV0pid::PreselectV0(AliESDv0 * const v0, Int_t idMC){
353 // Based on Armenteros plot preselet the possible identity of the V0 candidate
358 // momentum dependent armenteros plots
359 ArmenterosPlotMC(v0, idMC);
361 // comupte the armenteros variables
363 fV0cuts->Armenteros(v0, ar);
365 const Float_t alpha = ar[0];
366 const Float_t qt = ar[1];
367 //printf(" -D: Alpha: %f, QT: %f \n", alpha, qt);
369 if(TMath::Abs(alpha) > 1) return AliHFEV0cuts::kUndef;
371 fQA->Fill("h_AP_all_V0s", alpha, qt);
373 // fill a few MC tagged histograms - AP plots
376 case AliHFEV0cuts::kRecoGamma :
377 fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
378 fQA->Fill("h_AP_MC_Gamma", alpha, qt);
380 case AliHFEV0cuts::kRecoK0 :
381 fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
382 fQA->Fill("h_AP_MC_K0", alpha, qt);
384 case AliHFEV0cuts::kRecoLambda :
385 fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
386 fQA->Fill("h_AP_MC_Lambda", alpha, qt);
388 case AliHFEV0cuts::kRecoALambda :
389 fQA->Fill("h_AP_MC_all_V0s", alpha, qt);
390 fQA->Fill("h_AP_MC_ALambda", alpha, qt);
396 const Double_t cutAlphaG = 0.35;
397 const Double_t cutQTG = 0.05;
398 const Double_t cutAlphaG2[2] = {0.6, 0.8};
399 const Double_t cutQTG2 = 0.04;
402 const Float_t cutQTK0[2] = {0.1075, 0.215};
403 const Float_t cutAPK0[2] = {0.199, 0.8}; // parameters for curved QT cut
405 // Lambda & A-Lambda cuts
406 const Float_t cutQTL = 0.03;
407 const Float_t cutAlphaL[2] = {0.35, 0.7};
408 const Float_t cutAlphaAL[2] = {-0.7, -0.35};
409 const Float_t cutAPL[3] = {0.107, -0.69, 0.5}; // parameters fir curved QT cut
412 // Check for Gamma candidates
414 if( (TMath::Abs(alpha) < cutAlphaG) ){
415 fQA->Fill("h_AP_selected_V0s", alpha, qt);
416 return AliHFEV0cuts::kRecoGamma;
419 // additional region - should help high pT gammas
421 if( (TMath::Abs(alpha) > cutAlphaG2[0]) && (TMath::Abs(alpha) < cutAlphaG2[1]) ){
422 fQA->Fill("h_AP_selected_V0s", alpha, qt);
423 return AliHFEV0cuts::kRecoGamma;
428 // Check for K0 candidates
429 Float_t q = cutAPK0[0] * TMath::Sqrt(TMath::Abs(1 - alpha*alpha/(cutAPK0[1]*cutAPK0[1])));
430 if( (qt > cutQTK0[0]) && (qt < cutQTK0[1]) && (qt > q) ){
431 fQA->Fill("h_AP_selected_V0s", alpha, qt);
432 return AliHFEV0cuts::kRecoK0;
435 if( (alpha > 0) && (alpha > cutAlphaL[0]) && (alpha < cutAlphaL[1]) && (qt > cutQTL)){
436 q = cutAPL[0] * TMath::Sqrt(1 - ( (alpha + cutAPL[1]) * (alpha + cutAPL[1])) / (cutAPL[2]*cutAPL[2]) );
438 fQA->Fill("h_AP_selected_V0s", alpha, qt);
439 return AliHFEV0cuts::kRecoLambda;
443 // Check for A-Lambda candidates
444 if( (alpha < 0) && (alpha > cutAlphaAL[0]) && (alpha < cutAlphaAL[1]) && (qt > cutQTL)){
445 q = cutAPL[0] * TMath::Sqrt(1 - ( (alpha - cutAPL[1]) * (alpha - cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) );
447 fQA->Fill("h_AP_selected_V0s", alpha, qt);
448 return AliHFEV0cuts::kRecoLambda;
452 return AliHFEV0cuts::kUndef;
455 //____________________________________________________________
456 Bool_t AliHFEV0pid::IsGammaConv(TObject *v0){
461 if(!v0) return kFALSE;
463 AliVTrack* daughter[2];
464 Int_t pIndex = 0, nIndex = 0;
465 Double_t invMass = 0.;
470 AliESDv0 *esdV0 = static_cast<AliESDv0 *>(v0);
471 v0id = esdV0->GetLabel();
472 // apply FULL gamma cuts
473 if(!fV0cuts->GammaCuts(esdV0)) return kFALSE;
474 invMass = esdV0->GetEffMass(AliPID::kElectron, AliPID::kElectron);
475 pIndex = esdV0->GetPindex();
476 nIndex = esdV0->GetNindex();
479 // AOD Analysis - not possible to cut
480 AliAODv0 *aodV0 = static_cast<AliAODv0 *>(v0);
481 v0id = aodV0->GetID();
482 pIndex = aodV0->GetPosID();
483 nIndex = aodV0->GetNegID();
484 invMass = aodV0->InvMass2Prongs(0, 1, kElectron, kElectron);
486 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(pIndex));
487 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(nIndex));
488 if(!daughter[0] || !daughter[1]) return kFALSE;
491 AliDebug(1, Form("Gamma identified, daughter IDs: %d,%d", daughter[0]->GetID(), daughter[1]->GetID()));
493 // AFTER all gamma cuts
494 fQA->Fill("h_Pt_Gamma", mPt);
495 fQA->Fill("h_InvMassGamma", invMass);
497 // Identified gamma - store tracks in the electron containers
498 if(!fIndices->Find(daughter[0]->GetID())){
499 AliDebug(1, Form("Gamma identified, daughter IDs: %d,%d", daughter[0]->GetID(), daughter[1]->GetID()));
500 fElectrons->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
501 fIndices->Add(daughter[0]->GetID(), AliHFEV0cuts::kRecoElectron);
503 if(!fIndices->Find(daughter[1]->GetID())){
504 AliDebug(1, Form("Gamma identified, daughter IDs: %d,%d", daughter[1]->GetID(), daughter[1]->GetID()));
505 fElectrons->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
506 fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoElectron);
512 //____________________________________________________________
513 Bool_t AliHFEV0pid::IsK0s(TObject *v0){
518 if(!v0) return kFALSE;
520 AliVTrack* daughter[2];
521 Int_t pIndex = 0, nIndex = 0;
523 Double_t invMass = 0.;
527 AliESDv0 *esdV0 = static_cast<AliESDv0 *>(v0);
528 if(!fV0cuts->K0Cuts(esdV0)) return kFALSE;
529 v0id = esdV0->GetLabel();
530 pIndex = esdV0->GetPindex();
531 nIndex = esdV0->GetNindex();
532 invMass = esdV0->GetEffMass(AliPID::kPion, AliPID::kPion);
535 // AOD Analysis - not possible to cut
536 AliAODv0 *aodV0 = static_cast<AliAODv0 *>(v0);
538 pIndex = aodV0->GetPosID();
539 nIndex = aodV0->GetNegID();
540 invMass = aodV0->MassK0Short();
542 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(pIndex));
543 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(nIndex));
544 if(!daughter[0] || !daughter[1]) return kFALSE;
546 fQA->Fill("h_Pt_K0", mPt);
547 fQA->Fill("h_InvMassK0s", invMass);
549 AliDebug(1, Form("K0 identified, daughter IDs: %d,%d", daughter[0]->GetID(), daughter[1]->GetID()));
552 // Identified gamma - store tracks in the electron containers
553 if(!fIndices->Find(daughter[0]->GetID())){
554 AliDebug(1, Form("Adding K0 Pion track with ID %d", daughter[0]->GetID()));
555 fPionsK0->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
556 fIndices->Add(daughter[0]->GetID(), AliHFEV0cuts::kRecoPionK0);
558 if(!fIndices->Find(daughter[1]->GetID())){
559 AliDebug(1, Form("Adding K0 Pion track with ID %d", daughter[1]->GetID()));
560 fPionsK0->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
561 fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoPionK0);
567 //____________________________________________________________
568 Bool_t AliHFEV0pid::IsPhi(const TObject *v0) const {
570 // Identify Phi - very preliminary - requires diffrent approach (V0 fnder is not effective)
573 //const Double_t kPhiMass=TDatabasePDG::Instance()->GetParticle(333)->Mass(); // PDG phi mass
574 //AliVTrack* daughter[2];
575 //Double_t invMass = 0.;
577 if(!v0) return kFALSE;
579 //Int_t pIndex = 0, nIndex = 0;
582 //AliESDv0 *esdV0 = static_cast<AliESDv0 *>(v0);
583 //pIndex = esdV0->GetPindex();
584 //nIndex = esdV0->GetNindex();
587 // AOD Analysis - not possible to cut
593 //____________________________________________________________
594 Bool_t AliHFEV0pid::IsLambda(TObject *v0){
599 if(!v0) return kFALSE;
601 AliVTrack* daughter[2];
602 Int_t pIndex = 0, nIndex = 0;
603 Double_t invMass = 0.;
604 Bool_t isLambda = kTRUE; // Lambda - kTRUE, Anti Lambda - kFALSE
609 AliESDv0 *esdV0 = static_cast<AliESDv0 *>(v0);
610 v0id = esdV0->GetLabel();
611 if(!fV0cuts->LambdaCuts(esdV0,isLambda)) return kFALSE;
613 if(fV0cuts->CheckSigns(esdV0)){
614 pIndex = esdV0->GetPindex();
615 nIndex = esdV0->GetNindex();
618 pIndex = esdV0->GetNindex();
619 nIndex = esdV0->GetPindex();
623 // AOD Analysis - not possible to cut
625 // again - two cases as above
626 AliAODv0 *aodV0 = static_cast<AliAODv0 *>(v0);
627 v0id = aodV0->GetID();
628 pIndex = aodV0->GetPosID();
629 nIndex = aodV0->GetNegID();
630 invMass = aodV0->MassLambda();
633 daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(pIndex));
634 daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack(nIndex));
635 if(!daughter[0] || !daughter[1]) return kFALSE;
639 fQA->Fill("h_Pt_L", mPt);
640 fQA->Fill("h_InvMassL", invMass);
642 if(!fIndices->Find(daughter[0]->GetID())){
643 fProtons->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
644 fIndices->Add(daughter[0]->GetID(), AliHFEV0cuts::kRecoProton);
646 if(!fIndices->Find(daughter[1]->GetID())){
647 fPionsL->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
648 fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoPionL);
653 fQA->Fill("h_Pt_AL", mPt);
654 fQA->Fill("h_InvMassAL", invMass);
655 if(!fIndices->Find(daughter [1]->GetID())){
656 fProtons->Add(new AliHFEV0info(daughter[1], daughter[0]->GetID(), v0id));
657 fIndices->Add(daughter[1]->GetID(), AliHFEV0cuts::kRecoProton);
659 if(!fIndices->Find(daughter [0]->GetID())){
660 fPionsL->Add(new AliHFEV0info(daughter[0], daughter[1]->GetID(), v0id));
661 fIndices->Add(daughter [0]->GetID(), AliHFEV0cuts::kRecoPionL);
664 if(isLambda) fLambdas->Add(v0);
665 else fAntiLambdas->Add(v0);
670 //____________________________________________________________
671 Int_t AliHFEV0pid::IdentifyV0(TObject *esdV0, Int_t d[2]){
673 // for MC only, returns the V0 Id
677 // be carefull about changing the return values - they are used later selectively
678 // In particulra "-2" means that identity of either of the daughters could not be
682 AliESDv0 *v0 = dynamic_cast<AliESDv0 *>(esdV0);
685 AliESDtrack* dN, *dP;
688 iP = v0->GetPindex();
689 iN = v0->GetNindex();
690 if(iN < 0 || iP < 0) return -1;
691 if(iN >= fNtracks || iP >= fNtracks) return -1;
692 dP = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iP));
693 dN = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iN));
694 if(!dN || !dP) return -1;
697 // there is still a problem with wrong assignment of positive and negative
698 // V0 daughter in V0 finder - a check is necessary
699 // if the V0 daughters are miss-assigned - swap their labels
700 Bool_t sign = fV0cuts->CheckSigns(v0);
712 if(lN < 0 || lP < 0) return -2;
713 // get the associated MC particles
714 AliMCParticle *mcP, *mcN;
715 mcP = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lP));
716 mcN = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lN));
717 if(!mcP || !mcN) return -2;
719 // identify the daughter tracks and their mothers
721 pdgP = TMath::Abs(mcP->PdgCode());
722 pdgN = TMath::Abs(mcN->PdgCode());
723 // store the daughter ID for later use
726 //printf(" -D: pdgP: %i, pdgN: %i\n", pdgP, pdgN);
727 // lablel of the mother particle
728 // -1 may mean it was a primary particle
730 lPm = mcP->GetMother();
731 lNm = mcN->GetMother();
732 if(-1==lPm || -1==lNm) return -3;
734 // return if mothers are not the same particle
735 if(lPm != lNm) return -3;
736 // get MC mother particle - now we need only 1
737 AliMCParticle *m = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lPm));
740 Int_t pdgM = m->PdgCode();
742 // if(3122 == TMath::Abs(pdgM)){
743 // printf("-D: v0 signs: %i\n", fV0cuts->CheckSigns(v0));
744 // printf("-D: pdgM: %i, pdgN: %i, pdgP: %i \n", pdgM, pdgN, pdgP);
748 // now check the mother and daughters identity
749 if(22 == TMath::Abs(pdgM) && 11 == pdgN && 11 == pdgP) return AliHFEV0cuts::kRecoGamma;
750 if(310 == TMath::Abs(pdgM) && 211 == pdgN && 211 == pdgP) return AliHFEV0cuts::kRecoK0;
751 if(-3122 == pdgM && 2212 == pdgN && 211 == pdgP) return AliHFEV0cuts::kRecoALambda;
752 if(3122 == pdgM && 211 == pdgN && 2212 == pdgP) return AliHFEV0cuts::kRecoLambda;
754 return AliHFEV0cuts::kUndef;
757 //____________________________________________________________
758 void AliHFEV0pid::BenchmarkV0finder(){
760 // produce histograms for all findable V0s that are
761 // were selected byt the (oline) V0 finder and can
762 // be used to estimate the efficiency of teh V0 cuts
765 for(Int_t iv0 = 0; iv0 < fInputEvent->GetNumberOfV0s(); iv0++){
766 AliESDv0 *esdV0 = (static_cast<AliESDEvent *>(fInputEvent))->GetV0(iv0);
768 if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
769 // indetify the V0 candidate
770 Int_t idV0 = AliHFEV0cuts::kUndef;
771 Int_t idD[2] = {-1, -1};
772 idV0 = IdentifyV0(esdV0, idD);
775 //____________________________________________________________
776 void AliHFEV0pid::ArmenterosPlotMC(AliESDv0 * const v0, Int_t idMC){
778 // Armenteros plots as a function of Mohter Momentum
780 //const Float_t minP = 0.1;
781 //const Float_t maxP = 10.;
782 // approx log bins - over the 0.1 - 10 GeV/c
783 const Float_t bins[13] = {0.1, 0.1468, 0.2154, 0.3162, 0.4642, 0.6813, 1.0, 1.4678, 2.1544, 3.1623, 4.6416, 6.8129, 10.0};
786 fV0cuts->Armenteros(v0, ar);
789 if( (p <= bins[0]) || (p >= bins[12])) return;
792 Float_t tmp = bins[0];
799 if(AliHFEV0cuts::kRecoGamma == idMC) fQA->Fill("h_AP_MC_Gamma_p", pBin, ar[0], ar[1]);
800 if(AliHFEV0cuts::kRecoK0 == idMC) fQA->Fill("h_AP_MC_K0_p", pBin, ar[0], ar[1]);
801 if(AliHFEV0cuts::kRecoLambda == TMath::Abs(idMC)) fQA->Fill("h_AP_MC_Lambda_p", pBin, ar[0], ar[1]);
806 //____________________________________________________________
807 AliHFEV0pid::AliHFEV0pidTrackIndex::AliHFEV0pidTrackIndex():
813 , fIndexElectron(NULL)
820 // Default Constructor
824 //____________________________________________________________
825 AliHFEV0pid::AliHFEV0pidTrackIndex::~AliHFEV0pidTrackIndex(){
829 if(fIndexElectron) delete[] fIndexElectron;
830 if(fIndexPionK0) delete[] fIndexPionK0;
831 if(fIndexPionL) delete[] fIndexPionL;
832 if(fIndexProton) delete[] fIndexProton;
835 //____________________________________________________________
836 void AliHFEV0pid::AliHFEV0pidTrackIndex::Flush(){
841 if(fIndexElectron) delete[] fIndexElectron;
842 fIndexElectron = NULL;
843 if(fIndexPionK0) delete[] fIndexPionK0;
845 if(fIndexPionL) delete[] fIndexPionL;
847 if(fIndexKaon) delete[] fIndexKaon;
849 if(fIndexProton) delete[] fIndexProton;
859 //____________________________________________________________
860 void AliHFEV0pid::AliHFEV0pidTrackIndex::Init(Int_t capacity){
862 // Initialize container
864 fIndexElectron = new Int_t[capacity];
865 fIndexPionK0 = new Int_t[capacity];
866 fIndexPionL = new Int_t[capacity];
867 fIndexProton = new Int_t[capacity];
870 //____________________________________________________________
871 void AliHFEV0pid::AliHFEV0pidTrackIndex::Add(Int_t index, Int_t species){
873 // Add new index to the list of identified particles
876 case AliHFEV0cuts::kRecoElectron:
877 fIndexElectron[fNElectrons++] = index;
879 case AliHFEV0cuts::kRecoPionK0:
880 fIndexPionK0[fNPionsK0++] = index;
882 case AliHFEV0cuts::kRecoPionL:
883 fIndexPionL[fNPionsL++] = index;
885 case AliHFEV0cuts::kRecoProton:
886 fIndexProton[fNProtons++] = index;
891 //____________________________________________________________
892 Bool_t AliHFEV0pid::AliHFEV0pidTrackIndex::Find(Int_t index, Int_t species) const {
894 // Find track index in the specific sample of particles
897 Int_t *container = NULL; Int_t n = 0;
899 case AliHFEV0cuts::kRecoElectron:
900 container = fIndexElectron;
903 case AliHFEV0cuts::kRecoPionK0:
904 container = fIndexPionK0;
907 case AliHFEV0cuts::kRecoPionL:
908 container = fIndexPionL;
911 case AliHFEV0cuts::kRecoProton:
912 container = fIndexProton;
916 if(!container) return kFALSE;
917 if(n == 0) return kFALSE;
918 Bool_t found = kFALSE;
919 for(Int_t i = 0; i < n; i++){
920 if(container[i] == index){
928 //____________________________________________________________
929 Bool_t AliHFEV0pid::AliHFEV0pidTrackIndex::Find(Int_t index) const {
931 // Find index in all samples
933 if(Find(index, AliHFEV0cuts::kRecoElectron)) return kTRUE;
934 else if(Find(index, AliHFEV0cuts::kRecoPionK0)) return kTRUE;
935 else if(Find(index, AliHFEV0cuts::kRecoPionL)) return kTRUE;
936 else return Find(index, AliHFEV0cuts::kRecoProton);
939 //____________________________________________________________
940 TList *AliHFEV0pid::GetListOfQAhistograms(){
942 // Getter for V0 PID QA histograms
945 CLRBIT(fDestBits, 1);
947 TList *tmp = fV0cuts->GetList();
948 tmp->SetName("V0cuts");
952 tmp = fQA->GetList();
953 tmp->SetName("V0pid");
957 tmp = fV0cuts->GetListMC();
958 tmp->SetName("V0cutsMC");
959 //printf(" -D: adding MC V0 cuts stuff\n");