Cleanup the code. Fix memory leak. Now inherit from AliAnalysisTaskSE (Antoine, Phili...
[u/mrichter/AliRoot.git] / PWG3 / hfe / AliAnalysisTaskCheckV0tenderII.cxx
CommitLineData
3a72645a 1/**************************************************************************
2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercial purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
27de2dfb 15
16/* $Id$ */
17
3a72645a 18//
e3ae862b 19// Task for checking the performance of the V0 tender
3a72645a 20//
21//
6555e2ad 22// Authors:
3a72645a 23// Matus Kalisky <matus.kalisky@cern.ch>
24//
25
26#include <TH1F.h>
27#include <TList.h>
28
29#include "AliAnalysisManager.h"
30#include "AliMCEventHandler.h"
31#include "AliESDInputHandler.h"
32#include "AliESDv0.h"
33#include "AliESDtrack.h"
34#include "AliMCParticle.h"
35#include "AliMCEvent.h"
36#include "AliESDv0KineCuts.h"
37#include "AliKFVertex.h"
38
39#include "AliHFEtools.h"
40#include "AliHFEcollection.h"
41
42#include "AliAnalysisTaskCheckV0tenderII.h"
43
44ClassImp(AliAnalysisTaskCheckV0tenderII)
45
46//__________________________________________________________
47AliAnalysisTaskCheckV0tenderII::AliAnalysisTaskCheckV0tenderII():
48 AliAnalysisTaskSE("CheckV0tenderTask")
49 , fOutput(0x0)
50 , fColl(0x0)
51 , fCollMC(0x0)
52 , fV0cuts(0x0)
53 , fPrimaryVertex(0x0)
e3ae862b 54 , fpdgV0(-1)
55 , fpdgP(-1)
56 , fpdgN(-1)
3a72645a 57 , fEvents(0x0)
58{
59 //
60 // Default Constructor
61 //
62
63
64 DefineOutput(1, TH1F::Class());
65 DefineOutput(2, TList::Class());
66
67
68}
69//__________________________________________________________
70AliAnalysisTaskCheckV0tenderII::AliAnalysisTaskCheckV0tenderII(const Char_t *name):
71 AliAnalysisTaskSE(name)
72 , fOutput(0x0)
73 , fColl(0x0)
74 , fCollMC(0x0)
75 , fV0cuts(0x0)
76 , fPrimaryVertex(0x0)
77 , fpdgV0(0)
78 , fpdgP(0)
79 , fpdgN(0)
80 , fEvents(0x0)
81{
82 //
83 // Constructor
84 //
85
86 DefineOutput(1, TH1F::Class());
87 DefineOutput(2, TList::Class());
88
89}
90//__________________________________________________________
91AliAnalysisTaskCheckV0tenderII::~AliAnalysisTaskCheckV0tenderII(){
92 //
93 // Destructor
94 //
95
96 if (fOutput) delete fOutput;
e3ae862b 97 //if (fColl) delete fColl;
98 //if (fCollMC) delete fCollMC;
3a72645a 99 if (fV0cuts) delete fV0cuts;
100 if (fPrimaryVertex) delete fPrimaryVertex;
101 if (fEvents) delete fEvents;
102}
103//__________________________________________________________
104void AliAnalysisTaskCheckV0tenderII::UserCreateOutputObjects(){
105 //
106 // prepare output objects
107 //
108
109 fOutput = new TList();
110 //fOutput->SetOwner();
111 // Counter for number of events
112 fEvents = new TH1I("nEvents", "NumberOfEvents", 1, 1, 2);
113
114 fColl = new AliHFEcollection("V0_QA", "tender V0s for data");
115 fCollMC = new AliHFEcollection("V0_MC_QA", "tender V0s for MC");
116
117 fV0cuts = new AliESDv0KineCuts();
118 if(!fV0cuts){
119 AliError("Failed to initialize AliESDv0KineCuts instance");
120 return;
121 }
122
123 //
124 // Data histos
125 //
126
127 // total number of tagged V0s
128 fColl->CreateTH1F("h_NumberOf_V0s", "Number of tagged V0s; type; counts", 4, -0.5, 3.5);
129 // pT spectra of the tagged V0s
130 fColl->CreateTH1F("h_Gamma_pt", "p_{T} spectrum of tagged gammas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
131 fColl->CreateTH1F("h_K0_pt", "p_{T} spectrum of tagged K0s; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
132 fColl->CreateTH1F("h_Lambda_pt", "p_{T} spectrum of tagged Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
133 fColl->CreateTH1F("h_ALambda_pt", "p_{T} spectrum of tagged A-Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
134 // invariant mass of the V0s
135 fColl->CreateTH1F("h_Gamma_mass", "Inv. Mass of the ", 100, 0., 0.1);
136 fColl->CreateTH1F("h_K0_mass", "Inv. Mass of the ", 100, 0.45, 0.55);
137 fColl->CreateTH1F("h_Lambda_mass", "Inv. Mass of the ", 100, 1.08, 1.14);
138 fColl->CreateTH1F("h_ALambda_mass", "Inv. Mass of the ", 100, 1.08, 1.14);
139
140 // total number of tagged daughter particles (should correlate with number of V0s !
141 fColl->CreateTH1F("h_NumberOfDaughters", "Number of tagged daughters; type; counts", 3, -0.5, 2.5);
142 // pT spectra of tagged daughter particles
143 fColl->CreateTH1F("h_Electron_pt", "p_{T} spectrum of tagged; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
144 fColl->CreateTH1F("h_Pion_pt", "p_{T} spectrum of tagged; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
145 fColl->CreateTH1F("h_Proton_pt", "p_{T} spectrum of tagged; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
146
147 //
148 // MC histos
149 //
150
151 // pT spectra of the tagged V0s
152 fCollMC->CreateTH1F("h_Gamma_pt_S", "MC-S: p_{T} spectrum of tagged gammas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
153 fCollMC->CreateTH1F("h_K0_pt_S", "MC-S: p_{T} spectrum of tagged K0s; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
154 fCollMC->CreateTH1F("h_Lambda_pt_S", "MC-S: p_{T} spectrum of tagged Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
155 fCollMC->CreateTH1F("h_ALambda_pt_S", "MC-S: p_{T} spectrum of tagged A-Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
156 fCollMC->CreateTH1F("h_Gamma_pt_B", "MC-B: p_{T} spectrum of tagged gammas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
157 fCollMC->CreateTH1F("h_K0_pt_B", "MC-B: p_{T} spectrum of tagged K0s; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
158 fCollMC->CreateTH1F("h_Lambda_pt_B", "MC-B: p_{T} spectrum of tagged Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
159 fCollMC->CreateTH1F("h_ALambda_pt_B", "MC-B: p_{T} spectrum of tagged A-Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
160 // invariant mass of the V0s
161 fCollMC->CreateTH1F("h_Gamma_mass_S", "MC-S: Inv. Mass of the gamma; m (GeV/c^{2}); counts", 100, 0., 0.1);
162 fCollMC->CreateTH1F("h_K0_mass_S", "MC-S: Inv. Mass of the K0; m (GeV/c^{2}); counts", 100, 0.45, 0.55);
163 fCollMC->CreateTH1F("h_Lambda_mass_S", "MC-S: Inv. Mass of the Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
164 fCollMC->CreateTH1F("h_ALambda_mass_S", "MC-S: Inv. Mass of the A-Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
165 fCollMC->CreateTH1F("h_Gamma_mass_B", "MC-B: Inv. Mass of the gamma; m (GeV/c^{2}); counts", 100, 0., 0.1);
166 fCollMC->CreateTH1F("h_K0_mass_B", "MC-B: Inv. Mass of the K0; m (GeV/c^{2}); counts", 100, 0.45, 0.55);
167 fCollMC->CreateTH1F("h_Lambda_mass_B", "MC-B: Inv. Mass of the Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
168 fCollMC->CreateTH1F("h_ALambda_mass_B", "MC-B: Inv. Mass of the A-Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
169 // pT spectra of tagged daughter particles
170 fCollMC->CreateTH1F("h_Electron_pt_S", "MC-S: p_{T} spectrum of tagged electrons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
171 fCollMC->CreateTH1F("h_Pion_pt_S", "MC-S: p_{T} spectrum of tagged pions; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
172 fCollMC->CreateTH1F("h_Proton_pt_S", "MC-S: p_{T} spectrum of tagged protons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
173 fCollMC->CreateTH1F("h_Electron_pt_B", "MC-B: p_{T} spectrum of tagged electrons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
174 fCollMC->CreateTH1F("h_Pion_pt_B", "MC-B: p_{T} spectrum of tagged pions; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
175 fCollMC->CreateTH1F("h_Proton_pt_B", "MC-B: p_{T} spectrum of tagged protons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
e3ae862b 176 // background study - MC only
177 fCollMC->CreateTH1Fvector1(6, "h_Gamma_Bg", "MC gamma bg. - different sources; p_{T} (GeV/c), counts", 20, 0.1, 20, 0);
178 fCollMC->CreateTH1Fvector1(6, "h_K0_Bg", "MC K0 bg. - different sources; p_{T} (GeV/c), counts", 20, 0.1, 20, 0);
179 fCollMC->CreateTH1Fvector1(6, "h_Lambda_Bg", "MC gamma bg. - different sources; p_{T} (GeV/c), counts", 20, 0.1, 20, 0);
180 fCollMC->CreateTH1Fvector1(6, "h_ALambda_Bg", "MC gamma bg. - different sources; p_{T} (GeV/c), counts", 20, 0.1, 20, 0);
181
3a72645a 182
183 TList *tmp = fColl->GetList();
184 tmp->SetName(fColl->GetName());
185 fOutput->Add(tmp);
186 tmp = 0x0;
187 tmp = fCollMC->GetList();
188 tmp->SetName(fCollMC->GetName());
189 fOutput->Add(tmp);
190
191
192}
193//__________________________________________________________
194void AliAnalysisTaskCheckV0tenderII::UserExec(Option_t *){
195 //
196 // Event Loop
197 //
198 AliMCEventHandler* mcHandler = (dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()));
bf892a6a 199
200 //AliESDInputHandler *inh = dynamic_cast<AliESDInputHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
201 //AliESDpid *workingPID = NULL;
202 //if(inh && (workingPID = inh->GetESDpid())) workingPID = inh->GetESDpid();
203 //else workingPID = AliHFEtools::GetDefaultPID(mcHandler ? kTRUE : kFALSE);
3a72645a 204
205 // check the MC data
206 if(fMCEvent && !mcHandler ) return;
207 if(fMCEvent && !mcHandler->InitOk() ) return;
208 if(fMCEvent && !mcHandler->TreeK() ) return;
209 if(fMCEvent && !mcHandler->TreeTR() ) return;
210
e3ae862b 211 AliKFVertex *primaryVertex = new AliKFVertex(*(fInputEvent->GetPrimaryVertex()));
212 if(!primaryVertex) return;
213
214 fV0cuts->SetPrimaryVertex(primaryVertex);
3a72645a 215 fV0cuts->SetEvent(fInputEvent);
e3ae862b 216
3a72645a 217 ProcessV0s();
e3ae862b 218
3a72645a 219
e3ae862b 220 if (primaryVertex) delete primaryVertex;
3a72645a 221 fEvents->Fill(1.1);
222 PostData(1, fEvents);
223 PostData(2, fOutput);
224}
225//__________________________________________________________
226void AliAnalysisTaskCheckV0tenderII::Terminate(Option_t *){
227 //
228 // Do Post Processing
229 //
230}
231//__________________________________________________________
232void AliAnalysisTaskCheckV0tenderII::ProcessV0s(){
233 //
234 // loop over the V0s and extract the information about
235 // the V0s tagged by the V0 tender
236 //
237
e3ae862b 238 const TString type[4] = {"Gamma", "K0", "Lambda", "ALambda"};
239 TString name;
3a72645a 240
241 Int_t nV0s = fInputEvent->GetNumberOfV0s();
242 for(Int_t i=0; i<nV0s; ++i){
bf892a6a 243 AliESDv0 *esdV0 = (static_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
3a72645a 244 if(!esdV0) continue;
245 if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
246
e3ae862b 247 ResetPDGcodes();
248
3a72645a 249 // New standalone V0 selection software
e97c2edf 250 if( ! (fV0cuts->ProcessV0(esdV0, fpdgV0, fpdgP, fpdgN)) ) continue;
3a72645a 251
252 Int_t pid = PDGtoPIDv0(fpdgV0);
253 //printf(" -D: pdg: %i, pid: %i\n", fpdgV0, pid);
254 if(pid <= -1) continue;
255
256 fColl->Fill("h_NumberOf_V0s", pid);
e3ae862b 257 name = "h_" + type[pid] + "_pt";
3a72645a 258 Float_t pT = esdV0->Pt();
259 fColl->Fill(name, pT);
260 Float_t mass = MassV0(esdV0, pid);
e3ae862b 261 name = "h_" + type[pid] + "_mass";
3a72645a 262 fColl->Fill(name, mass);
263
264 ProcessDaughters(esdV0);
265 if(fMCEvent){
266 ProcessV0sMC(esdV0);
267 ProcessDaughtersMC(esdV0);
e3ae862b 268 ProcessBackground(esdV0);
3a72645a 269 }
270
271 }
272
273}
274//__________________________________________________________
275void AliAnalysisTaskCheckV0tenderII::ProcessDaughters(AliESDv0 * const v0){
276 //
277 // produce some check plots for V0 tender tagged single tracks
278 //
279
e3ae862b 280 const TString type[3] = {"Electron", "Pion", "Proton"};
281 TString name;
3a72645a 282
283 if(!v0) return;
284
285 // daughter tracks
286 const Int_t nTracks = fInputEvent->GetNumberOfTracks();
287 AliESDtrack *d[2];
288 Int_t iN, iP;
289 iN = iP = -1;
290 iP = v0->GetPindex();
291 iN = v0->GetNindex();
292 if(iN < 0 || iP < 0) return;
293 if(iN >= nTracks || iP >= nTracks) return;
294 d[0] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iP));
295 d[1] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iN));
296 if(!d[0] || !d[1]) return;
297
298 for(Int_t i=0; i<2; ++i){
299 Int_t pid = (i == 0) ? PDGtoPID(fpdgP) : PDGtoPID(fpdgN);
300 if(pid < 0) continue;
301 fColl->Fill("h_NumberOfDaughters", pid*1.0);
302 Float_t pT = d[i]->Pt();
e3ae862b 303 name = "h_" + type[pid] + "_pt";
3a72645a 304 fColl->Fill(name, pT);
3a72645a 305 }
306}
307//__________________________________________________________
308void AliAnalysisTaskCheckV0tenderII::ProcessV0sMC(AliESDv0 * const v0){
309 //
310 // check all V0tender selected V0 on their true identity
311 //
312
e3ae862b 313 const TString type[4] = {"Gamma", "K0", "Lambda", "ALambda"};
314 TString name;
3a72645a 315
316 Int_t pid = PDGtoPIDv0(fpdgV0);
317 if(pid < 0) return;
318
319 // true if fpdgV0 agrees with MC pdg of the mother
320 Bool_t id = kFALSE;
321 const Int_t nTracks = fInputEvent->GetNumberOfTracks();
322 // both ESD daughtr tracks
323 Int_t iN, iP;
324 iN = iP = -1;
325 iP = v0->GetPindex();
326 iN = v0->GetNindex();
327 if(iN < 0 || iP < 0) return;
328 if(iN >= nTracks || iP >= nTracks) return;
329 AliESDtrack *dP, *dN;
330 dP = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iP));
331 dN = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iN));
332 if(!dN || !dP) return;
333
334 // MC labels of the daughter tracks
335 Int_t lN, lP;
336 lN = dN->GetLabel();
337 lP = dP->GetLabel();
338 if(lN < 0 || lP < 0) return;
339
340 // MC daughter particles
341 AliMCParticle *mcP, *mcN;
342 mcP = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lP));
343 mcN = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lN));
344 if(!mcP || !mcN) return;
345
346 // labels of the mother particles
347 Int_t lPm, lNm;
348 lPm = mcP->GetMother();
349 lNm = mcN->GetMother();
350 AliMCParticle *m = 0x0;
351 // particles without mother particle are most probably
352 // primary particles
353 if(lPm >= 0){
354 m = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lPm));
355 if(!m) return;
356 }
357 if(m){
358 Int_t pdg = m->PdgCode();
359 //if(lPm == lNm) printf(" -D: pdg: %i, fpdgV0: %i \n", pdg, fpdgV0);
360 if((lPm == lNm) && (pdg == fpdgV0)) id = kTRUE;
361 }
362
e3ae862b 363 if(id) name = "h_" + type[pid] + "_pt_S";
364 else name = "h_" + type[pid] + "_pt_B";
3a72645a 365 Float_t pT = v0->Pt();
366 fCollMC->Fill(name, pT);
367
e3ae862b 368 if(id) name = "h_" + type[pid] + "_mass_S";
369 else name = "h_" + type[pid] + "_mass_B";
3a72645a 370 Float_t mass = MassV0(v0, pid);
371 fCollMC->Fill(name, mass);
372
373
374}
375//__________________________________________________________
376void AliAnalysisTaskCheckV0tenderII::ProcessDaughtersMC(AliESDv0 * const v0){
377 //
378 // check the identity of the V0tender selected V0 daughters
379 // !!! for positive check only the true identity plays a role here,
380 // not the true V0 mother identity (e.g. selected electron could come
bf892a6a 381 // from primary vertex or pion dalitz deca or true gamma conversion) !!!!
3a72645a 382 //
383
e3ae862b 384 const TString type[3] = {"Electron", "Pion", "Proton"};
385 TString name;
3a72645a 386
387 if(!v0) return;
388
389 // daughter tracks
390 const Int_t nTracks = fInputEvent->GetNumberOfTracks();
391 AliESDtrack *d[2];
392 Int_t iN, iP;
393 iN = iP = -1;
394 iP = v0->GetPindex();
395 iN = v0->GetNindex();
396 if(iN < 0 || iP < 0) return;
397 if(iN >= nTracks || iP >= nTracks) return;
398 d[0] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iP));
399 d[1] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iN));
400 if(!d[0] || !d[1]) return;
401
402 //printf(" *** fpdgV0: %i, fpdgP: %i, fpdgN: %i \n", fpdgV0, fpdgP, fpdgN);
403
404 for(Int_t i=0; i<2; ++i){
405 Bool_t id = kFALSE;
406 Int_t pid = (i == 0) ? PDGtoPID(fpdgP) : PDGtoPID(fpdgN);
407 Int_t pdg = (i == 0) ? fpdgP : fpdgN;
408 if(pid < 0) continue;
409 Float_t pT = d[i]->Pt();
410 Int_t label = d[i]->GetLabel();
411 if(label < 0) continue;
412 AliMCParticle *mcp = 0x0;
413 mcp = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(label));
414 if(!mcp) continue;
415 Int_t pdgMC = mcp->PdgCode();
416 //printf(" * pid: %i, pdg: %i, pdgMC: %i \n", pid, pdg, pdgMC);
417 if(pdgMC == pdg) id = kTRUE;
e3ae862b 418 if(id) name = "h_" + type[pid] + "_pt_S";
419 else name = "h_" + type[pid] + "_pt_B";
3a72645a 420 fCollMC->Fill(name, pT);
421 }
422}
423//__________________________________________________________
e3ae862b 424void AliAnalysisTaskCheckV0tenderII::ProcessBackground(AliESDv0 * const v0){
425 //
426 // look at the compostition and the properties of the
427 // miss-identified V0
428 //
429
430 if(!v0) return;
431 Float_t pt = v0->Pt();
432 const TString type[4] = {"Gamma", "K0", "Lambda", "ALambda"};
433 TString name;
434 // daughter tracks
435 const Int_t nTracks = fInputEvent->GetNumberOfTracks();
436 AliESDtrack *d[2];
437 Int_t index[2] = {-1, -1};
438 index[0] = v0->GetPindex();
439 index[1] = v0->GetNindex();
440 if(index[0] < 0 || index[1] < 0) return;
441 if(index[0] >= nTracks || index[1] >= nTracks) return;
442 d[0] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(index[0]));
443 d[1] = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(index[1]));
444 if(!d[0] || !d[1]) return;
445
446 // daughter MC particles
447 Int_t label[2] = {d[0]->GetLabel(), d[1]->GetLabel()};
448 if(label[0] < 0 || label[1] < 0) return;
449 AliMCParticle *dmc[2] = {0x0, 0x0};
450 dmc[0] = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(label[0]));
451 dmc[1] = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(label[1]));
452 if(!dmc[0] || !dmc[1]) return;
453
454 Bool_t primary[2] = {dmc[0]->Particle()->IsPrimary(), dmc[1]->Particle()->IsPrimary()};
455
456 // mother MC particles
457 Int_t labelM[2] = {dmc[0]->GetMother(), dmc[1]->GetMother()};
458 AliMCParticle *mmc[2] = {0x0, 0x0};
459 if(!primary[0]) mmc[0] = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(labelM[0]));
460 if(!primary[1]) mmc[1] = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(labelM[1]));
461
462 // if particle is not primary it MUST have a mother particle
463 if(!primary[0] && !mmc[0]) return;
464 if(!primary[1] && !mmc[1]) return;
465
466 Int_t pdgM[2] = {0, 0};
467 if(mmc[0]) pdgM[0] = mmc[0]->PdgCode();
468 if(mmc[1]) pdgM[1] = mmc[1]->PdgCode();
469
470 // select only miss-identified V0s
471 if(labelM[0] == labelM[1] && pdgM[0] == fpdgV0) return;
472
473 // indices for V0 background histograms
474 // 0 - random match
475 // 1 - gamma
476 // 2 - K0
477 // 3 - lambda
478 // 4 - a-lambda
479 // 5 - other dacay or interaction
480
481 Int_t ixM = -1;
482
483 Int_t typeV0 = (labelM[0] != labelM[1]) ? 0 : 1;
484 if(0 == typeV0) ixM = 0;
485 else{
486 ixM = PDGtoPIDv0(pdgM[0]) + 1;
487 if(ixM < 0) ixM = 5;
488 }
0817c20b 489 Int_t ix = PDGtoPIDv0(fpdgV0);
490 if(0 <= ix ){
491 name = "h_" + type[PDGtoPIDv0(fpdgV0)] + "_Bg";
492 fCollMC->Fill(name, ixM, pt);
493 }
e3ae862b 494
495 // now look at the daughter tracks
496
497
498}
499//__________________________________________________________
3a72645a 500Float_t AliAnalysisTaskCheckV0tenderII::MassV0(AliESDv0 * const v0, Int_t id){
501 //
502 // Get the V0 effective mass
503 //
504
505 Float_t mass = -0.1;
506 Bool_t sign = CheckSigns(v0);
507 if(0 == id){
508 mass = v0->GetEffMass(0, 0);
509 }
510 else if(1 == id){
511 mass = v0->GetEffMass(2, 2);
512 }
513 else if(2 == id){
514 mass = (sign) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);
515 }
516 else if(3 == id){
517 mass = (sign) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
518 }
519 else{
520 AliWarning(Form("Unrecognized V0 id: %i", id));
521 }
522
523 return mass;
524
525}
526//__________________________________________________________
527Bool_t AliAnalysisTaskCheckV0tenderII::CheckSigns(AliESDv0 * const v0){
528 //
529 // check wheter the sign was correctly applied to
530 // V0 daughter tracks
531 // This function should become obsolete once the V0 finder will be updated (fixed)
532 //
533
534 Bool_t correct = kFALSE;
535
536 Int_t pIndex = 0, nIndex = 0;
537 pIndex = v0->GetPindex();
538 nIndex = v0->GetNindex();
539
540 AliESDtrack* d[2];
541 d[0] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(pIndex));
542 d[1] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(nIndex));
543
544 Int_t sign[2];
545 sign[0] = (int)d[0]->GetSign();
546 sign[1] = (int)d[1]->GetSign();
547
548 if(-1 == sign[0] && 1 == sign[1]){
549 correct = kFALSE;
550 //v0->SetIndex(0, pIndex); // set the index of the negative v0 track
551 //v0->SetIndex(1, nIndex); // set the index of the positive v0 track
552 }
553 else{
554 correct = kTRUE;
555 }
556
557 //pIndex = v0->GetPindex();
558 //nIndex = v0->GetNindex();
559 //printf("-D2: P: %i, N: %i\n", pIndex, nIndex);
560
561 return correct;
562}
563//__________________________________________________________
6555e2ad 564Int_t AliAnalysisTaskCheckV0tenderII::PDGtoPIDv0(Int_t pdgV0) const {
3a72645a 565 //
566 // convert thereconstructed V0 pdg to local pid
567 //
568
569 switch(pdgV0){
570 case 22: return 0; break;
571 case 310: return 1; break;
572 case 3122: return 2; break;
573 case -3122: return 3; break;
574 }
575
576 return -1;
577
578}
579//__________________________________________________________
6555e2ad 580Int_t AliAnalysisTaskCheckV0tenderII::PDGtoPID(Int_t pdg) const {
3a72645a 581 //
582 // convert daughter pdg code to local pid
583 //
584 switch(TMath::Abs(pdg)){
585 case 11: return 0; break;
586 case 211: return 1; break;
587 case 2212: return 2; break;
588 }
589 return -1;
590
591}
e3ae862b 592//__________________________________________________________
593void AliAnalysisTaskCheckV0tenderII::ResetPDGcodes(){
594 //
595 // reset the PDG codes of the V0 and daughter
596 // particles to default values
597 //
598
599 fpdgV0 = -1;
600 fpdgP = -1;
601 fpdgN = -1;
602
603}