1 //-------------------------------------------------------------------------
3 // This is the PROOF-enabled version of ITS/AliITSComparisonV2.C macro.
4 // Origin: Andrei.Zalite@cern.ch
6 //-------------------------------------------------------------------------
12 #include "TClonesArray.h"
20 #include "AliVEvent.h"
21 #include "AliESDEvent.h"
22 #include "AliMCEvent.h"
23 #include "AliESDtrack.h"
24 #include "AliTrackReference.h"
25 #include "AliMCComparisonTrack.h"
27 #include "AliITSComparisonTask.h"
29 ClassImp(AliITSComparisonTask)
31 extern TStyle *gStyle;
33 AliITSComparisonTask::AliITSComparisonTask()
34 : AliAnalysisTaskSE("AliITSComaprisonTask"),
49 // Default constructor
50 AliInfo("Default constructor AliITSComparisonTask");
51 // Define input and output slots here
52 // Input slot #0 works with a TChain
53 DefineInput(0, TChain::Class());
54 // Output slot #1 TList
55 DefineOutput(1, TList::Class());
60 AliITSComparisonTask::AliITSComparisonTask(const char* name)
61 : AliAnalysisTaskSE(name),
77 AliInfo("Constructor AliITSComparisonTask");
78 // Define input and output slots here
79 // Input slot #0 works with a TChain
80 DefineInput(0, TChain::Class());
81 // Output slot #1 TList
82 DefineOutput(1, TList::Class());
87 void AliITSComparisonTask::UserCreateOutputObjects()
91 AliInfo("AliITSComparisonTask::UserCreateOutputObjects");
92 // Create output container
93 fListOfHistos = new TList();
95 fGood = new TH1F("fGood", "Pt for good tracks", 34, 0.2, 7.0);
96 fFound = new TH1F("fFound", "Pt for found tracks", 34, 0.2, 7.0);
97 fFake = new TH1F("fFake", "Pt for fake tracks", 34, 0.2, 7.0);
98 fP = new TH1F("fP", "PHI resolution", 50, -20., 20.);
99 fL = new TH1F("fL", "LAMBDA resolution", 50, -20., 20.);
100 fPt = new TH1F("fPt", "Relative Pt resolution", 30, -10., 10.);
101 fTip = new TH1F("fTip", "Transverse impact parameter", 40, -800., 800.);
102 fE = new TH1F("fE", "dE/dX for pions with 0.4<p<0.5 GeV/c", 50, 0., 100.);
103 fEp = new TH2F("fEp", "dE/dX vs momentum", 50, 0., 2., 50, 0., 400.);
104 fGoodPhi = new TH1F("fGoodPhi", "Phi for good tracks", 90, 0., 2.*TMath::Pi());
105 fFoundPhi = new TH1F("fFoundPhi", "Phi for found tracks", 90, 0., 2.*TMath::Pi());
106 fLip = new TH1F("fLip","Longitudinal impact parameter", 40, -800., 800.);
109 fListOfHistos->Add(fGood);
110 fListOfHistos->Add(fFound);
111 fListOfHistos->Add(fFake);
112 fListOfHistos->Add(fP);
113 fListOfHistos->Add(fL);
114 fListOfHistos->Add(fPt);
115 fListOfHistos->Add(fTip);
116 fListOfHistos->Add(fE);
117 fListOfHistos->Add(fEp);
118 fListOfHistos->Add(fGoodPhi);
119 fListOfHistos->Add(fFoundPhi);
120 fListOfHistos->Add(fLip);
124 void AliITSComparisonTask::UserExec(Option_t *)
127 // Called for each event
130 AliMCEvent* mcEvent = MCEvent();
132 Printf("ERROR: Could not retrieve MC event");
138 TClonesArray dummy("AliMCComparisonTrack",1000), *refs=&dummy;
140 // Loop over all MC tracks and select "good" tracks
141 for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfTracks(); iTracks++) {
142 AliMCParticle* track = (AliMCParticle*)mcEvent->GetTrack(iTracks);
144 Printf("ERROR: Could not receive track %d (mc loop)", iTracks);
149 if (track->Pt() < 0.2) continue;
150 if (track->Pt() > 7.0) continue;
151 if (TMath::Abs(track->Pz()/track->Pt()) > 0.999) continue;
153 Double_t vx = track->Xv(), vy = track->Yv(), vz = track->Zv();
154 if (TMath::Sqrt(vx*vx+vy*vy) > 3.5) continue;
155 if (TMath::Abs(vz) > 50.) continue;
157 // Loop over Track References
158 Bool_t labelTPC = kFALSE;
159 AliTrackReference* trackRef = 0;
160 UInt_t iTSLayerMap = 0;
161 for (Int_t iTrackRef = 0; iTrackRef < track->GetNumberOfTrackReferences(); iTrackRef++) {
162 trackRef = track->GetTrackReference(iTrackRef);
164 Int_t detectorId = trackRef->DetectorId();
165 if (detectorId == AliTrackReference::kITS) {
166 Float_t radius = trackRef->R();
167 if (radius > 2.5 && radius < 5.5)
169 else if (radius > 5.5 && radius < 8.5)
170 iTSLayerMap |= 1 << 1;
171 else if (radius > 13.5 && radius < 16.5)
172 iTSLayerMap |= 1 << 2;
173 else if (radius > 22. && radius < 26.)
174 iTSLayerMap |= 1 << 3;
175 else if (radius > 36. && radius < 41.)
176 iTSLayerMap |= 1 << 4;
177 else if (radius > 42. && radius < 46.)
178 iTSLayerMap |= 1 << 5;
180 Printf("Wrong radius %f ", radius);
185 if (detectorId == AliTrackReference::kTPC) {
190 } // track references loop
192 // Skip tracks that passed not all ITS layers
193 if (iTSLayerMap != 0x3F) continue;
197 AliMCComparisonTrack* ref = new((*refs)[nt]) AliMCComparisonTrack();
198 ref->SetLabel(iTracks);
199 TParticle* particle = track->Particle();
200 Int_t pdg = particle->GetPdgCode();
202 ref->SetPz(track->Pz());
203 Float_t pt = track->Pt();
206 Float_t phig = track->Phi();
208 fGoodPhi->Fill(phig);
215 AliVEvent* event = InputEvent();
217 Printf("ERROR: Could not retrieve event");
221 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
228 Int_t mcGoods = refs->GetEntriesFast();
230 // Loop over all "good" MC tracks
231 for (Int_t k = 0; k < mcGoods; k++) {
232 AliMCComparisonTrack* ref = (AliMCComparisonTrack*)refs->UncheckedAt(k);
233 Int_t mcLabel = ref->GetLabel();
234 Float_t ptg = ref->GetPt();
235 Float_t phiG = ref->GetPhi();
237 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
238 AliESDtrack* track = esd->GetTrack(iTrack);
240 Printf("ERROR: Could not receive track %d", iTrack);
244 // if (! track->IsOn(AliESDtrack::kTPCrefit)) continue;
245 if (! track->IsOn(AliESDtrack::kITSrefit)) continue;
247 Int_t label = track->GetLabel();
249 if (mcLabel == TMath::Abs(label)) {
250 if (mcLabel == label) {
253 fFoundPhi->Fill(phiG);
261 Float_t phi = track->Phi();
262 Double_t pt = track->Pt();
264 Float_t lam = TMath::ATan2(track->Pz(),pt);
267 Float_t phig = ref->GetPhi();
268 fP->Fill((phi - phig)*1000.);
270 Float_t lamg = TMath::ATan2(ref->GetPz(),ptg);
271 fL->Fill((lam - lamg)*1000.);
273 fPt->Fill((pt1 - 1/ptg)/(1/ptg)*100.);
275 Float_t d,z; track->GetImpactParameters(d,z);
276 fTip->Fill(10000.*d);
277 fLip->Fill(10000.*z);
279 Float_t mom = pt/TMath::Cos(lam);
280 Float_t dedx = track->GetTPCsignal();
281 fEp->Fill(mom, dedx, 1.);
283 Int_t pdg = ref->GetPDG();
284 if (TMath::Abs(pdg )== 211) { // pions
285 if (mom > 0.4 && mom < 0.5) {
298 Printf(" Results: " );
299 Printf(" Found %d Fake %d Lost %d ", nfound, nfake, nlost);
304 PostData(1, fListOfHistos);
307 //________________________________________________________________________
308 void AliITSComparisonTask::Terminate(Option_t *)
311 // Draw result to the screen
312 // Called once at the end of the query
313 fListOfHistos = dynamic_cast<TList*>(GetOutputData(1));
314 if (!fListOfHistos) {
315 Printf("ERROR: fListOfHistos not available");
319 fGood = dynamic_cast<TH1F*>(fListOfHistos->At(0));
320 fFound = dynamic_cast<TH1F*>(fListOfHistos->At(1));
321 fFake = dynamic_cast<TH1F*>(fListOfHistos->At(2));
322 fP = dynamic_cast<TH1F*>(fListOfHistos->At(3));
323 fL = dynamic_cast<TH1F*>(fListOfHistos->At(4));
324 fPt = dynamic_cast<TH1F*>(fListOfHistos->At(5));
325 fTip = dynamic_cast<TH1F*>(fListOfHistos->At(6));
326 fE = dynamic_cast<TH1F*>(fListOfHistos->At(7));
327 fEp = dynamic_cast<TH2F*>(fListOfHistos->At(8));
328 fGoodPhi = dynamic_cast<TH1F*>(fListOfHistos->At(9));
329 fFoundPhi = dynamic_cast<TH1F*>(fListOfHistos->At(10));
330 fLip = dynamic_cast<TH1F*>(fListOfHistos->At(11));
332 gStyle->SetOptStat(111110);
333 gStyle->SetOptFit(1);
335 TCanvas* c1 = new TCanvas("c1", "", 0, 0, 700, 850);
339 TPad* p1 = new TPad("p1", "", 0., 0.3, 0.5, 0.6); p1->Draw();
340 p1->cd(); p1->SetFillColor(42); p1->SetFrameFillColor(10);
341 fP->SetFillColor(4); fP->SetXTitle("(mrad)");
342 if (fP->GetEntries() < minc) fP->Draw(); else fP->Fit("gaus"); c1->cd();
344 TPad* p2 = new TPad("p2", "", 0.5, 0.3, 1., 0.6); p2->Draw();
345 p2->cd(); p2->SetFillColor(42); p2->SetFrameFillColor(10);
346 fL->SetFillColor(4); fL->SetXTitle("(mrad)");
347 if (fL->GetEntries() < minc) fL->Draw(); else fL->Fit("gaus"); c1->cd();
349 TPad* p3 = new TPad("p3", "", 0., 0., 0.5, 0.3); p3->Draw();
350 p3->cd(); p3->SetFillColor(42); p3->SetFrameFillColor(10);
351 fPt->SetFillColor(2); fPt->SetXTitle("(%)");
352 if (fPt->GetEntries() < minc) fPt->Draw(); else fPt->Fit("gaus"); c1->cd();
354 TPad* p4 = new TPad("p4", "", 0.5, 0., 1., 0.3); p4->Draw();
355 p4->cd(); p4->SetFillColor(42); p4->SetFrameFillColor(10);
356 fTip->SetFillColor(6); fTip->SetXTitle("(micron)");
357 if (fTip->GetEntries() < minc) fTip->Draw(); else fTip->Fit("gaus");
358 fLip->Draw("same"); c1->cd();
360 TPad* p5 = new TPad("p5", "", 0., 0.6, 1., 1.); p5->Draw(); p5->cd();
361 p5->SetFillColor(41); p5->SetFrameFillColor(10);
362 fFound->Sumw2(); fGood->Sumw2(); fFake->Sumw2();
363 TH1F* hg = new TH1F("hg", "Efficiency for good tracks", 34, 0.2, 7.0);
364 TH1F* hf = new TH1F("hf", "Efficiency for fake tracks", 34, 0.2, 7.0);
365 hg->Divide(fFound,fGood,1.,1.,"B");
366 hf->Divide(fFake,fGood,1.,1.,"B");
367 hg->SetLineColor(4); hg->SetLineWidth(2);
369 hg->SetYTitle("Tracking efficiency");
370 hg->SetXTitle("Pt (GeV/c)");
373 TLine* line1 = new TLine(0.2, 1.0, 7.0, 1.0); line1->SetLineStyle(4);
375 TLine* line2 = new TLine(0.2, 0.9, 7.0, 0.9); line2->SetLineStyle(4);
379 hf->SetFillStyle(3013);
382 hf->Draw("histsame");
383 TText* text = new TText(0.461176, 0.248448, "Fake tracks");
384 text->SetTextSize(0.05);
386 text = new TText(0.453919, 1.11408, "Good tracks");
387 text->SetTextSize(0.05);
390 TCanvas* c2 = new TCanvas("c2", "", 320, 32, 530, 590);
392 TPad* p6 = new TPad("p6", "", 0., 0., 1., .5); p6->Draw();
393 p6->cd(); p6->SetFillColor(42); p6->SetFrameFillColor(10);
394 fE->SetFillColor(2); fE->SetFillStyle(3005);
395 fE->SetXTitle("Arbitrary Units");
396 if (fE->GetEntries() < minc) fE->Draw(); else fE->Fit("gaus"); c2->cd();
398 TPad* p7 = new TPad("p7", "", 0., 0.5, 1., 1.); p7->Draw();
399 p7->cd(); p7->SetFillColor(42); p7->SetFrameFillColor(10);
400 fEp->SetFillColor(2); fEp->SetFillStyle(3005);
401 fEp->SetMarkerStyle(8);
402 fEp->SetMarkerSize(0.4);
403 fEp->SetXTitle("p (GeV/c)"); fEp->SetYTitle("dE/dX (Arb. Units)");
404 fEp->Draw(); c1->cd();
406 TCanvas* c3 = new TCanvas("c3", "", 10, 10, 510, 510);
407 c3->SetFillColor(42); c3->SetFrameFillColor(10);
408 fFoundPhi->Sumw2(); fGoodPhi->Sumw2();
409 TH1F* hgphi = new TH1F("hgphi", "Efficiency for good tracks (Phi)",
410 90, 0., 2.*TMath::Pi());
411 hgphi->Divide(fFoundPhi, fGoodPhi, 1., 1., "B");
412 hgphi->SetYTitle("Tracking efficiency");
413 hgphi->SetXTitle("Phi (rad)");
416 TFile fc("AliITSComparison.root","RECREATE");