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