]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/KINK/AliAnalysisKinkESDMC.cxx
Files and updates for PWG2kink.par as required by Mihaela
[u/mrichter/AliRoot.git] / PWG2 / KINK / AliAnalysisKinkESDMC.cxx
CommitLineData
10eaad41 1/**************************************************************************
2 * Authors: Martha Spyropoulou-Stassinaki and the members
3 * of the Greek group at Physics Department of Athens University
4 * Paraskevi Ganoti, Anastasia Belogianni and Filimon Roukoutakis
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 **************************************************************************/
15
16//-----------------------------------------------------------------
17// AliAnalysisKinkESDMC class
18// Example of an analysis task for kink topology study
19// Kaons from kink topology are 'identified' in this code
20//-----------------------------------------------------------------
21
22#include "TChain.h"
23#include "TTree.h"
24#include "TH1F.h"
25#include "TH2F.h"
26#include "TH1D.h"
27#include "TH2D.h"
28#include "TParticle.h"
29#include <TVector3.h>
30#include "TF1.h"
31
32#include "AliAnalysisTaskSE.h"
33#include "AliAnalysisManager.h"
34
35#include "AliESDEvent.h"
36#include "AliESDInputHandler.h"
37
38#include "AliMCEventHandler.h"
39#include "AliMCEvent.h"
40
41#include "AliAnalysisKinkESDMC.h"
42
43#include "AliStack.h"
44#include "AliESDpid.h"
45#include "AliESDkink.h"
46
47ClassImp(AliAnalysisKinkESDMC)
48
49//________________________________________________________________________
50AliAnalysisKinkESDMC::AliAnalysisKinkESDMC()
51 : AliAnalysisTaskSE(), fESD(0),fHistPtESD(),fHistPt(0),fHistQtAll(0),fHistQt1(0),fHistQt2(0)
52 , fHistPtKaon(0),fHistPtKPDG(0),fHistEta(0),fHistEtaK(0),fptKMC(0),fMultiplMC(0),fESDMult(0),fgenpt(0),frad(0),
53 fKinkKaon(0), fKinkKaonBg(0), fM1kaon(0), fgenPtEtR(0),fPtKink(0), fptKink(0),
54 fcodeH(0), fdcodeH(0), fAngMomK(0),fAngMomPi(0),
55 fRpr(0),fZpr(0),
56 fAngMomKC(0),
57 fZvXv(0),fZvYv(0), fXvYv(0), fPtPrKink(0), f1(0), f2(0),
58 fListOfHistos(0)
59
60{
61 // Constructor
62
63}
64
65
66//________________________________________________________________________
67AliAnalysisKinkESDMC::AliAnalysisKinkESDMC(const char *name)
68 : AliAnalysisTaskSE(name), fESD(0),fHistPtESD(0),fHistPt(0),fHistQtAll(0),fHistQt1(0),fHistQt2(0)
69 , fHistPtKaon(0),fHistPtKPDG(0),fHistEta(0),fHistEtaK(0),fptKMC(0),fMultiplMC(0),fESDMult(0),fgenpt(0),frad(0),
70 fKinkKaon(0), fKinkKaonBg(0), fM1kaon(0), fgenPtEtR(0),fPtKink(0), fptKink(0),
71 fcodeH(0), fdcodeH(0), fAngMomK(0),fAngMomPi(0),
72 fRpr(0),fZpr(0),
73 fAngMomKC(0),
74 fZvXv(0),fZvYv(0), fXvYv(0), fPtPrKink(0), f1(0), f2(0),
75 fListOfHistos(0)
76 // , faod(0), fSignalTree(0)
77
78{
79 // Constructor
80
81 // Define input and output slots here
82 // Input slot #0 works with a TChain
83 DefineInput(0, TChain::Class());
84 DefineOutput(1, TList::Class());
85}
86
87//________________________________________________________________________
88void AliAnalysisKinkESDMC::ConnectInputData(Option_t *)
89{
90 // Connect ESD or AOD here
91 // Called once
92
93 TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
94 if (!tree) {
95 Printf("ERROR: Could not read chain from input slot 0");
96 } else {
97 // Disable all branches, we want to process only MC
98 tree->SetBranchStatus("*", kTRUE);
99 // tree->SetBranchStatus("fTracks.*", kTRUE);
100 tree->SetBranchStatus("*Calo*", kFALSE);
101 tree->SetBranchStatus("FMD", kFALSE);
102
103 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
104
105 if (!esdH) {
106 Printf("ERROR: Could not get ESDInputHandler");
107 } else
108 fESD = esdH->GetEvent();
109 }
110}
111
112//________________________________________________________________________
113void AliAnalysisKinkESDMC::CreateOutputObjects()
114{
115 // Create histograms
116 // Called once
117
118 f1=new TF1("f1","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",1.1,10.0);
119 f1->SetParameter(0,0.493677);
120 f1->SetParameter(1,0.9127037);
121 f1->SetParameter(2,TMath::Pi());
122
123
124 f2=new TF1("f2","((atan([0]*[1]*(1.0/(sqrt((x^2)*(1.0-([1]^2))-([0]^2)*([1]^2))))))*180.)/[2]",0.1,10.0);
125 f2->SetParameter(0,0.13957018);
126 f2->SetParameter(1,0.2731374);
127 f2->SetParameter(2,TMath::Pi());
128 //Open file 1= CAF
129 //OpenFile(1);
130
131 fHistPtESD = new TH1F("fHistPtESD", "P_{T} distribution", 60, 0.0, 6.0);
132 fHistPtESD->GetXaxis()->SetTitle("P_{T} (GeV/c)");
133 fHistPtESD->GetYaxis()->SetTitle("dN/dP_{T} (c/GeV)");
134 fHistPtESD->SetMarkerStyle(kFullCircle);
135 fHistPt = new TH1F("fHistPt", "P_{T} distribution", 60, 0.0, 6.0);
136 fHistQtAll = new TH1F("fHistQtAll", "Q_{T} distr All Kinks ",100, 0.0,.250);
137 fHistQt1= new TH1F("fHistQt1", "Q_{T} distribution",100, 0.0,.300);
138 fHistQt2= new TH1F("fHistQt2", "Q_{T} distribution",100, 0.0,.300);
139 fHistPtKaon = new TH1F("fHistPtKaon", "P_{T}Kaon distribution", 60, 0.0, 6.0);
140 fHistPtKPDG = new TH1F("fHistPtKPDG", "P_{T}Kaon distribution", 60, 0.0, 6.0);
141 fHistEta= new TH1F("fHistEta", "Eta distribution", 26,-1.3, 1.3);
142 fHistEtaK= new TH1F("fHistEtaK", "EtaK distribution", 26,-1.3, 1.3);
143 fptKMC= new TH1F("fptKMC", "P_{T}Kaon generated", 60, 0.0, 6.0);
144 fMultiplMC= new TH1F("fMultiplMC", "charge multiplicity MC", 60, 0.5,120.5);
145 fESDMult= new TH1F("fESDMult", "charge multipliESD", 60, 0.5,120.5);
146 fgenpt= new TH1F("fgenpt", "genpt K distribution", 60, 0.0, 6.0);
147 frad= new TH1F("frad", "radius K generated",100,0.,400.);
148 fKinkKaon= new TH1F("fKinkKaon", "P_{T}Kaon kinks identi", 60, 0.0, 6.0);
149 fKinkKaonBg= new TH1F("fKinkKaonBg", "P_{T}Kaon kinks backgr", 60, 0.0, 6.0);
150 fM1kaon= new TH1F("fM1kaon","Invar m(kaon) from kink->mu+netrino decay",80,0.0, 0.8);
151 fgenPtEtR= new TH1F("fgenPtEtR", "P_{T}Kaon distribution", 60, 0.0, 6.0);
152 fPtKink= new TH1F("fPtKink", "P_{T}Kaon Kink bution", 60, 0.0, 6.0);
153 fptKink= new TH1F("fptKink", "P_{T}Kaon Kink bution", 60, 0.0, 6.0);
154 fcodeH = new TH2F("fcodeH", "code vrs dcode dist. kinks,K",100,0.,2500.,100,0.,2500.);
155 fdcodeH = new TH2F("fdcodeH", "code vrs dcode dist. kinks,K",100,0.,2500.,100,0.,2500.);
156 fAngMomK= new TH2F("fAngMomK","Decay angle vrs Mother Mom,K",100,0.0,5.0,80,0.,80.);
157 fAngMomPi= new TH2F("fAngMomPi","Decay angle vrs Mother Mom,Pi",100,0.0,5.0,80,0.,80.);
158
159 fRpr = new TH1D("fRpr", "rad distribution PID pr",50,0.0, 2.5);
160 fZpr = new TH1D("fZpr", "z distribution PID pr ",60,-15.,15.);
161 fAngMomKC= new TH2F("fAngMomKC","Decay angle vrs Mother Mom,K",100,0.0,5.0,80,0.,80.);
162 fZvXv= new TH2F("fZvXv","Xv-Zv main vtx",60,-0.5,0.5,60, -15., 15.0);
163 fZvYv= new TH2F("fZvYv","Yv-Zv main vtx",60,-0.5,0.5, 60, -15., 15.);
164 fXvYv= new TH2F("fXvYv","Xv-Yv main vtx", 60,-1.5,1.5, 60, -1.5, 1.5);
165 fPtPrKink=new TH1F("fPtPrKink","pt of ESD kaonKink tracks", 60, 0.0, 6.0);
166
167 fListOfHistos=new TList();
168
169 fListOfHistos->Add(fHistPtESD);
170 fListOfHistos->Add(fHistPt);
171 fListOfHistos->Add(fHistQtAll);
172 fListOfHistos->Add(fHistQt1);
173 fListOfHistos->Add(fHistQt2);
174 fListOfHistos->Add(fHistPtKaon);
175 fListOfHistos->Add(fHistPtKPDG);
176 fListOfHistos->Add(fHistEta);
177 fListOfHistos->Add(fHistEtaK);
178 fListOfHistos->Add(fptKMC);
179 fListOfHistos->Add(fMultiplMC);
180 fListOfHistos->Add(fESDMult);
181 fListOfHistos->Add(fgenpt);
182 fListOfHistos->Add(frad);
183 fListOfHistos->Add(fKinkKaon);
184 fListOfHistos->Add(fKinkKaonBg);
185 fListOfHistos->Add(fM1kaon);
186 fListOfHistos->Add(fgenPtEtR);
187 fListOfHistos->Add(fPtKink);
188 fListOfHistos->Add(fptKink);
189 fListOfHistos->Add(fcodeH);
190 fListOfHistos->Add(fdcodeH);
191 fListOfHistos->Add(fdcodeH);
192 fListOfHistos->Add(fAngMomK);
193 fListOfHistos->Add(fAngMomPi);
194 fListOfHistos->Add(fRpr);
195 fListOfHistos->Add(fZpr);
196 fListOfHistos->Add(fAngMomKC);
197 fListOfHistos->Add(fZvXv);
198 fListOfHistos->Add(fZvYv);
199 fListOfHistos->Add(fXvYv);
200 fListOfHistos->Add(fPtPrKink);
201}
202
203//________________________________________________________________________
204void AliAnalysisKinkESDMC::Exec(Option_t *)
205{
206 // Main loop
207 // Called for each event
208
209 // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
210 // This handler can return the current MC event
211
212 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
213 if (!eventHandler) {
214 Printf("ERROR: Could not retrieve MC event handler");
215 return;
216 }
217
218 AliMCEvent* mcEvent = eventHandler->MCEvent();
219 if (!mcEvent) {
220 Printf("ERROR: Could not retrieve MC event");
221 return;
222 }
223
224 Printf("MC particles: %d", mcEvent->GetNumberOfTracks());
225
226 if (!fESD) {
227 Printf("ERROR: fESD not available");
228 return;
229 }
230 AliStack* stack=mcEvent->Stack();
231
232//primary tracks in MC
233 Int_t nPrim = stack->GetNprimary();
234//
235 const AliESDVertex *vertex=GetEventVertex(fESD);
236 if(!vertex) return;
237//
238 Double_t vpos[3];
239 vertex->GetXYZ(vpos);
240 fZpr->Fill(vpos[2]);
241 fZvXv->Fill( vpos[1], vpos[2]);
242 fZvYv->Fill( vpos[0], vpos[2]);
243
244 Double_t vtrack[3], ptrack[3];
245
246
247
248 Int_t nGoodTracks = fESD->GetNumberOfTracks();
249 fESDMult->Fill(nGoodTracks);
250//
251// track loop
252 for (Int_t iTracks = 0; iTracks < fESD->GetNumberOfTracks(); iTracks++) {
253 AliESDtrack* track = fESD->GetTrack(iTracks);
254 if (!track) {
255 Printf("ERROR: Could not receive track %d", iTracks);
256 continue;
257 }
258
259
260 fHistPt->Fill(track->Pt());
261
262 // pt cut at 300 MeV
263 if ( (track->Pt())<.300)continue;
264
265 UInt_t status=track->GetStatus();
266
267 // if((status&AliESDtrack::kITSrefit)==0) continue;
268 if((status&AliESDtrack::kTPCrefit)==0) continue;
269 if((track->GetTPCchi2()/track->GetTPCclusters(0))>3.5) continue;
270
271 Double_t extCovPos[15];
272 track->GetExternalCovariance(extCovPos);
273 if(extCovPos[0]>2) continue;
274 if(extCovPos[2]>2) continue;
275 if(extCovPos[5]>0.5) continue;
276 if(extCovPos[9]>0.5) continue;
277 if(extCovPos[14]>2) continue;
278
279
280 track->GetXYZ(vtrack);
281 fXvYv->Fill(vtrack[0],vtrack[1]);
282
283// track momentum
284 track->GetPxPyPz(ptrack);
285
286 TVector3 trackMom(ptrack[0],ptrack[1],ptrack[2]);
287
288 Double_t trackEta=trackMom.Eta();
289
290
291
292 Float_t nSigmaToVertex = GetSigmaToVertex(track);
293
294 Float_t bpos[2];
295 Float_t bCovpos[3];
296 track->GetImpactParameters(bpos,bCovpos);
297
298 if (bCovpos[0]<=0 || bCovpos[2]<=0) {
299 Printf("Estimated b resolution lower or equal zero!");
300 bCovpos[0]=0; bCovpos[2]=0;
301 }
302
303 Float_t dcaToVertexXYpos = bpos[0];
304 Float_t dcaToVertexZpos = bpos[1];
305
306 fRpr->Fill(dcaToVertexZpos);
307
308 if(nSigmaToVertex>=4) continue;
309 if((dcaToVertexXYpos>3.0)||(dcaToVertexZpos>3.0)) continue;
310
311
312 Int_t label = track->GetLabel();
313 label = TMath::Abs(label);
314 TParticle * part = stack->Particle(label);
315 if (!part) continue;
316
317 // cut on eta
318 if( (TMath::Abs(trackEta )) > 1.1 ) continue;
319 fHistPtESD->Fill(track->Pt());
320
321 // Add Kink analysis
322
323 Int_t indexKinkPos=track->GetKinkIndex(0);
324// loop on kinks
325 if(indexKinkPos<0){
326 fPtKink->Fill(track->Pt()); /// pt from track
327
328 // select kink class
329
330 AliESDkink *kink=fESD->GetKink(TMath::Abs(indexKinkPos)-1);
331//
332
333 Int_t eSDfLabel1=kink->GetLabel(0);
334 TParticle *particle1= stack->Particle(TMath::Abs(eSDfLabel1));
335 Int_t code1= particle1->GetPdgCode();
336
337 Int_t eSDfLabeld=kink->GetLabel(1);
338 TParticle *particled= stack->Particle(TMath::Abs(eSDfLabeld));
339 Int_t dcode1= particled->GetPdgCode();
340
341 const TVector3 motherMfromKink(kink->GetMotherP());
342 const TVector3 daughterMKink(kink->GetDaughterP());
343 Float_t qT=kink->GetQt();
344
345 fHistEta->Fill(trackEta) ; // Eta distr
346 fHistQtAll->Fill(qT) ; // Qt distr
347
348 fptKink->Fill(motherMfromKink.Pt()); /// pt from kink
349
350// fake kinks, small Qt and small kink angle
351 if ( qT<0.012) fcodeH->Fill( TMath::Abs(code1), TMath::Abs(dcode1) );
352 if(qT<0.012) continue; // remove fake kinks
353
354 Float_t kinkAngle=TMath::RadToDeg()*kink->GetAngle(2);
355
356//remove the double taracks
357
358 if( (kinkAngle<1.) ) continue;
359//
360 if( (qT>0.05)&& ( qT<0.25) ) fHistQt2->Fill(qT); // candidate kaon kinks
361
362 if (TMath::Abs(code1)==321) fHistPtKPDG->Fill(track->Pt()); // ALL KAONS (pdg) inside kink sample
363
364// maximum decay angle at a given mother momentum
365 Double_t maxDecAngKmu=f1->Eval(motherMfromKink.Mag(),0.,0.,0.);
366 Double_t maxDecAngpimu=f2->Eval(motherMfromKink.Mag(),0.,0.,0.);
367// two dimensional plot
368 if(TMath::Abs(code1)==321) fAngMomK->Fill(motherMfromKink.Mag(), kinkAngle);
369 if(TMath::Abs(code1)==211) fAngMomPi->Fill(motherMfromKink.Mag(), kinkAngle);
370//
371// invariant mass of mother track decaying to mu
372 Float_t energyDaughterMu=TMath::Sqrt(daughterMKink.Mag()*daughterMKink.Mag()+0.105658*0.105658);
373 Float_t p1XM= motherMfromKink.Px();
374 Float_t p1YM= motherMfromKink.Py();
375 Float_t p1ZM= motherMfromKink.Pz();
376 Float_t p2XM= daughterMKink.Px();
377 Float_t p2YM= daughterMKink.Py();
378 Float_t p2ZM= daughterMKink.Pz();
379 Float_t p3Daughter=TMath::Sqrt(((p1XM-p2XM)*(p1XM-p2XM))+((p1YM-p2YM)*(p1YM-p2YM))+((p1ZM-p2ZM)*(p1ZM-p2ZM)));
380 Double_t invariantMassKmu= TMath::Sqrt((energyDaughterMu+p3Daughter)*(energyDaughterMu+p3Daughter)-motherMfromKink.Mag()*motherMfromKink.Mag());
381
382 fM1kaon->Fill(invariantMassKmu);
383
384// kaon selection from kinks
385if((kinkAngle>maxDecAngpimu)&&(qT>0.05)&&(qT<0.25)&&((kink->GetR()>110.)&&(kink->GetR()<230.))&&(TMath::Abs(trackEta)<1.1)&&(invariantMassKmu<0.6)) {
386 fHistEtaK->Fill(trackEta);
387
388 if(TMath::Abs(code1)==321) fAngMomKC->Fill(motherMfromKink.Mag(), kinkAngle); // real kaons
389
390 fHistQt1 ->Fill(qT) ; // Qt distr
391 if ( (kinkAngle>maxDecAngKmu) && ( motherMfromKink.Mag()> 1.1 ) ) continue; // maximum angle selection
392
393 fHistPtKaon->Fill(track->Pt()); //all PID kink-kaon
394
395// background inside the identified kaons, e.g KK , Kp
396 if((TMath::Abs(code1)==321)&&(( TMath::Abs(dcode1)) >=(TMath::Abs(code1)) )) fdcodeH->Fill( TMath::Abs(code1), TMath::Abs(dcode1));
397
398 if((TMath::Abs(code1)==321) ) {
399
400 if ( label< nPrim ) fPtPrKink->Fill(track->Pt());
401 fKinkKaon->Fill(track->Pt()); }
402 else {
403 fKinkKaonBg->Fill(track->Pt());
404 } // primary and secondary
405
406 } // kink selection
407
408
409 } //End Kink Information
410
411
412 } //track loop
413
414 // loop over mc particles
415
416 // variable to count tracks
417 Int_t nMCTracks = 0;
418
419for (Int_t iMc = 0; iMc < nPrim; ++iMc)
420 {
421 // AliDebug(AliLog::kDebug+1, Form("MC Loop: Processing particle %d.", iMc));
422
423 TParticle* particle = stack->Particle(iMc);
424
425 if (!particle)
426 {
427 // AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (mc loop).", iMc));
428 continue;
429 }
430
431 //
432 Float_t eta = particle->Eta();
433
434 if (TMath::Abs(eta) > 1.1)
435 continue;
436
437 Double_t ptK = particle->Pt();
438
439 if( ptK <0.300) continue;
440
441 Float_t code = particle->GetPdgCode();
442 if (code==321 || code==-321){
443
444 fptKMC->Fill(ptK);
445
446
447 Int_t firstD=particle->GetFirstDaughter();
448 Int_t lastD=particle->GetLastDaughter();
449//loop on secondaries
450 for (Int_t k=firstD;k<lastD;k++) {
451 TParticle* daughter1=stack->Particle(k+1);
452 Float_t dcode = daughter1->GetPdgCode();
453
454
455
456 frad->Fill(daughter1->R());
457
458 if (((daughter1->R())>110)&&((daughter1->R())<230) ){
459 if (( ( code==321 )&& ( dcode ==-13 ))|| (( code == -321 )&& ( dcode == 13))) fgenPtEtR->Fill( ptK );//to muon
460 if (( ( code==321 )&& ( dcode ==211 ))|| (( code == -321 )&& ( dcode ==-211))) fgenPtEtR->Fill( ptK );//to pion
461 // if (( ( code==321 )&& ( dcode ==-11 ))|| (( code == -321 )&& ( dcode ==11))) fgenPtEtR->Fill( ptK );// to electr
462
463 }
464// next only to mu decay
465 if (((code==321)&&(dcode==-13))||((code==-321)&&(dcode==13))){
466
467
468 if (((daughter1->R())>110)&&((daughter1->R())<230)){
469 fgenpt->Fill(ptK);
470 }
471 }
472 }
473
474 }
475
476
477
478 nMCTracks++;
479 }// end of mc particle
480
481// printf("hello mc");
482 fMultiplMC->Fill(nMCTracks);
483
484 PostData(1, fListOfHistos);
485
486}
487
488//________________________________________________________________________
489void AliAnalysisKinkESDMC::Terminate(Option_t *)
490{
491 // Draw result to the screen
492 // Called once at the end of the query
493
494}
495//____________________________________________________________________//
496
497Float_t AliAnalysisKinkESDMC::GetSigmaToVertex(AliESDtrack* esdTrack) const {
498 // Calculates the number of sigma to the vertex.
499
500 Float_t b[2];
501 Float_t bRes[2];
502 Float_t bCov[3];
503
504 esdTrack->GetImpactParameters(b,bCov);
505
506 if (bCov[0]<=0 || bCov[2]<=0) {
507 //AliDebug(1, "Estimated b resolution lower or equal zero!");
508 bCov[0]=0; bCov[2]=0;
509 }
510 bRes[0] = TMath::Sqrt(bCov[0]);
511 bRes[1] = TMath::Sqrt(bCov[2]);
512
513 if (bRes[0] == 0 || bRes[1] ==0) return -1;
514
515 Float_t d = TMath::Sqrt(TMath::Power(b[0]/bRes[0],2) + TMath::Power(b[1]/bRes[1],2));
516
517 if (TMath::Exp(-d * d / 2) < 1e-10) return 1000;
518
519 d = TMath::ErfInverse(1 - TMath::Exp(-d * d / 2)) * TMath::Sqrt(2);
520
521 return d;
522 }
523//____________________________________________________________________//
524
525const AliESDVertex* AliAnalysisKinkESDMC::GetEventVertex(const AliESDEvent* esd) const
526 // Get the vertex from the ESD and returns it if the vertex is valid
527
528{
529 // Get the vertex
530
531 const AliESDVertex* vertex = esd->GetPrimaryVertex();
532
533 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
534 else
535 {
536 vertex = esd->GetPrimaryVertexSPD();
537 if((vertex->GetStatus()==kTRUE)&&(vertex->GetNContributors()>2)) return vertex;
538 else
539 return 0;
540 }
541}