]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/global/AliAnalysisTaskVertexESD.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Added method to find the local maximum (Giacomo)
[u/mrichter/AliRoot.git] / PWG1 / global / AliAnalysisTaskVertexESD.cxx
CommitLineData
388ca814 1/**************************************************************************
2 * Copyright(c) 1998-2009, 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 **************************************************************************/
15
16//*************************************************************************
17// Class AliAnalysisTaskVertexESD
18// AliAnalysisTask to extract from ESD the information for the analysis
19// of the primary vertex reconstruction efficiency and resolution
20// (for MC events) and distributions (for real data). Three vertices:
21// - SPD tracklets
22// - ITS+TPC tracks
23// - TPC-only tracks
24//
25// Author: A.Dainese, andrea.dainese@pd.infn.it
26//*************************************************************************
27
28#include <TChain.h>
29#include <TTree.h>
30#include <TNtuple.h>
31#include <TBranch.h>
32#include <TClonesArray.h>
33#include <TObjArray.h>
34#include <TH1F.h>
35#include <TH2F.h>
36#include <TCanvas.h>
a874bb01 37#include <TGraphAsymmErrors.h>
38#include <TFile.h>
388ca814 39
40#include "AliAnalysisTask.h"
41#include "AliAnalysisManager.h"
42
43#include "AliMultiplicity.h"
44#include "AliVertexerTracks.h"
45#include "AliESDtrack.h"
46#include "AliExternalTrackParam.h"
47#include "AliESDVertex.h"
48#include "AliESDEvent.h"
388ca814 49#include "AliESDInputHandler.h"
388ca814 50#include "AliTrackReference.h"
38d12819 51//#include "AliTriggerAnalysis.h"
388ca814 52
53#include "AliMCEventHandler.h"
54#include "AliMCEvent.h"
55#include "AliStack.h"
56#include "AliLog.h"
57
58#include "AliGenEventHeader.h"
59#include "AliAnalysisTaskVertexESD.h"
60
61
62ClassImp(AliAnalysisTaskVertexESD)
63
64//________________________________________________________________________
65AliAnalysisTaskVertexESD::AliAnalysisTaskVertexESD(const char *name) :
74917c1f 66AliAnalysisTaskSE(name),
a874bb01 67 fCheckEventType(kTRUE),
68 fReadMC(kFALSE),
d0bbaf02 69 fRecoVtxTPC(kTRUE),
a874bb01 70 fRecoVtxITSTPC(kTRUE),
71 fRecoVtxITSTPCHalfEvent(kFALSE),
72 fOnlyITSTPCTracks(kFALSE),
73 fOnlyITSSATracks(kFALSE),
74 fFillNtuple(kFALSE),
75 fFillTreeBeamSpot(kFALSE),
76 fESD(0),
77 fOutput(0),
78 fNtupleVertexESD(0),
79 fhSPDVertexX(0),
80 fhSPDVertexY(0),
81 fhSPDVertexZ(0),
82 fhTRKVertexX(0),
83 fhTRKVertexY(0),
84 fhTRKVertexZ(0),
85 fhTPCVertexX(0),
86 fhTPCVertexY(0),
87 fhTPCVertexZ(0),
88 fhTrackRefs(0),
89 fTreeBeamSpot(0),
90 fhTriggeredTrklets(0),
91 fhSPD3DTrklets(0),
92 fhSPDZTrklets(0),
93 fhTRKTrklets(0),
94 fhTRKcTrklets(0),
95 fhTRKncTrklets(0),
d0bbaf02 96 fhTPCTrklets(0),
97 fhTPCcTrklets(0),
98 fhTPCncTrklets(0),
a874bb01 99 fhSPD3DZreco(0),
100 fhSPDZZreco(0),
101 fhSPDVertexXPile(0),
102 fhSPDVertexYPile(0),
103 fhSPDVertexZPile(0),
8d3b93cb 104 fhSPDVertexDiffZPileContr2(0),
105 fhSPDVertexDiffZPileContr3(0),
106 fhSPDVertexDiffZPileContr4(0),
107 fhSPDVertexDiffZPileContr5(0),
a874bb01 108 fhSPDContributorsPile(0),
109 fhSPDDispContributors(0)
388ca814 110{
111 // Constructor
112
113 // Define input and output slots here
388ca814 114 // Output slot #0 writes into a TList container
41deb8e3 115 DefineOutput(1, TList::Class()); //My private output
388ca814 116}
117//________________________________________________________________________
118AliAnalysisTaskVertexESD::~AliAnalysisTaskVertexESD()
119{
120 // Destructor
121
122 // histograms are in the output list and deleted when the output
123 // list is deleted by the TSelector dtor
124
125 if (fOutput) {
126 delete fOutput;
127 fOutput = 0;
128 }
129}
130
388ca814 131
132//________________________________________________________________________
74917c1f 133void AliAnalysisTaskVertexESD::UserCreateOutputObjects()
388ca814 134{
135 // Create histograms
136 // Called once
137
138 // Several histograms are more conveniently managed in a TList
139 fOutput = new TList;
140 fOutput->SetOwner();
141
19076e2d 142 fNtupleVertexESD = new TNtuple("fNtupleVertexESD","vertices","run:tstamp:bunchcross:triggered:dndygen:xtrue:ytrue:ztrue:xSPD:xerrSPD:ySPD:yerrSPD:zSPD:zerrSPD:ntrksSPD:SPD3D:dphiSPD:xTPC:xerrTPC:yTPC:yerrTPC:zTPC:zerrTPC:ntrksTPC:constrTPC:xTRK:xerrTRK:yTRK:yerrTRK:zTRK:zerrTRK:ntrksTRK:constrTRK:ntrklets:nESDtracks:nITSrefit5or6:nTPCin:nTPCinEta09:SPD0cls:xSPDp:xerrSPDp:ySPDp:yerrSPDp:zSPDp:zerrSPDp:ntrksSPDp:xTPCnc:xerrTPCnc:yTPCnc:yerrTPCnc:zTPCnc:zerrTPCnc:ntrksTPCnc:xTPCc:xerrTPCc:yTPCc:yerrTPCc:zTPCc:zerrTPCc:ntrksTPCc:xTRKnc:xerrTRKnc:yTRKnc:yerrTRKnc:zTRKnc:zerrTRKnc:ntrksTRKnc:xTRKc:xerrTRKc:yTRKc:yerrTRKc:zTRKc:zerrTRKc:ntrksTRKc:deltaxTRKnc:deltayTRKnc:deltazTRKnc:deltaxerrTRKnc:deltayerrTRKnc:deltazerrTRKnc:ntrksEvenTRKnc:ntrksOddTRKnc");
388ca814 143
144 fOutput->Add(fNtupleVertexESD);
145
b588553c 146 fhSPDVertexX = new TH1F("fhSPDVertexX","SPDVertex x; x vertex [cm]; events",200,-1,1);
147 fOutput->Add(fhSPDVertexX);
148 fhSPDVertexY = new TH1F("fhSPDVertexY","SPDVertex y; y vertex [cm]; events",200,-1,1);
149 fOutput->Add(fhSPDVertexY);
150 fhSPDVertexZ = new TH1F("fhSPDVertexZ","SPDVertex z; z vertex [cm]; events",200,-20,20);
151 fOutput->Add(fhSPDVertexZ);
152 fhTRKVertexX = new TH1F("fhTRKVertexX","TRKVertex x; x vertex [cm]; events",200,-1,1);
153 fOutput->Add(fhTRKVertexX);
154 fhTRKVertexY = new TH1F("fhTRKVertexY","TRKVertex y; y vertex [cm]; events",200,-1,1);
155 fOutput->Add(fhTRKVertexY);
156 fhTRKVertexZ = new TH1F("fhTRKVertexZ","TRKVertex z; z vertex [cm]; events",200,-20,20);
157 fOutput->Add(fhTRKVertexZ);
158 fhTPCVertexX = new TH1F("fhTPCVertexX","TPCVertex x; x vertex [cm]; events",200,-3,3);
159 fOutput->Add(fhTPCVertexX);
160 fhTPCVertexY = new TH1F("fhTPCVertexY","TPCVertex y; y vertex [cm]; events",200,-3,3);
161 fOutput->Add(fhTPCVertexY);
162 fhTPCVertexZ = new TH1F("fhTPCVertexZ","TPCVertex z; z vertex [cm]; events",200,-20,20);
163 fOutput->Add(fhTPCVertexZ);
a874bb01 164
165
166 fhSPDVertexXPile = new TH1F("fhSPDVertexXPile","SPDVertexPile x; x vertex [cm]; events",200,-20,20);
167 fOutput->Add(fhSPDVertexXPile);
168 fhSPDVertexYPile = new TH1F("fhSPDVertexYPile","SPDVertexPile y; y vertex [cm]; events",200,-20,20);
169 fOutput->Add(fhSPDVertexYPile);
170 fhSPDVertexZPile = new TH1F("fhSPDVertexZPile","SPDVertexPile z; z vertex [cm]; events",200,-40,40);
171 fOutput->Add(fhSPDVertexZPile);
8d3b93cb 172
173 fhSPDVertexDiffZPileContr2 = new TH1F("fhSPDVertexDiffZPileContr2","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
174 fOutput->Add(fhSPDVertexDiffZPileContr2);
175 fhSPDVertexDiffZPileContr3 = new TH1F("fhSPDVertexDiffZPileContr3","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
176 fOutput->Add(fhSPDVertexDiffZPileContr3);
177 fhSPDVertexDiffZPileContr4 = new TH1F("fhSPDVertexDiffZPileContr4","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
178 fOutput->Add(fhSPDVertexDiffZPileContr4);
179 fhSPDVertexDiffZPileContr5 = new TH1F("fhSPDVertexDiffZPileContr5","SPDVertexDiff z; zmain - zdiff [cm]; events",2000,-100,100);
180 fOutput->Add(fhSPDVertexDiffZPileContr5);
b588553c 181
388ca814 182 fhTrackRefs = new TH2F("fhTrackRefs","Track references; x; y",1000,-4,4,1000,-4,4);
183 fOutput->Add(fhTrackRefs);
184
4f0574ad 185 fTreeBeamSpot = new TTree("fTreeBeamSpot", "BeamSpotTree");
186 UShort_t triggered, ntrkletsS, ntrksTRKnc;
187 UInt_t run, bx;
188 Float_t cetTimeLHC,xTRKnc, yTRKnc, zTRKnc;
189 fTreeBeamSpot->Branch("run", &run, "run/i");
190 fTreeBeamSpot->Branch("cetTimeLHC", &cetTimeLHC, "cetTimeLHC/F");
191 fTreeBeamSpot->Branch("bx", &bx, "bx/i");
192 fTreeBeamSpot->Branch("triggered", &triggered, "triggered/s");
193 fTreeBeamSpot->Branch("ntrkletsS", &ntrkletsS, "ntrkletsS/s");
194 fTreeBeamSpot->Branch("xTRKnc", &xTRKnc, "xTRKnc/F");
195 fTreeBeamSpot->Branch("yTRKnc", &yTRKnc, "yTRKnc/F");
196 fTreeBeamSpot->Branch("zTRKnc", &zTRKnc, "zTRKnc/F");
197 fTreeBeamSpot->Branch("ntrksTRKnc", &ntrksTRKnc, "ntrksTRKnc/s");
198 fOutput->Add(fTreeBeamSpot);
aea5af08 199
d0bbaf02 200 Int_t nbinTrklets=29;
201 Float_t lowTrklets[30]={-0.5,0.5,1.5,2.5,3.5,4.5,5.5,7.5,10.5,25.5,50.5,75.5,100.5,150.5,200.5,300.5,400.5,500.5,600.5,800.5,1000.5,1500.5,2000.5,2500.5,3000.5,4000.5,5000.5,6000.5,8000.5,10000.5};
a874bb01 202 fhTriggeredTrklets = new TH1F("fhTriggeredTrklets","trklets dist for triggered ev.; ntrklets; entries",nbinTrklets,lowTrklets);
203 fOutput->Add(fhTriggeredTrklets);
204 fhSPD3DTrklets = new TH1F("fhSPD3DTrklets","trklets dist for SPD3D ev.; ntrklets; entries",nbinTrklets,lowTrklets);
205 fOutput->Add(fhSPD3DTrklets);
206 fhSPDZTrklets = new TH1F("fhSPDZTrklets","trklets dist for SPDZ ev.; ntrklets; entries",nbinTrklets,lowTrklets);
207 fOutput->Add(fhSPDZTrklets);
208 fhTRKTrklets = new TH1F("fhTRKTrklets","trklets dist for TRK ev.; ntrklets; entries",nbinTrklets,lowTrklets);
209 fOutput->Add(fhTRKTrklets);
210 fhTRKcTrklets = new TH1F("fhTRKcTrklets","trklets dist for TRKc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
211 fOutput->Add(fhTRKcTrklets);
212 fhTRKncTrklets = new TH1F("fhTRKncTrklets","trklets dist for TRKnc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
213 fOutput->Add(fhTRKncTrklets);
d0bbaf02 214 fhTPCTrklets = new TH1F("fhTPCTrklets","trklets dist for TPC ev.; ntrklets; entries",nbinTrklets,lowTrklets);
215 fOutput->Add(fhTPCTrklets);
216 fhTPCcTrklets = new TH1F("fhTPCcTrklets","trklets dist for TPCc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
217 fOutput->Add(fhTPCcTrklets);
218 fhTPCncTrklets = new TH1F("fhTPCncTrklets","trklets dist for TPCnc ev.; ntrklets; entries",nbinTrklets,lowTrklets);
219 fOutput->Add(fhTPCncTrklets);
a874bb01 220
221 Int_t nbinZreco = 16;
222 Double_t lowZreco[17]={-15.0,-10.0,-7.0,-5,-4,-3,-2,-1,0,1,2,3,4,5,7,10,15};
223 fhSPD3DZreco = new TH1F("fhSPD3DZreco","Zreco dist for SPD3D ev.; Zreco [cm]; entries",nbinZreco,lowZreco);
224 fOutput->Add(fhSPD3DZreco);
225 fhSPDZZreco = new TH1F("fhSPDZZreco","Zreco dist for SPDZ ev.; Zreco [cm]; entries",nbinZreco,lowZreco);
226 fOutput->Add(fhSPDZZreco);
227
228 fhSPDContributorsPile = new TH1F("fhSPDContributorsPile","ncontributors pile up vertex; ncontributors; entries",200,-0.5,199.5);
229 fOutput->Add(fhSPDContributorsPile);
230
231 fhSPDDispContributors = new TH2F("fhSPDDispContributors","ncontributors main-pile; ncontributors main; ncontributors pile",200,-0.5,199.5,200,-0.5,199.5);
232 fOutput->Add(fhSPDDispContributors);
233
997f86cb 234 PostData(1, fOutput);
30b05a14 235
f6232e6e 236 return;
388ca814 237}
238
239//________________________________________________________________________
74917c1f 240void AliAnalysisTaskVertexESD::UserExec(Option_t *)
388ca814 241{
242 // Main loop
243 // Called for each event
244
74917c1f 245 if (!InputEvent()) {
388ca814 246 Printf("ERROR: fESD not available");
247 return;
248 }
74917c1f 249 AliESDEvent* esdE = (AliESDEvent*) InputEvent();
250
129c7f70 251 // Select PHYSICS events (type=7, for data) and MC events (type=0)
252 // fCheckEventType is kFALSE if fReadMC is kTRUE, hence check is skipped
253 if(fCheckEventType) {
254 if(esdE->GetEventType()!=7 && esdE->GetEventType()!=0) return;
255 }
38d12819 256
257
388ca814 258 TArrayF mcVertex(3);
259 mcVertex[0]=9999.; mcVertex[1]=9999.; mcVertex[2]=9999.;
260 Float_t dNchdy=-999.;
388ca814 261 // *********** MC info ***************
262 if (fReadMC) {
263 AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
264 if (!eventHandler) {
265 Printf("ERROR: Could not retrieve MC event handler");
266 return;
267 }
268
269 AliMCEvent* mcEvent = eventHandler->MCEvent();
270 if (!mcEvent) {
271 Printf("ERROR: Could not retrieve MC event");
272 return;
273 }
274
275 AliStack* stack = mcEvent->Stack();
276 if (!stack) {
277 AliDebug(AliLog::kError, "Stack not available");
278 return;
279 }
280
281 AliHeader* header = mcEvent->Header();
282 if (!header) {
283 AliDebug(AliLog::kError, "Header not available");
284 return;
285 }
286 AliGenEventHeader* genHeader = header->GenEventHeader();
287 genHeader->PrimaryVertex(mcVertex);
288
289 Int_t ngenpart = (Int_t)stack->GetNtrack();
290 //printf("# generated particles = %d\n",ngenpart);
291 dNchdy=0;
292 for(Int_t ip=0; ip<ngenpart; ip++) {
293 TParticle* part = (TParticle*)stack->Particle(ip);
294 // keep only electorns, muons, pions, kaons and protons
295 Int_t apdg = TMath::Abs(part->GetPdgCode());
296 if(apdg!=11 && apdg!=13 && apdg!=211 && apdg!=321 && apdg!=2212) continue;
297 // reject secondaries
298 if(TMath::Sqrt((part->Vx()-mcVertex[0])*(part->Vx()-mcVertex[0])+(part->Vy()-mcVertex[1])*(part->Vy()-mcVertex[1]))>0.0010) continue;
299 // reject incoming protons
300 Double_t energy = part->Energy();
301 if(energy>900.) continue;
302 Double_t pz = part->Pz();
303 Double_t y = 0.5*TMath::Log((energy+pz+1.e-13)/(energy-pz+1.e-13));
304 if(TMath::Abs(y)<1.0) dNchdy += 0.5; // count 1/2 of particles in |y|<1
305 // tracks refs
306 TClonesArray *trefs=0;
307 Int_t ntrefs = mcEvent->GetParticleAndTR(ip,part,trefs);
308 if(ntrefs<0) continue;
309 for(Int_t iref=0; iref<ntrefs; iref++) {
310 AliTrackReference *tref = (AliTrackReference*)trefs->At(iref);
311 if(tref->R()>10.) continue;
312 fhTrackRefs->Fill(tref->X(),tref->Y());
313 }
314 }
315 //printf("# primary particles = %7.1f\n",dNchdy);
316 }
317 // *********** MC info ***************
318
388ca814 319
320 // Trigger
38d12819 321 //ULong64_t triggerMask;
322 //ULong64_t spdFO = (1 << 14);
323 //ULong64_t v0left = (1 << 10);
324 //ULong64_t v0right = (1 << 11);
388ca814 325
74917c1f 326 //triggerMask=esdE->GetTriggerMask();
388ca814 327 // MB1: SPDFO || V0L || V0R
38d12819 328 //Bool_t eventTriggered = (triggerMask & spdFO || ((triggerMask & v0left) || (triggerMask & v0right)));
388ca814 329 //MB2: GFO && V0R
330 //triggerMask & spdFO && ((triggerMask&v0left) || (triggerMask&v0right))
38d12819 331 //Bool_t eventTriggered = (triggerMask & spdFO);
388ca814 332
38d12819 333 //static AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis();
74917c1f 334 Bool_t eventTriggered = 0;//triggerAnalysis->IsTriggerFired(esdE, AliTriggerAnalysis::kSPDGFO /*| AliTriggerAnalysis::kOfflineFlag*/);
388ca814 335
8c7764e0 336 // use response of AliPhysicsSelection
337 eventTriggered = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
338
4f0574ad 339
340 Int_t nInFile = esdE->GetEventNumberInFile();
341
30b05a14 342 const AliMultiplicity *alimult = esdE->GetMultiplicity();
343 Int_t ntrklets=0,spd0cls=0;
344 if(alimult) {
345 ntrklets = alimult->GetNumberOfTracklets();
346
30b05a14 347 for(Int_t l=0;l<alimult->GetNumberOfTracklets();l++){
348 if(alimult->GetDeltaPhi(l)<-9998.) ntrklets--;
349 }
350 spd0cls = alimult->GetNumberOfSingleClusters()+ntrklets;
351 }
352
353
4f0574ad 354 UShort_t triggered, ntrkletsS, ntrksTRKnc;
355 UInt_t run, bx;
356 Float_t cetTimeLHC,xTRKnc, yTRKnc, zTRKnc;
357 fTreeBeamSpot->SetBranchAddress("run", &run);
358 fTreeBeamSpot->SetBranchAddress("cetTimeLHC", &cetTimeLHC);
359 fTreeBeamSpot->SetBranchAddress("bx", &bx);
360 fTreeBeamSpot->SetBranchAddress("triggered", &triggered);
361 fTreeBeamSpot->SetBranchAddress("ntrkletsS", &ntrkletsS);
362 fTreeBeamSpot->SetBranchAddress("xTRKnc", &xTRKnc);
363 fTreeBeamSpot->SetBranchAddress("yTRKnc", &yTRKnc);
364 fTreeBeamSpot->SetBranchAddress("zTRKnc", &zTRKnc);
365 fTreeBeamSpot->SetBranchAddress("ntrksTRKnc", &ntrksTRKnc);
366
367
aea5af08 368 Double_t tstamp = esdE->GetTimeStamp();
369 Float_t cetTime =(tstamp-1262304000.)+7200.;
370
371 Float_t cetTime1h =(tstamp-1262304000.)+3600.;
30b05a14 372
4f0574ad 373 cetTimeLHC = (Float_t)cetTime1h;
374
74917c1f 375 Int_t ntracks = esdE->GetNumberOfTracks();
388ca814 376 Int_t nITS5or6=0,nTPCin=0,nTPCinEta09=0;
74917c1f 377 //printf("Tracks # = %d\n",esdE->GetNumberOfTracks());
388ca814 378 for(Int_t itr=0; itr<ntracks; itr++) {
74917c1f 379 AliESDtrack *t = esdE->GetTrack(itr);
388ca814 380 if(t->GetNcls(0)>=5) nITS5or6++;
74917c1f 381 Double_t z0; t->GetZAt(0,esdE->GetMagneticField(),z0);
382 if(t->GetNcls(1)>0 && TMath::Abs(t->GetD(0,0,esdE->GetMagneticField()))<2.8 && TMath::Abs(z0)<20) {
388ca814 383 nTPCin++;
384 if(TMath::Abs(t->GetTgl())<1.5) nTPCinEta09++;
385 }
388ca814 386 }
387
4f0574ad 388
74917c1f 389 const AliESDVertex *spdv=esdE->GetPrimaryVertexSPD();
39258194 390 const AliESDVertex *spdvp=esdE->GetPileupVertexSPD(0);
74917c1f 391 const AliESDVertex *tpcv=esdE->GetPrimaryVertexTPC();
392 const AliESDVertex *trkv=esdE->GetPrimaryVertexTracks();
b588553c 393
394 // fill histos
395
396 if(spdv) {
397 if(spdv->GetNContributors()>0) {
398 TString title=spdv->GetTitle();
399 if(title.Contains("3D")) {
400 fhSPDVertexX->Fill(spdv->GetXv());
401 fhSPDVertexY->Fill(spdv->GetYv());
402 }
403 fhSPDVertexZ->Fill(spdv->GetZv());
404 }
405 }
406
407 if(trkv) {
408 if(trkv->GetNContributors()>0) {
409 fhTRKVertexX->Fill(trkv->GetXv());
410 fhTRKVertexY->Fill(trkv->GetYv());
411 fhTRKVertexZ->Fill(trkv->GetZv());
412 }
413 }
414
415 if(tpcv) {
416 if(tpcv->GetNContributors()>0) {
417 fhTPCVertexX->Fill(tpcv->GetXv());
418 fhTPCVertexY->Fill(tpcv->GetYv());
419 fhTPCVertexZ->Fill(tpcv->GetZv());
420 }
421 }
39258194 422
a874bb01 423
424
425 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
a874bb01 426 // filling Vertex reco efficiency plots
d0bbaf02 427
a874bb01 428 if(eventTriggered ? 1. : 0.){
429 fhTriggeredTrklets->Fill(ntrklets);
d0bbaf02 430
431 if(spdv){
432 if(spdv->GetNContributors()>0.5){
433 TString spdtitle = spdv->GetTitle();
434 if(spdtitle.Contains("vertexer: 3D") ? 1. : 0.){
435 fhSPD3DTrklets->Fill(ntrklets);
436 fhSPD3DZreco->Fill(spdv->GetZv());
437 }
438 else{
439 fhSPDZTrklets->Fill(ntrklets);
440 fhSPDZZreco->Fill(spdv->GetZv());
441 }
a874bb01 442 }
443 }
d0bbaf02 444
445 if(trkv){
446 if(trkv->GetNContributors()>0.5)fhTRKTrklets->Fill(ntrklets);
447 }
a874bb01 448 if(fRecoVtxITSTPC) {
449 AliESDVertex *trkvc = ReconstructPrimaryVertexITSTPC(kTRUE);
450 if(trkvc->GetNContributors()>0.5)fhTRKcTrklets->Fill(ntrklets);
451 delete trkvc; trkvc=0;
452 AliESDVertex *trkvnc = ReconstructPrimaryVertexITSTPC(kFALSE);
453 if(trkvnc->GetNContributors()>0.5)fhTRKncTrklets->Fill(ntrklets);
454 delete trkvnc; trkvnc=0;
d0bbaf02 455 }
456
457 if(tpcv){
458 if(tpcv->GetNContributors()>0.5)fhTPCTrklets->Fill(ntrklets);
459 }
460 if(fRecoVtxTPC) {
461 AliESDVertex *tpcvc = ReconstructPrimaryVertexTPC(kTRUE);
462 if(tpcvc->GetNContributors()>0.5)fhTPCcTrklets->Fill(ntrklets);
463 delete tpcvc; tpcvc=0;
464 AliESDVertex *tpcvnc = ReconstructPrimaryVertexTPC(kFALSE);
465 if(tpcvnc->GetNContributors()>0.5)fhTPCncTrklets->Fill(ntrklets);
466 delete tpcvnc; tpcvnc=0;
467 }
a874bb01 468 }
469 /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
470
39258194 471 Float_t xpile=-999.;
472 Float_t ypile=-999.;
473 Float_t zpile=-999.;
474 Float_t expile=-999.;
475 Float_t eypile=-999.;
476 Float_t ezpile=-999.;
477 Int_t ntrkspile=-1;
478
b6292968 479 if(esdE->GetNumberOfPileupVerticesSPD()>0 && spdvp && spdv){
a874bb01 480
481 if(spdvp->GetNContributors()>0) {
482
483 fhSPDVertexXPile->Fill(spdvp->GetXv());
484 fhSPDVertexYPile->Fill(spdvp->GetYv());
485 fhSPDContributorsPile->Fill(spdvp->GetNContributors());
486 fhSPDDispContributors->Fill(spdv->GetNContributors(),spdvp->GetNContributors());
487 fhSPDVertexZPile->Fill(spdvp->GetZv());
8d3b93cb 488 if(spdvp->GetNContributors()>=2) {fhSPDVertexDiffZPileContr2->Fill(spdv->GetZv()-spdvp->GetZv());}
489 if(spdvp->GetNContributors()>=3) {fhSPDVertexDiffZPileContr3->Fill(spdv->GetZv()-spdvp->GetZv());}
490 if(spdvp->GetNContributors()>=4) {fhSPDVertexDiffZPileContr4->Fill(spdv->GetZv()-spdvp->GetZv());}
491 if(spdvp->GetNContributors()>=5) {fhSPDVertexDiffZPileContr5->Fill(spdv->GetZv()-spdvp->GetZv());}
492
a874bb01 493 }
494
39258194 495 xpile=spdvp->GetXv();
496 expile=spdvp->GetXRes();
497 ypile=spdvp->GetYv();
498 eypile=spdvp->GetYRes();
499 zpile=spdvp->GetZv();
500 ezpile=spdvp->GetZRes();
501 ntrkspile=spdvp->GetNContributors();
502 }
388ca814 503
b588553c 504 // fill ntuple
19076e2d 505 Int_t isize=82;
506 Float_t xnt[82]; for(Int_t iii=0;iii<isize;iii++) xnt[iii]=0.;
30b05a14 507
aea5af08 508 Int_t isizeSecNt=9;
509 Float_t secnt[9]; for(Int_t iii=0;iii<isizeSecNt;iii++) secnt[iii]=0.;
30b05a14 510
97805f7d 511 Int_t index=0;
30b05a14 512 Int_t indexSecNt=0;
97805f7d 513
74917c1f 514 xnt[index++]=(Float_t)esdE->GetRunNumber();
30b05a14 515 secnt[indexSecNt++]=(Float_t)esdE->GetRunNumber();
4f0574ad 516 run = (Int_t)esdE->GetRunNumber();
8c7764e0 517
30b05a14 518 xnt[index++]=cetTime; //(Float_t)esdE->GetTimeStamp();
aea5af08 519 //secnt[indexSecNt++]=cetTime;
520 secnt[indexSecNt++]=cetTime1h;
8c7764e0 521
30b05a14 522 xnt[index++]=(Float_t)esdE->GetBunchCrossNumber();
4f0574ad 523 secnt[indexSecNt++]=(Float_t)esdE->GetBunchCrossNumber();
524 bx = (Int_t)esdE->GetBunchCrossNumber();
97805f7d 525
30b05a14 526 xnt[index++]=(eventTriggered ? 1. : 0.);
527 secnt[indexSecNt++]=(eventTriggered ? 1. : 0.);
4f0574ad 528 triggered = (UShort_t)(eventTriggered ? 1 : 0);
529
30b05a14 530 xnt[index++]=(Float_t)dNchdy;
531
97805f7d 532 xnt[index++]=mcVertex[0];
533 xnt[index++]=mcVertex[1];
534 xnt[index++]=mcVertex[2];
388ca814 535
97805f7d 536 xnt[index++]=spdv->GetXv();
537 xnt[index++]=spdv->GetXRes();
538 xnt[index++]=spdv->GetYv();
539 xnt[index++]=spdv->GetYRes();
540 xnt[index++]=spdv->GetZv();
541 xnt[index++]=spdv->GetZRes();
542 xnt[index++]=spdv->GetNContributors();
8c7764e0 543 TString spdtitle = spdv->GetTitle();
544 xnt[index++]=(spdtitle.Contains("vertexer: 3D") ? 1. : 0.);
545 xnt[index++]=spdv->GetDispersion();
388ca814 546
97805f7d 547 xnt[index++]=tpcv->GetXv();
548 xnt[index++]=tpcv->GetXRes();
549 xnt[index++]=tpcv->GetYv();
550 xnt[index++]=tpcv->GetYRes();
551 xnt[index++]=tpcv->GetZv();
552 xnt[index++]=tpcv->GetZRes();
553 xnt[index++]=tpcv->GetNContributors();
8c7764e0 554 TString tpctitle = tpcv->GetTitle();
555 xnt[index++]=(tpctitle.Contains("WithConstraint") ? 1. : 0.);
388ca814 556
97805f7d 557 xnt[index++]=trkv->GetXv();
558 xnt[index++]=trkv->GetXRes();
559 xnt[index++]=trkv->GetYv();
560 xnt[index++]=trkv->GetYRes();
561 xnt[index++]=trkv->GetZv();
562 xnt[index++]=trkv->GetZRes();
563 xnt[index++]=trkv->GetNContributors();// tpccontrorig;
8c7764e0 564 TString trktitle = trkv->GetTitle();
565 xnt[index++]=(trktitle.Contains("WithConstraint") ? 1. : 0.);
388ca814 566
97805f7d 567 xnt[index++]=float(ntrklets);
30b05a14 568 secnt[indexSecNt++]=float(ntrklets);
4f0574ad 569 ntrkletsS = (UShort_t)ntrklets;
30b05a14 570
97805f7d 571 xnt[index++]=float(ntracks);
572 xnt[index++]=float(nITS5or6);
30b05a14 573
97805f7d 574 xnt[index++]=float(nTPCin);
575 xnt[index++]=float(nTPCinEta09);
388ca814 576
97805f7d 577 xnt[index++]=spd0cls;
39258194 578
579 xnt[index++]=xpile;
580 xnt[index++]=expile;
581 xnt[index++]=ypile;
582 xnt[index++]=eypile;
583 xnt[index++]=zpile;
584 xnt[index++]=ezpile;
585 xnt[index++]=(Float_t)ntrkspile;
586
587
588
8c7764e0 589 // add recalculated vertices TRK and TPC
97805f7d 590
8c7764e0 591 if(fRecoVtxTPC) {
592 AliESDVertex *tpcvnc = ReconstructPrimaryVertexTPC(kFALSE);
593 xnt[index++]=tpcvnc->GetXv();
594 xnt[index++]=tpcvnc->GetXRes();
595 xnt[index++]=tpcvnc->GetYv();
596 xnt[index++]=tpcvnc->GetYRes();
597 xnt[index++]=tpcvnc->GetZv();
598 xnt[index++]=tpcvnc->GetZRes();
599 xnt[index++]=tpcvnc->GetNContributors();
600 delete tpcvnc; tpcvnc=0;
601
602 AliESDVertex *tpcvc = ReconstructPrimaryVertexTPC(kTRUE);
603 xnt[index++]=tpcvc->GetXv();
604 xnt[index++]=tpcvc->GetXRes();
605 xnt[index++]=tpcvc->GetYv();
606 xnt[index++]=tpcvc->GetYRes();
607 xnt[index++]=tpcvc->GetZv();
608 xnt[index++]=tpcvc->GetZRes();
609 xnt[index++]=tpcvc->GetNContributors();
610 delete tpcvc; tpcvc=0;
39258194 611 } else index+=14;
8c7764e0 612
613 if(fRecoVtxITSTPC) {
614 AliESDVertex *trkvnc = ReconstructPrimaryVertexITSTPC(kFALSE);
615 xnt[index++]=trkvnc->GetXv();
616 xnt[index++]=trkvnc->GetXRes();
617 xnt[index++]=trkvnc->GetYv();
618 xnt[index++]=trkvnc->GetYRes();
619 xnt[index++]=trkvnc->GetZv();
620 xnt[index++]=trkvnc->GetZRes();
621 xnt[index++]=trkvnc->GetNContributors();
30b05a14 622
623 secnt[indexSecNt++]=trkvnc->GetXv();
624 secnt[indexSecNt++]=trkvnc->GetYv();
625 secnt[indexSecNt++]=trkvnc->GetZv();
626 secnt[indexSecNt++]=trkvnc->GetNContributors();
4f0574ad 627
628 xTRKnc = (Float_t)trkvnc->GetXv();
629 yTRKnc = (Float_t)trkvnc->GetYv();
630 zTRKnc = (Float_t)trkvnc->GetZv();
ee98f8eb 631 ntrksTRKnc = (trkvnc->GetNContributors()<0 ? 0 : (UShort_t)trkvnc->GetNContributors());
4f0574ad 632
30b05a14 633 delete trkvnc; trkvnc=0;
634
635
8c7764e0 636 AliESDVertex *trkvc = ReconstructPrimaryVertexITSTPC(kTRUE);
637 xnt[index++]=trkvc->GetXv();
638 xnt[index++]=trkvc->GetXRes();
639 xnt[index++]=trkvc->GetYv();
640 xnt[index++]=trkvc->GetYRes();
641 xnt[index++]=trkvc->GetZv();
642 xnt[index++]=trkvc->GetZRes();
643 xnt[index++]=trkvc->GetNContributors();
644 delete trkvc; trkvc=0;
39258194 645 } else index+=14;
646
4d95fd6c 647 if(fRecoVtxITSTPCHalfEvent) {
648 AliESDVertex *trkvncodd = ReconstructPrimaryVertexITSTPC(kFALSE,1);
649 AliESDVertex *trkvnceven = ReconstructPrimaryVertexITSTPC(kFALSE,2);
650 if(trkvncodd->GetNContributors()>0 &&
651 trkvnceven->GetNContributors()>0) {
652 xnt[index++]=trkvncodd->GetXv()-trkvnceven->GetXv();
4d95fd6c 653 xnt[index++]=trkvncodd->GetYv()-trkvnceven->GetYv();
4d95fd6c 654 xnt[index++]=trkvncodd->GetZv()-trkvnceven->GetZv();
007ec307 655 xnt[index++]=TMath::Sqrt(trkvncodd->GetXRes()*trkvncodd->GetXRes()+trkvnceven->GetXRes()*trkvnceven->GetXRes());
656 xnt[index++]=TMath::Sqrt(trkvncodd->GetYRes()*trkvncodd->GetYRes()+trkvnceven->GetYRes()*trkvnceven->GetYRes());
4d95fd6c 657 xnt[index++]=TMath::Sqrt(trkvncodd->GetZRes()*trkvncodd->GetZRes()+trkvnceven->GetZRes()*trkvnceven->GetZRes());
19076e2d 658 xnt[index++]=trkvnceven->GetNContributors();
659 xnt[index++]=trkvncodd->GetNContributors();
4d95fd6c 660 } else {
661 xnt[index++]=0.;
662 xnt[index++]=0.;
663 xnt[index++]=0.;
664 xnt[index++]=0.;
665 xnt[index++]=0.;
666 xnt[index++]=0.;
667 xnt[index++]=-1;
19076e2d 668 xnt[index++]=-1;
4d95fd6c 669 }
670 delete trkvncodd; trkvncodd=0;
671 delete trkvnceven; trkvnceven=0;
19076e2d 672 } else index+=8;
4d95fd6c 673
674
8c7764e0 675 if(index>isize) printf("AliAnalysisTaskVertexESD: ERROR, index!=isize\n");
b588553c 676 if(fFillNtuple) fNtupleVertexESD->Fill(xnt);
30b05a14 677
678 if(indexSecNt>isizeSecNt) printf("AliAnalysisTaskVertexESD: ERROR, indexSecNt!=isizeSecNt\n");
4f0574ad 679
680 // if(indexTree>isizeTree) printf("AliAnalysisTaskVertexESD: ERROR, indexTree!=isizeTree\n");
681 // only every second event (to reduce output size)
682 if(fFillTreeBeamSpot && (nInFile%2 == 0)) fTreeBeamSpot->Fill();
683
f6232e6e 684
388ca814 685 // Post the data already here
41deb8e3 686 PostData(1, fOutput);
388ca814 687
688 return;
689}
690
691//________________________________________________________________________
692void AliAnalysisTaskVertexESD::Terminate(Option_t *)
693{
694 // Draw result to the screen
695 // Called once at the end of the query
2ac70243 696 fOutput = dynamic_cast<TList*> (GetOutputData(1));
388ca814 697 if (!fOutput) {
698 Printf("ERROR: fOutput not available");
699 return;
700 }
30b05a14 701
a874bb01 702 //////////////////////////////////////////////////////
703 /*
704 TH1F *fhTriggeredTrklets=(TH1F*)fOutput->FindObject("fhTriggeredTrklets");
705 TH1F *fhSPDZTrklets=(TH1F*)fOutput->FindObject("fhSPDZTrklets");
706 TH1F *fhSPD3DTrklets=(TH1F*)fOutput->FindObject("fhSPD3DTrklets");
707 TH1F *fhTRKTrklets=(TH1F*)fOutput->FindObject("fhTRKTrklets");
708 TH1F *fhTRKcTrklets=(TH1F*)fOutput->FindObject("fhTRKcTrklets");
709 TH1F *fhTRKncTrklets=(TH1F*)fOutput->FindObject("fhTRKncTrklets");
710 TH1F *fhSPDZZreco=(TH1F*)fOutput->FindObject("fhSPDZZreco");
711 TH1F *fhSPD3DZreco=(TH1F*)fOutput->FindObject("fhSPD3DZreco");
712
713 TGraphAsymmErrors *fhSPDZEffTrklets=new TGraphAsymmErrors(fhSPDZTrklets,fhTriggeredTrklets,"w");
714 fhSPDZEffTrklets->SetName("fhSPDZEffTrklets");
715 fhSPDZEffTrklets->SetDrawOption("AP");
716 TGraphAsymmErrors *fhSPD3DEffTrklets=new TGraphAsymmErrors(fhSPD3DTrklets,fhTriggeredTrklets,"w");
717 fhSPD3DEffTrklets->SetName("fhSPD3DEffTrklets");
718 TH1F * fhSPDOverallTrklets=(TH1F*)fhSPDZTrklets->Clone("fhSPDOverallTrklets");
719 fhSPDOverallTrklets->Add(fhSPD3DTrklets);
720 TGraphAsymmErrors *fhSPDOverallEffTrklets=new TGraphAsymmErrors(fhSPDOverallTrklets,fhTriggeredTrklets,"w");
721 fhSPDOverallEffTrklets->SetName("fhSPDOverallEffTrklets");
722 TGraphAsymmErrors *fhTRKEffTrklets=new TGraphAsymmErrors(fhTRKTrklets,fhTriggeredTrklets,"w");
723 fhTRKEffTrklets->SetName("fhTRKEffTrklets");
724 TGraphAsymmErrors *fhTRKcEffTrklets=new TGraphAsymmErrors(fhTRKcTrklets,fhTriggeredTrklets,"w");
725 fhTRKcEffTrklets->SetName("fhTRKcEffTrklets");
726 TGraphAsymmErrors *fhTRKncEffTrklets=new TGraphAsymmErrors(fhTRKncTrklets,fhTriggeredTrklets,"w");
727 fhTRKncEffTrklets->SetName("fhTRKncEffTrklets");
728 TH1F * fhSPDOverallZreco=(TH1F*)fhSPDZZreco->Clone("fhSPDOverallZreco");
729 fhSPDOverallZreco->Add(fhSPD3DZreco);
730 TGraphAsymmErrors *fhSPDEffZreco=new TGraphAsymmErrors(fhSPD3DZreco,fhSPDOverallZreco,"w");
731 fhSPDEffZreco->SetName("fhSPDEffZreco");
732
733 TH1F *fhEff = new TH1F("hEff","hEff",6,0.5,6.5);
734 Int_t count=1;
735 if(fhSPDZTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
736 fhEff->Fill(count,fhSPDZTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
737 fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhSPDZTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
738 }
739 fhEff->GetXaxis()->SetBinLabel(count,"SPDZ");
740
741 count++;
742 if(fhSPD3DTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
743 fhEff->Fill(count,fhSPD3DTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
744 fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhSPD3DTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
745 }
746 fhEff->GetXaxis()->SetBinLabel(count,"SPD3D");
747
748 count++;
749 if(fhSPDOverallTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
750 fhEff->Fill(count,fhSPDOverallTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
751 fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhSPDOverallTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
752 }
753 fhEff->GetXaxis()->SetBinLabel(count,"SPD Overall");
754
755 count++;
756 if(fhTRKTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
757 fhEff->Fill(count,fhTRKTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
758 fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhTRKTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
759 }
760 fhEff->GetXaxis()->SetBinLabel(count,"TRK");
761
762 count++;
763 if(fhTRKcTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
764 fhEff->Fill(count,fhTRKcTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
765 fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhTRKcTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
766 }
767 fhEff->GetXaxis()->SetBinLabel(count,"TRKc");
768
769 count++;
770 if(fhTRKncTrklets->GetEntries()!=0 && fhTriggeredTrklets->GetEntries()!=0){
771 fhEff->Fill(count,fhTRKncTrklets->GetEntries()/fhTriggeredTrklets->GetEntries());
772 fhEff->SetBinError(count,fhEff->GetBinContent(count)*TMath::Sqrt(1/fhTRKncTrklets->GetEntries()+1/fhTriggeredTrklets->GetEntries()));
773 }
774 fhEff->GetXaxis()->SetBinLabel(count,"TRKnc");
775
776 count++;
777 fhEff->Print("all");
778
779 TFile* fileEff = new TFile("VtxEff.root","recreate");
780 fhSPDZEffTrklets->Write();
781 fhSPD3DEffTrklets->Write();
782 fhSPDOverallEffTrklets->Write();
783 fhTRKEffTrklets->Write();
784 fhTRKcEffTrklets->Write();
785 fhTRKncEffTrklets->Write();
786 fhSPDEffZreco->Write();
787 fhEff->Write();
788 fileEff->Close();
789 delete fileEff;
790
791 /////////////////////////////////////////
792 */
793
794
30b05a14 795 if (!fNtupleVertexESD){
796 Printf("ERROR: fNtuple not available");
797 return;
798 }
799
388ca814 800 fNtupleVertexESD = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleVertexESD"));
801
30b05a14 802
388ca814 803 return;
804}
805
806//_________________________________________________________________________
8c7764e0 807AliESDVertex* AliAnalysisTaskVertexESD::ReconstructPrimaryVertexTPC(Bool_t constr) const {
388ca814 808 // On the fly reco of TPC vertex from ESD
74917c1f 809 AliESDEvent* evt = (AliESDEvent*) fInputEvent;
810 AliVertexerTracks vertexer(evt->GetMagneticField());
0faeb7ba 811 if(evt->GetNumberOfTracks()<500) {
812 vertexer.SetTPCMode(); // defaults
813 } else {
814 vertexer.SetTPCMode(0.1,1.0,5.0,10,1,3.,0.1,1.5,3.,30.,1,1);// PbPb
815 }
8c7764e0 816 Float_t diamondcovxy[3]; evt->GetDiamondCovXY(diamondcovxy);
817 Double_t pos[3]={evt->GetDiamondX(),evt->GetDiamondY(),0};
818 Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
819 AliESDVertex *initVertex = new AliESDVertex(pos,cov,1.,1);
388ca814 820 vertexer.SetVtxStart(initVertex);
821 delete initVertex;
8c7764e0 822 if(!constr) vertexer.SetConstraintOff();
388ca814 823
74917c1f 824 return vertexer.FindPrimaryVertex(evt);
388ca814 825}
94f15f9b 826
827//_________________________________________________________________________
4d95fd6c 828AliESDVertex* AliAnalysisTaskVertexESD::ReconstructPrimaryVertexITSTPC(Bool_t constr,Int_t mode) const {
94f15f9b 829 // On the fly reco of ITS+TPC vertex from ESD
4d95fd6c 830 // mode=0 use all tracks
831 // mode=1 use odd-number tracks
832 // mode=2 use even-number tracks
833
74917c1f 834 AliESDEvent* evt = (AliESDEvent*) fInputEvent;
835 AliVertexerTracks vertexer(evt->GetMagneticField());
0faeb7ba 836 if(evt->GetNumberOfTracks()<500) {
837 vertexer.SetITSMode(); // defaults
838 vertexer.SetMinClusters(4); // default is 5
839 } else {
840 vertexer.SetITSMode(0.1,0.1,0.5,5,1,3.,100.,1000.,3.,30.,1,1);// PbPb
841 }
30b05a14 842 //vertexer.SetITSpureSA();
8c7764e0 843 Float_t diamondcovxy[3]; evt->GetDiamondCovXY(diamondcovxy);
844 Double_t pos[3]={evt->GetDiamondX(),evt->GetDiamondY(),0};
845 Double_t cov[6]={diamondcovxy[0],diamondcovxy[1],diamondcovxy[2],0.,0.,10.*10.};
846 AliESDVertex *initVertex = new AliESDVertex(pos,cov,1.,1);
94f15f9b 847 vertexer.SetVtxStart(initVertex);
848 delete initVertex;
8c7764e0 849 if(!constr) vertexer.SetConstraintOff();
94f15f9b 850
cb393b39 851 // use only ITS-TPC or only ITS-SA tracks
5457d525 852 if(fOnlyITSTPCTracks || fOnlyITSSATracks || mode>0) {
cb393b39 853 Int_t iskip=0;
74917c1f 854 Int_t *skip = new Int_t[evt->GetNumberOfTracks()];
855 for(Int_t itr=0;itr<evt->GetNumberOfTracks(); itr++) {
856 AliESDtrack* track = evt->GetTrack(itr);
cb393b39 857 if(fOnlyITSTPCTracks && track->GetNcls(1)==0) { // skip ITSSA
858 skip[iskip++]=itr;
4d95fd6c 859 continue;
cb393b39 860 }
861 if(fOnlyITSSATracks && track->GetNcls(1)>0) { // skip ITSTPC
862 skip[iskip++]=itr;
4d95fd6c 863 continue;
cb393b39 864 }
4d95fd6c 865 if(mode==1 && itr%2==0) skip[iskip++]=itr;
866 if(mode==2 && itr%2==1) skip[iskip++]=itr;
cb393b39 867 }
868 vertexer.SetSkipTracks(iskip,skip);
869 delete [] skip; skip=NULL;
870 }
871
74917c1f 872 return vertexer.FindPrimaryVertex(evt);
94f15f9b 873}