10 #include "TParticlePDG.h"
13 #include "AliAnalysisTask.h"
14 #include "AliAnalysisManager.h"
16 #include "AliESDEvent.h"
18 #include "AliESDVertex.h"
19 #include "AliESDInputHandler.h"
20 #include "AliESDtrackCuts.h"
21 #include "AliMultiplicity.h"
23 #include "AliMCParticle.h"
24 #include "AliMCEvent.h"
25 #include "AliAnalysisTaskEfficiency.h"
26 #include "AliExternalTrackParam.h"
27 #include "AliTrackReference.h"
29 #include "AliHeader.h"
30 #include "AliGenEventHeader.h"
31 #include "AliGenCocktailEventHeader.h"
32 #include "AliGenPythiaEventHeader.h"
33 #include "AliGenDPMjetEventHeader.h"
35 // Analysis Task for basic QA on the ESD
36 // Authors: Jan Fiete Grosse-Oetringhaus, Christian Klein-Boesing, Veronica Canoa, AM, Eva Sicking
38 ClassImp(AliAnalysisTaskEfficiency)
40 //________________________________________________________________________
41 AliAnalysisTaskEfficiency::AliAnalysisTaskEfficiency(const char *name)
42 : AliAnalysisTaskSE(name)
73 ,fh2EtaCorrelationShift(0)
75 ,fh2PhiCorrelationShift(0)
77 ,fh2PtCorrelationShift(0)
78 ,fHistDeltaphiprimaries(0)
79 ,fHistDeltaphisecondaries(0)
80 ,fHistDeltaphireject(0)
89 ,fHistTPCTRDdeltaphi(0)
109 for(int i = 0;i< 2;++i){
110 fh2MultSpdChips[i] = 0;
113 for(int i = 0;i < kTotalAC;++i){
114 fh2PhiPadRow[i] = fh2PhiLayer[i] = 0;
117 for(int i = 0; i < kTotalVtx; ++i){
118 fh2VertexCorrelation[i] =
119 fh2VertexCorrelationShift[i] = 0;
121 fh1VertexShiftNorm[i] = 0;
124 for(int i = 0;i< 8;++i){
125 fHistRECptCharge[i]=0;
126 fHistMCptCharge[i]=0;
127 fHistMCNRptCharge[i]=0;
128 fHistFAKEptCharge[i]=0;
130 fHistRECetaCharge[i]=0;
131 fHistMCetaCharge[i]=0;
132 fHistMCNRetaCharge[i]=0;
133 fHistFAKEetaCharge[i]=0;
135 fHistRECphiCharge[i]=0;
136 fHistMCphiCharge[i]=0;
137 fHistMCNRphiCharge[i]=0;
138 fHistFAKEphiCharge[i]=0;
142 DefineOutput(1, TList::Class());
146 //________________________________________________________________________
147 void AliAnalysisTaskEfficiency::UserCreateOutputObjects()
152 fHists = new TList();
153 fHistRECpt = new TH1F("fHistRECpt", " p_{T} distribution: all reconstructed", 100, 0., 20.);
154 fHistMCpt = new TH1F("fHistMCpt", " p_{T} distribution: all MC", 100, 0., 20.);
155 fHistMCNRpt = new TH1F("fHistMCNRpt", " p_{T} distribution: all not-reconstructed MC", 100, 0., 20.);
156 fHistFAKEpt = new TH1F("fHistFAKEpt", " p_{T} distribution: all Fake", 100, 0., 20.);
157 fHistMULTpt = new TH1F("fHistMULTpt", " p_{T} distribution: multiply rec.", 100, 0., 20.);
159 fHistRECeta = new TH1F("fHistRECeta", " #eta distribution: all reconstructed", 100,-1.0,1.0);
160 fHistMCeta = new TH1F("fHistMCeta", " #eta distribution: all MC", 100,-1.0,1.0);
161 fHistMCNReta = new TH1F("fHistMCNReta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0);
162 fHistFAKEeta = new TH1F("fHistFAKEeta", " #eta distribution: all Fake", 100,-1.0,1.0);
163 fHistMULTeta = new TH1F("fHistMULTeta", " #eta distribution: multiply rec.", 100,-1.0,1.0);
165 fHistRECphi = new TH1F("fHistRECphi", " #phi distribution: all reconstructed", 314, 0., 6.28);
166 fHistMCphi = new TH1F("fHistMCphi", " #phi distribution: all MC", 314, 0., 6.28);
167 fHistMCNRphi = new TH1F("fHistMCNRphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28);
168 fHistFAKEphi = new TH1F("fHistFAKEphi", " #phi distribution: all Fake", 314, 0., 6.28);
169 fHistMULTphi = new TH1F("fHistMULTphi", " #phi distribution: multipli rec.", 314, 0., 6.28);
171 fHistRECHPTeta = new TH1F("fHistRECHPTeta", " #eta distribution: all reconstructed", 100,-1.0,1.0);
172 fHistMCHPTeta = new TH1F("fHistMCHPTeta", " #eta distribution: all MC", 100,-1.0,1.0);
173 fHistMCNRHPTeta = new TH1F("fHistMCNRHPTeta", " #eta distribution: all not-reconstructed MC", 100,-1.0,1.0);
174 fHistFAKEHPTeta = new TH1F("fHistFAKEHPTeta", " #eta distribution: all Fake", 100,-1.0,1.0);
176 fHistRECHPTphi = new TH1F("fHistRECHPTphi", " #phi distribution: all reconstructed", 314, 0., 6.28);
177 fHistMCHPTphi = new TH1F("fHistMCHPTphi", " #phi distribution: all MC", 314, 0., 6.28);
178 fHistMCNRHPTphi = new TH1F("fHistMCNRHPTphi", " #phi distribution: all not-reconstructed MC", 314, 0., 6.28);
179 fHistFAKEHPTphi = new TH1F("fHistFAKEHPTphi", " #phi distribution: all Fake", 314, 0., 6.28);
181 fHistRecMult = new TH1F("fHistRecMult", "Multiple reconstructed tracks", 50, 0., 50.);
182 fHistNCluster = new TH1F("fHistNCluster", "Number of clusters for suspicious tracks", 300, 0., 300.);
186 fh1VtxEff = new TH1F("fh1VtxEff","VtxEff TPC",4,-0.5,3.5);
187 fh1VtxEff->GetXaxis()->SetBinLabel(1,"NO TPC VTX");
188 fh1VtxEff->GetXaxis()->SetBinLabel(2,"TPC VTX");
189 fh1VtxEff->GetXaxis()->SetBinLabel(3,"NO SPD VTX");
190 fh1VtxEff->GetXaxis()->SetBinLabel(4,"SPD VTX");
192 TString labels[kTotalAC];
193 labels[kNegA]="NegA";
194 labels[kPosA]="PosA";
195 labels[kNegC]="NegC";
196 labels[kPosC]="PosC";
198 for(int i = 0;i< kTotalAC;++i){
199 fh2PhiPadRow[i] = new TH2F(Form("fh2PhiPadRow_%d",i),Form("Padrow vs phi AC %s;phi;padrow",labels[i].Data()),360,0.,360.,159,-0.5,158.5);
200 fh2PhiLayer[i] = new TH2F(Form("fh2PhiLayer_%d",i),Form("layer vs phi AC %s;phi;layer",labels[i].Data()),360,0.,360.,6,-0.5,5.5);
202 for(int i = 0;i<2;++i){
203 fh2MultSpdChips[i] = new TH2F(Form("fh2ChipsSpdMult_%d",i),"mult SPD vs chips;chips;Mult",201,-1.5,199.5,201,-1.5,199.5);
205 fh2VtxTpcSpd = new TH2F("fh2VtxTpcSpd","SPD vtx vs TPC;tpc-z;spd-z",200,-20.,20.,200,-20,20.);
207 TString labelsVtx[kTotalVtx];
208 labelsVtx[kTPC] = "TPC";
209 labelsVtx[kSPD] = "SPD";
210 for(int i = 0; i < kTotalVtx; ++i){
211 fh2VertexCorrelation[i] = new TH2F(Form("fh2VertexCorrelation_%d",i),Form("VertexCorrelation %s;MC z-vtx;ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 120, -30, 30);
212 fh2VertexCorrelationShift[i] = new TH2F(Form("fh2VertexCorrelationShift_%d",i), Form("VertexCorrelationShift %s;MC z-vtx;MC z-vtx - ESD z-vtx",labelsVtx[i].Data()), 120, -30, 30, 100, -1, 1);
213 fh1VertexShift[i] = new TH1F(Form("fh1VertexShift_%d",i), Form("VertexShift %s;(MC z-vtx - ESD z-vtx);Entries",labelsVtx[i].Data()), 201, -2, 2);
214 fh1VertexShiftNorm[i] = new TH1F(Form("fh1VertexShiftNorm_%d",i), Form("VertexShiftNorm %s;(MC z-vtx - ESD z-vtx) / #sigma_{ESD z-vtx};Entries",labelsVtx[i].Data()), 200, -100, 100);
218 fh2EtaCorrelation = new TH2F("fh2EtaCorrelation", "EtaCorrelation;MC #eta;ESD #eta", 120, -3, 3, 120, -3, 3);
219 fh2EtaCorrelationShift = new TH2F("fh2EtaCorrelationShift", "EtaCorrelationShift;ESD #eta;MC #eta - ESD #eta", 120, -3, 3, 100, -0.1, 0.1);
220 fh2PhiCorrelation = new TH2F("fh2PhiCorrelation", "PhiCorrelation;MC #phi;ESD #phi", 120, 0, 360, 120, 0, 360);
221 fh2PhiCorrelationShift = new TH2F("fh2PhiCorrelationShift", "PhiCorrelationShift;ESD #phi;MC #phi - ESD #phi", 120, 0, 360, 100, -10, 10);
222 fh2PtCorrelation = new TH2F("fh2PtCorrelation", "PtCorrelation;MC p_T;ESD p_T", 200, 0., 200., 200, 0, 200);
223 fh2PtCorrelationShift = new TH2F("fh2PtCorrelationShift", "PtCorrelationShift;ESD p_T;MC p_T - ESD #p_T", 120, 0, 10, 100, -2, 2);
225 fHistDeltaphiprimaries =new TH1F("fHistDeltaphiprimaries", " #Delta#phi distribution: primaries", 314, -0.2,0.2);
226 fHistDeltaphisecondaries=new TH1F("fHistDeltaphisecondaries", "#Delta#phi distribution: secondaries", 314, -0.2,0.2);
227 fHistDeltaphireject =new TH1F("fHistDeltaphireject", " #Delta#phi distribution: reject trackled", 314, -0.2,0.2);
229 fHistTPCITSdeltax =new TH1F("fHistTPCITSdeltax", "TPC-ITS matching:#Delta x", 100,-20,20);
230 fHistTPCITSdeltay =new TH1F("fHistTPCITSdeltay", "TPC-ITS matching:#Delta y", 100,-20,20);
231 fHistTPCITSdeltaz =new TH1F("fHistTPCITSdeltaz", "TPC-ITS matching:#Delta z", 100,-20,20);
232 fHistTPCITSdeltar =new TH1F("fHistTPCITSdeltar", "TPC-ITS matching:#Delta r", 100,0,20);
234 fHistTPCTRDdeltax = new TH1F("fHistTPCTRDdeltax", "TPC-TRD matching:#Delta x", 400,-400, 400);
235 fHistTPCTRDdeltay = new TH1F("fHistTPCTRDdeltay", "TPC-TRD matching:#Delta y", 400,-400, 400);
236 fHistTPCTRDdeltaz = new TH1F("fHistTPCTRDdeltaz", "TPC-TRD matching:#Delta z", 400,-400, 400);
237 fHistTPCTRDdeltar = new TH1F("fHistTPCTRDdeltar", "TPC-TRD matching:#Delta r", 100,0,20);
238 fHistTPCTRDdeltaphi = new TH1F("fHistTPCTRDdeltarphi", "TPC-TRD matching:#Delta #phi", 128,-3.14 , 3.14);
240 fPtPullsPos = new TProfile("fPtPullsPos", "Pt Pulls for pos. primaries", 50 , 0., 10., -5., 5.,"S");
241 fPtPullsNeg = new TProfile("fPtPullsNeg", "Pt Pulls for neg. primaries", 50 , 0., 10., -5., 5.,"S");
242 fPtShiftPos = new TProfile("fPtShiftPos", "Pt Shift for pos. primaries", 50 , 0., 10., -0.2, 0.2,"S");
243 fPtShiftNeg = new TProfile("fPtShiftNeg", "Pt Shift for neg. primaries", 50 , 0., 10., -0.2, 0.2,"S");
245 fEtaMultiFluc = new TProfile("fEtaMultiFluc", "eta multiplicity fluctuations", 10, -2., 2., 0., 600., "S");
246 fEtaMulti = new TH1F("fEtaMulti", "eta multiplicity fluctuations", 10, -2., 2.);
247 fEtaMultiH = new TH1F("fEtaMultiH", "eta multiplicity fluctuations", 600, 0., 600.);
249 fPhiMultiFluc = new TProfile("fPhiMultiFluc", "phi multiplicity fluctuations", 64, 0., 6.4, 0., 100., "S");
250 fPhiMulti = new TH1F("fPhiMulti", "phi multiplicity fluctuations", 64, 0., 6.4);
251 fPhiMultiH = new TH1F("fPhiMultiH", "phi multiplicity fluctuations", 100, 0., 100.);
254 fVertexRvsZ = new TH2F("fVertexRvsZ", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.);
255 fVertexRvsZC = new TH2F("fVertexRvsZC", "SPD Vertex R vs Z", 200, -1., 1., 200, -1., 1.);
257 TString charge[8]={"Deu", "AntiDeu", "Tri", "AntiTri", "He3", "AntiHe3", "He4", "AntiHe4"};
258 for(Int_t i=0;i<8;i++){
259 fHistRECptCharge[i] = new TH1F(Form("fHistRECptCharge%s",charge[i].Data()),
260 "p_{T} distribution: all reconstructed",
262 fHistMCptCharge[i] = new TH1F(Form("fHistMCptCharge%s",charge[i].Data()),
263 "p_{T} distribution: all MC",
265 fHistMCNRptCharge[i] = new TH1F(Form("fHistMCNRptCharge%s",charge[i].Data()),
266 "p_{T} distribution: all not-reconstructed MC",
268 fHistFAKEptCharge[i] = new TH1F( Form("fHistFAKEptCharge%s",charge[i].Data()),
269 "p_{T} distribution: all Fake",
272 fHistRECetaCharge[i] = new TH1F(Form("fHistRECetaCharge%s",charge[i].Data()),
273 "p_{T} distribution: all reconstructed",
275 fHistMCetaCharge[i] = new TH1F(Form("fHistMCetaCharge%s",charge[i].Data()),
276 "p_{T} distribution: all MC",
278 fHistMCNRetaCharge[i] = new TH1F(Form("fHistMCNRetaCharge%s",charge[i].Data()),
279 "p_{T} distribution: all not-reconstructed MC",
281 fHistFAKEetaCharge[i] = new TH1F( Form("fHistFAKEetaCharge%s",charge[i].Data()),
282 "p_{T} distribution: all Fake",
285 fHistRECphiCharge[i] = new TH1F(Form("fHistRECphiCharge%s",charge[i].Data()),
286 "p_{T} distribution: all reconstructed",
288 fHistMCphiCharge[i] = new TH1F(Form("fHistMCphiCharge%s",charge[i].Data()),
289 "p_{T} distribution: all MC",
291 fHistMCNRphiCharge[i] = new TH1F(Form("fHistMCNRphiCharge%s",charge[i].Data()),
292 "p_{T} distribution: all not-reconstructed MC",
294 fHistFAKEphiCharge[i] = new TH1F( Form("fHistFAKEphiCharge%s",charge[i].Data()),
295 "p_{T} distribution: all Fake",
304 fHists->Add(fHistRECpt);
305 fHists->Add(fHistMCpt);
306 fHists->Add(fHistMCNRpt);
307 fHists->Add(fHistFAKEpt);
308 fHists->Add(fHistMULTpt);
310 fHists->Add(fHistRECeta);
311 fHists->Add(fHistMCeta);
312 fHists->Add(fHistMCNReta);
313 fHists->Add(fHistFAKEeta);
314 fHists->Add(fHistMULTeta);
316 fHists->Add(fHistRECphi);
317 fHists->Add(fHistMCphi);
318 fHists->Add(fHistMCNRphi);
319 fHists->Add(fHistFAKEphi);
320 fHists->Add(fHistMULTphi);
322 fHists->Add(fHistRECHPTeta);
323 fHists->Add(fHistMCHPTeta);
324 fHists->Add(fHistMCNRHPTeta);
325 fHists->Add(fHistFAKEHPTeta);
327 fHists->Add(fHistRECHPTphi);
328 fHists->Add(fHistMCHPTphi);
329 fHists->Add(fHistMCNRHPTphi);
330 fHists->Add(fHistFAKEHPTphi);
334 fHists->Add(fh1VtxEff);
335 for(int i = 0;i < kTotalAC;++i){
336 fHists->Add(fh2PhiPadRow[i]);
337 fHists->Add(fh2PhiLayer[i]);
339 for(int i = 0;i<2;++i){
340 fHists->Add(fh2MultSpdChips[i]);
342 fHists->Add(fh2VtxTpcSpd);
344 for(int i = 0; i < kTotalVtx; ++i){
345 fHists->Add(fh2VertexCorrelation[i]);
346 fHists->Add(fh2VertexCorrelationShift[i]);
347 fHists->Add(fh1VertexShift[i]);
348 fHists->Add(fh1VertexShiftNorm[i]);
352 fHists->Add(fh2EtaCorrelation);
353 fHists->Add(fh2EtaCorrelationShift);
354 fHists->Add(fh2PhiCorrelation);
355 fHists->Add(fh2PhiCorrelationShift);
356 fHists->Add(fh2PtCorrelation);
357 fHists->Add(fh2PtCorrelationShift);
359 fHists->Add(fHistDeltaphiprimaries);
360 fHists->Add(fHistDeltaphisecondaries);
361 fHists->Add(fHistDeltaphireject);
363 fHists->Add(fHistTPCITSdeltax);
364 fHists->Add(fHistTPCITSdeltay);
365 fHists->Add(fHistTPCITSdeltaz);
366 fHists->Add(fHistTPCITSdeltar);
368 fHists->Add(fHistTPCTRDdeltax);
369 fHists->Add(fHistTPCTRDdeltay);
370 fHists->Add(fHistTPCTRDdeltaz);
371 fHists->Add(fHistTPCTRDdeltar);
372 fHists->Add(fHistTPCTRDdeltaphi);
374 fHists->Add(fHistRecMult);
375 fHists->Add(fHistNCluster);
378 fHists->Add(fPtPullsPos);
379 fHists->Add(fPtPullsNeg);
380 fHists->Add(fPtShiftPos);
381 fHists->Add(fPtShiftNeg);
383 fHists->Add(fVertexRvsZ);
384 fHists->Add(fVertexRvsZC);
385 fHists->Add(fEtaMultiFluc);
386 fHists->Add(fEtaMulti);
387 fHists->Add(fEtaMultiH);
388 fHists->Add(fPhiMultiFluc);
389 fHists->Add(fPhiMulti);
390 fHists->Add(fPhiMultiH);
392 for(Int_t i=0;i<8;i++){
393 fHists->Add(fHistRECptCharge[i]);
394 fHists->Add(fHistMCptCharge[i]);
395 fHists->Add(fHistMCNRptCharge[i]);
396 fHists->Add(fHistFAKEptCharge[i]);
398 fHists->Add(fHistRECetaCharge[i]);
399 fHists->Add(fHistMCetaCharge[i]);
400 fHists->Add(fHistMCNRetaCharge[i]);
401 fHists->Add(fHistFAKEetaCharge[i]);
403 fHists->Add(fHistRECphiCharge[i]);
404 fHists->Add(fHistMCphiCharge[i]);
405 fHists->Add(fHistMCNRphiCharge[i]);
406 fHists->Add(fHistFAKEphiCharge[i]);
409 for (Int_t i=0; i<fHists->GetEntries(); ++i) {
410 TH1 *h1 = dynamic_cast<TH1*>(fHists->At(i));
423 //________________________________________________________________________
424 void AliAnalysisTaskEfficiency::UserExec(Option_t *)
427 // Analysis of MC true particles and reconstructed tracks
428 // Different track types (Global, TPC, ITS) can be selected
431 Printf("ERROR: fESD not available");
434 static Int_t nfc = 0;
436 //AliESDInputHandler* esdI = (AliESDInputHandler*)
437 //(AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler();
438 //AliESDEvent* hltEvent = esdI->GetHLTEvent();
440 AliStack* stack = MCEvent()->Stack();
441 AliESDEvent* esdE = (AliESDEvent*) fInputEvent;
444 // Fetch the SPD vertex:
446 const AliESDVertex* spdVertex = esdE->GetPrimaryVertexSPD();
448 if(spdVertex->GetNContributors()<=0){
452 spdVertex->GetXYZ(vtxSPD);
456 // Fetch the TPC vertex
458 const AliESDVertex* tpcVertex = esdE->GetPrimaryVertexTPC();
460 if(tpcVertex->GetNContributors()<=0){
464 tpcVertex->GetXYZ(vtxTPC);
470 x = spdVertex->GetX() + 0.07;
471 y = spdVertex->GetY() - 0.25;
472 z = spdVertex->GetZ();
473 if (TMath::Abs(z) < 10.) {
474 fVertexRvsZ->Fill(x,y);
475 if (TMath::Sqrt(x*x + y*y) > 5. * 0.028)
476 fVertexRvsZC->Fill(x,y);
480 //Printf("%s:%d %5d",(char*)__FILE__,__LINE__, esdE->GetNumberOfTracks());
481 TArrayI labels(esdE->GetNumberOfTracks());
483 // Track loop to fill a pT spectrum
484 //printf("ESD track loop \n");
486 AliESDtrack *tpcP = 0x0;
488 for (Int_t iTracks = 0; iTracks < esdE->GetNumberOfTracks(); iTracks++) {
490 //prevent mem leak for TPConly track
491 if(fTrackType==2&&tpcP){
496 AliVParticle *track = esdE->GetTrack(iTracks);
497 AliESDtrack *esdtrack = dynamic_cast<AliESDtrack*>(track);
498 esdtrack->PropagateToDCA(esdE->GetPrimaryVertex(),
499 esdE->GetMagneticField(), 10000.);
502 Printf("ERROR: Could not receive track %d", iTracks);
506 // run Task for global tracks or ITS tracks or TPC tracks
507 const AliExternalTrackParam *tpcPin = 0x0;
512 Int_t indexAC = GetIndexAC(esdtrack);
514 //select in the steering macro, which track type should be analysed.
515 // Four track types are selectable:
517 // 1. ITS tracks (SA or Pure SA)
518 // 2. TPC only tracks
519 // Track selection copied from PWGPP/AliAnalysisTaskQAsym
522 //Fill all histograms with global tracks
525 if (!fCuts->AcceptTrack(tpcP)) continue;
530 else if(fTrackType==1){
531 //Fill all histograms with ITS tracks
534 if (!fCuts->AcceptTrack(tpcP)) continue;
539 else if(fTrackType==2){
540 //Fill all histograms with TPC track information
541 tpcPin = esdtrack->GetInnerParam();
542 if (!tpcPin) continue;
543 tpcP = AliESDtrackCuts::
544 GetTPCOnlyTrack(dynamic_cast<AliESDEvent*>(esdE),
547 if (!fCuts->AcceptTrack(tpcP)) continue;
548 if(tpcP->GetNcls(1)>160)continue;//jacek's track cut
549 if(tpcP->GetConstrainedChi2TPC()<0)continue; // jacek's track cut
555 Printf("ERROR: wrong track type \n");
562 labels.AddAt(-1, iTracks);
566 if (TMath::Abs(eta) > 0.9) {
569 //Int_t charge=track->Charge();
571 Int_t ind = TMath::Abs(esdtrack->GetLabel());
572 AliMCParticle* mcPtest = (AliMCParticle*) MCEvent()->GetTrack(ind);
573 if (stack->IsPhysicalPrimary(mcPtest->Label())){
574 if(TMath::Abs(mcPtest->PdgCode())>=1000020030){//all helium (+-1000020030,+-1000020040)
575 pt*=2;// reconstruction takes charge=1 -> for particles with charge=2, pT,rec need to be corrected
578 Int_t index=ConvertePDG(mcPtest->PdgCode());
579 if(fDebug>1)Printf("PDG=%d, index=%d", mcPtest->PdgCode(), index);
580 //fill tracks comming from d,t,3He, 4He
582 fHistRECptCharge [index] ->Fill(pt);
583 fHistRECetaCharge[index]->Fill(eta);
584 fHistRECphiCharge[index]->Fill(phi);
588 fHistRECpt ->Fill(pt);
589 fHistRECeta->Fill(eta);
590 fHistRECphi->Fill(phi);
593 fHistRECHPTeta->Fill(eta);
594 fHistRECHPTphi->Fill(phi);
597 // Int_t ind = TMath::Abs(esdtrack->GetLabel());
598 AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(ind);
599 if (!(stack->IsPhysicalPrimary(mcP->Label()))) {
600 fHistFAKEpt ->Fill(pt);
601 fHistFAKEeta->Fill(eta);
602 fHistFAKEphi->Fill(phi);
605 fHistFAKEHPTeta->Fill(eta);
606 fHistFAKEHPTphi->Fill(phi);
611 labels.AddAt(TMath::Abs(esdtrack->GetLabel()), iTracks);
615 Float_t phiDeg = TMath::RadToDeg()*track->Phi();
621 Int_t labeltpcits = esdtrack->GetLabel();
622 AliExternalTrackParam * tpcPCopy = new AliExternalTrackParam(*tpcP);
623 Double_t xk = 43.6; // cm
624 Double_t bz = 5.0; // kG
625 if(tpcPCopy->PropagateTo(xk,bz)){
626 Double_t alpha=tpcPCopy->GetAlpha();
627 // if(tpcPCopy->Rotate(0.)){
628 Float_t xtpc=tpcPCopy->GetX();
629 Float_t ytpc=tpcPCopy->GetY();
630 Float_t ztpc=tpcPCopy->GetZ();
632 xpr=xtpc*TMath::Cos(alpha)-ytpc*TMath::Sin(alpha);
633 ypr=xtpc*TMath::Sin(alpha)+ytpc*TMath::Cos(alpha);
635 AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpcits));
636 Int_t ntref = mcPart->GetNumberOfTrackReferences();
638 for (Int_t k = 0; k < ntref; k++) {
639 AliTrackReference* ref = mcPart->GetTrackReference(k);
640 Float_t xhits = ref->X();
641 Float_t yhits = ref->Y();
642 Float_t radio = TMath::Sqrt(xhits*xhits+yhits*yhits);
644 if(ref->DetectorId() == 0 && (radio > 42.8 && radio < 43.7)) {
646 Float_t xits = ref->X();
647 Float_t yits = ref->Y();
648 Float_t zits = ref->Z();
650 Float_t difx=(xits-xpr);
651 fHistTPCITSdeltax->Fill(difx);
652 Float_t dify=(yits-ypr);
653 fHistTPCITSdeltay->Fill(dify);
654 Float_t difz=(zits-ztpc);
655 fHistTPCITSdeltaz->Fill(difz);
656 Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz);
657 fHistTPCITSdeltar->Fill(deltar);
666 const AliExternalTrackParam *trd = esdtrack->GetOuterParam();
668 if (trd){Int_t labeltpctrd = track->GetLabel();
670 AliExternalTrackParam * trdCopy = new AliExternalTrackParam(*trd);
671 Double_t xktrd=296.0;
673 if(trdCopy->PropagateTo(xktrd,bztrd)){
674 Float_t xtpc2=trdCopy->GetX();
675 Float_t ytpc2=trdCopy->GetY();
676 Float_t ztpc2=trdCopy->GetZ();
677 Double_t alpha=trdCopy->GetAlpha();
679 xpr=xtpc2*TMath::Cos(alpha)-ytpc2*TMath::Sin(alpha);
680 ypr=xtpc2*TMath::Sin(alpha)+ytpc2*TMath::Cos(alpha);
681 AliMCParticle* mcPart = (AliMCParticle*) MCEvent()->GetTrack(abs(labeltpctrd));
682 Int_t ntreftrd = mcPart->GetNumberOfTrackReferences();
684 for (Int_t k = 0; k < ntreftrd; k++) {
685 AliTrackReference* ref = mcPart->GetTrackReference(k);
686 if(ref->DetectorId() == 3 && ref->Label() == labeltpctrd){
688 Float_t xtrd = ref->X();
689 Float_t ytrd = ref->Y();
690 Float_t ztrd = ref->Z();
692 Float_t difx=(xtrd-xpr);
693 fHistTPCTRDdeltax->Fill(difx);
694 Float_t dify=(ytrd-ypr);
695 fHistTPCTRDdeltay->Fill(dify);
696 Float_t difz=(ztrd-ztpc2);
697 fHistTPCTRDdeltaz->Fill(difz);
698 Float_t deltar = TMath::Sqrt(difx * difx + dify * dify + difz * difz);
699 fHistTPCTRDdeltar->Fill(deltar);
700 Float_t phi_tpc = TMath::ATan2(ypr, xpr);
701 Float_t phi_trd = TMath::ATan2(ytrd, xtrd);
702 Float_t dphi = phi_trd - phi_tpc;
703 if (dphi > TMath::Pi()) dphi -= 2. * TMath::Pi();
704 if (dphi < - TMath::Pi()) dphi += 2. * TMath::Pi();
706 fHistTPCTRDdeltaphi->Fill(dphi);
717 TBits bits(esdtrack->GetTPCClusterMap());
718 for(unsigned int ip = 0;ip < bits.GetNbits();++ip){
720 fh2PhiPadRow[indexAC]->Fill(phiDeg,ip);
723 for(int i = 0;i < 6;++i){
724 if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i);
727 else if(!fFieldOn){ // field not on all momenta are set to MPV
728 TBits bits(esdtrack->GetTPCClusterMap());
729 for(unsigned int ip = 0;ip < bits.GetNbits();++ip){
731 fh2PhiPadRow[indexAC]->Fill(phiDeg,ip);
733 for(int i = 0;i < 6;++i){
734 if(esdtrack->HasPointOnITSLayer(i))fh2PhiLayer[indexAC]->Fill(phiDeg,i);
739 // Fill the correlation
741 TParticle *part = MCEvent()->Stack()->Particle(TMath::Abs(esdtrack->GetLabel()));
743 Float_t mcEta = part->Eta();
744 Float_t mcPhi = TMath::RadToDeg()*part->Phi();
745 Float_t mcPt = part->Pt();
746 // if (pt - mcPt > 20.) {
747 fh2EtaCorrelation->Fill(mcEta,eta);
748 fh2EtaCorrelationShift->Fill(eta,mcEta-eta);
749 fh2PhiCorrelation->Fill(mcPhi,phiDeg);
750 fh2PhiCorrelationShift->Fill(phiDeg,mcPhi-phiDeg);
751 fh2PtCorrelation->Fill(mcPt,pt);
752 fh2PtCorrelationShift->Fill(pt,mcPt-pt);
754 Double_t sigmaPt = TMath::Sqrt(esdtrack->GetSigma1Pt2());
755 if (track->Charge() > 0.) {
756 fPtPullsPos->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt);
757 fPtShiftPos->Fill(mcPt, (1./pt - 1./mcPt));
759 fPtPullsNeg->Fill(mcPt, (1./pt - 1./mcPt) / sigmaPt);
760 fPtShiftNeg->Fill(mcPt, (1./pt - 1./mcPt));
769 //prevent mem leak for TPConly track
770 if(fTrackType==2&&tpcP){
775 //Printf("%s:%d",(char*)__FILE__,__LINE__);
779 fh1VtxEff->Fill("TPC VTX",1);
782 fh1VtxEff->Fill("NO TPC VTX",1);
786 fh1VtxEff->Fill("SPD VTX",1);
789 fh1VtxEff->Fill("NO SPD VTX",1);
793 // Vertex correlation SPD vs. TPC
794 if(tpcVertex&&spdVertex){
795 fh2VtxTpcSpd->Fill(vtxTPC[2],vtxSPD[2]);
797 // Printf("%s:%d",(char*)__FILE__,__LINE__);
798 // Multiplicity checks in the SPD
799 const AliMultiplicity *spdMult = esdE->GetMultiplicity();
800 // Multiplicity Analysis
801 Int_t nt = spdMult->GetNumberOfTracklets();
805 for (Int_t ib = 1; ib <= 10; ib++) {
806 Int_t ic = Int_t(fEtaMulti->GetBinContent(ib));
807 Float_t eta = -1.8 + Float_t(ib-1) * 0.4;
808 fEtaMultiFluc->Fill(eta, Float_t(ic));
809 if (ib == 5) fEtaMultiH->Fill(Float_t(ic));
812 for (Int_t ib = 1; ib <= 64; ib++) {
813 Int_t ic = Int_t(fPhiMulti->GetBinContent(ib));
814 Float_t phi = 0.05 + Float_t(ib-1) * 0.1;
815 fPhiMultiFluc->Fill(phi, Float_t(ic));
816 if (ib == 2) fPhiMultiH->Fill(Float_t(ic));
823 for (Int_t j = 0; j < nt; j++) {
824 fEtaMulti->Fill(spdMult->GetEta(j));
825 fPhiMulti->Fill(spdMult->GetPhi(j));
828 Short_t chips0 = spdMult->GetNumberOfFiredChips(0);
829 Short_t chips1 = spdMult->GetNumberOfFiredChips(1);
830 //Printf("%s:%d",(char*)__FILE__,__LINE__);
831 Int_t inputCount = 0;
834 // get multiplicity from ITS tracklets
835 for (Int_t i=0; i<spdMult->GetNumberOfTracklets(); ++i){
836 Float_t deltaPhi = spdMult->GetDeltaPhi(i);
837 // prevent values to be shifted by 2 Pi()
838 if (deltaPhi < -TMath::Pi())
839 deltaPhi += TMath::Pi() * 2;
840 if (deltaPhi > TMath::Pi())
841 deltaPhi -= TMath::Pi() * 2;
842 if (TMath::Abs(deltaPhi) > 1)
843 printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, spdMult->GetTheta(i), spdMult->GetPhi(i), deltaPhi);
845 Int_t label1=spdMult->GetLabel(i,0);
846 Int_t label2=spdMult->GetLabel(i,1);
847 if ( label1==label2){
848 if(stack->IsPhysicalPrimary(label1) == kTRUE)
849 fHistDeltaphiprimaries->Fill(deltaPhi);
851 fHistDeltaphisecondaries->Fill(deltaPhi);
854 fHistDeltaphireject->Fill(deltaPhi);
858 // Printf("%s:%d",(char*)__FILE__,__LINE__);
859 fh2MultSpdChips[0]->Fill(chips0,inputCount);
860 fh2MultSpdChips[1]->Fill(chips1,inputCount);
861 //Printf("%s:%d",(char*)__FILE__,__LINE__);
866 Int_t nmc = MCEvent()->GetNumberOfTracks();
867 AliGenEventHeader* header = MCEvent()->GenEventHeader();
869 AliGenCocktailEventHeader* cH = dynamic_cast<AliGenCocktailEventHeader*> (header);
870 AliGenPythiaEventHeader* pH;
874 pH = dynamic_cast<AliGenPythiaEventHeader*>(header);
876 TObject* entry = cH->GetHeaders()->FindObject("Pythia");
877 pH = dynamic_cast<AliGenPythiaEventHeader*>(entry);
881 if (pH) iproc = pH->ProcessType();
884 header->PrimaryVertex(mcV);
885 Float_t vzMC = mcV[2];
888 // Fill the correlation with MC
889 if(spdVertex&&MCEvent()){
890 Float_t diff = vzMC-vtxSPD[2];
891 fh2VertexCorrelation[kSPD]->Fill(vzMC,vtxSPD[2]);
892 fh2VertexCorrelationShift[kSPD]->Fill(vzMC,diff);
893 fh1VertexShift[kSPD]->Fill(diff);
894 if(spdVertex->GetZRes()>0)fh1VertexShiftNorm[kSPD]->Fill(diff/spdVertex->GetZRes());
896 if(tpcVertex&&MCEvent()){
897 Float_t diff = vzMC-vtxTPC[2];
898 fh2VertexCorrelation[kTPC]->Fill(vzMC,vtxTPC[2]);
899 fh2VertexCorrelationShift[kTPC]->Fill(vzMC,diff);
900 fh1VertexShift[kTPC]->Fill(diff);
901 if(tpcVertex->GetZRes()>0)fh1VertexShiftNorm[kTPC]->Fill(diff/tpcVertex->GetZRes());
907 //printf("MC particle loop \n");
908 for (Int_t i = 0; i < nmc; i++) {
909 AliMCParticle* mcP = (AliMCParticle*) MCEvent()->GetTrack(i);
910 //printf("MC particle loop %5d \n", i);
912 if (!(stack->IsPhysicalPrimary(mcP->Label()))) continue;
913 if (mcP->Particle()->GetPDG()->Charge() == 0) continue;
914 // Int_t charge=Int_t(mcP->Particle()->GetPDG()->Charge());
915 Float_t phi = mcP->Phi();
916 Float_t eta = mcP->Eta();
917 Float_t pt = mcP->Pt();
918 if (TMath::Abs(eta) > 0.9) continue;
919 Int_t ntr = mcP->GetNumberOfTrackReferences();
926 for (Int_t j = 0; j < ntr; j++) {
928 AliTrackReference* ref = mcP->GetTrackReference(j);
930 if (ref->DetectorId() == 0) nITS++;
931 if (ref->DetectorId() == 1) nTPC++;
932 if (ref->DetectorId() == 2) nFRA++;
940 fHistMCpt ->Fill(pt);
941 fHistMCeta->Fill(eta);
942 fHistMCphi->Fill(phi);
943 Int_t index=ConvertePDG(mcP->PdgCode());
945 fHistMCptCharge [index] ->Fill(pt);
946 fHistMCetaCharge[index]->Fill(eta);
947 fHistMCphiCharge[index]->Fill(phi);
952 fHistMCHPTeta->Fill(eta);
953 fHistMCHPTphi->Fill(phi);
956 Bool_t reco = kFALSE;
960 for (Int_t j = 0; j < esdE->GetNumberOfTracks(); j++) {
961 if (i == labels.At(j)) {
964 //AliESDtrack* test = esdE->GetTrack(j);
967 AliESDtrack* track = esdE->GetTrack(jold);
968 nclus = track->GetTPCclusters(0);
969 fHistNCluster->Fill(nclus);
972 track = esdE->GetTrack(j);
973 nclus = track->GetTPCclusters(0);
974 fHistNCluster->Fill(nclus);
982 fHistRecMult->Fill(Float_t(multR), 1.);
984 fHistMCNRpt ->Fill(pt);
985 fHistMCNReta->Fill(eta);
986 fHistMCNRphi->Fill(phi);
988 fHistMCNRptCharge [index] ->Fill(pt);
989 fHistMCNRetaCharge[index]->Fill(eta);
990 fHistMCNRphiCharge[index]->Fill(phi);
995 fHistMCNRHPTeta->Fill(eta);
996 fHistMCNRHPTphi->Fill(phi);
1001 fHistMULTpt ->Fill(pt);
1002 fHistMULTeta->Fill(eta);
1003 fHistMULTphi->Fill(phi);
1005 } // MC particle loop
1006 //printf("End of MC particle loop \n");
1011 //________________________________________________________________________
1012 void AliAnalysisTaskEfficiency::Terminate(Option_t *)
1018 Bool_t AliAnalysisTaskEfficiency::SelectJouri(AliESDtrack* track)
1021 Bool_t selected = kTRUE;
1022 AliESDEvent* esdE = (AliESDEvent*) fInputEvent;
1024 if (track->GetTPCclusters(0) < 50) selected = kFALSE;
1026 const AliESDVertex *vtx=esdE->GetPrimaryVertexSPD();
1027 if (!vtx->GetStatus()) selected = kFALSE;
1029 Double_t zv=vtx->GetZv();
1032 const AliExternalTrackParam *ct=track->GetTPCInnerParam();
1033 if (!ct) selected = kFALSE;
1037 Float_t rv = TMath::Sqrt(xyz[0] * xyz[0] + xyz[1] * xyz[1]);
1038 if (rv > 3.0) selected = kFALSE;
1039 if (TMath::Abs(xyz[2] - zv)>2.5) selected = kFALSE;
1044 Int_t AliAnalysisTaskEfficiency::GetIndexAC(AliESDtrack *track){
1045 if(!track)return -1;
1047 // Crude selection for A and C side
1049 if (track->Eta() > 0) {
1050 if (track->Charge() > 0)
1056 if (track->Charge() > 0)
1067 Int_t AliAnalysisTaskEfficiency::ConvertePDG(Int_t pdg)
1070 // function converts the pdg values for nuclei (d, t, 3He,
1071 // 4He and their antiparticles) into indices for histo-array
1072 // filled in UserExec
1074 if ( pdg == 1000010020 )return 0;
1075 else if ( pdg == -1000010020 )return 1;
1076 else if ( pdg == 1000010030 )return 2;
1077 else if ( pdg == -1000010030 )return 3;
1078 else if ( pdg == 1000020030 )return 4;
1079 else if ( pdg == -1000020030 )return 5;
1080 else if ( pdg == 1000020040 )return 6;
1081 else if ( pdg == -1000020040 )return 7;