]>
Commit | Line | Data |
---|---|---|
1 | /************************************************************************** | |
2 | * Copyright(c) 1998-1999, 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 | /* $Id: AliTRDefficiencyMC.cxx 27496 2008-07-22 08:35:45Z cblume $ */ | |
17 | ||
18 | //////////////////////////////////////////////////////////////////////////// | |
19 | // // | |
20 | // Reconstruction QA // | |
21 | // // | |
22 | // Authors: // | |
23 | // Markus Fasel <M.Fasel@gsi.de> // | |
24 | // // | |
25 | //////////////////////////////////////////////////////////////////////////// | |
26 | ||
27 | #include <TObjArray.h> | |
28 | #include <TClonesArray.h> | |
29 | #include <TPad.h> | |
30 | #include <TLegend.h> | |
31 | #include <TProfile.h> | |
32 | #include <TMath.h> | |
33 | #include <TDatabasePDG.h> | |
34 | #include "TTreeStream.h" | |
35 | ||
36 | #include "AliMagF.h" | |
37 | #include "AliPID.h" | |
38 | #include "AliESDtrack.h" | |
39 | #include "AliMathBase.h" | |
40 | #include "AliTrackReference.h" | |
41 | #include "AliAnalysisManager.h" | |
42 | ||
43 | #include "AliTRDcluster.h" | |
44 | #include "AliTRDseedV1.h" | |
45 | #include "AliTRDtrackV1.h" | |
46 | #include "Cal/AliTRDCalPID.h" | |
47 | #include "info/AliTRDtrackInfo.h" | |
48 | #include "AliTRDinfoGen.h" | |
49 | #include "AliTRDefficiencyMC.h" | |
50 | ||
51 | ClassImp(AliTRDefficiencyMC) | |
52 | Float_t AliTRDefficiencyMC::fgPCut = 0.2; //[GeV/c] | |
53 | Float_t AliTRDefficiencyMC::fgPhiCut = 50.; //[deg] | |
54 | Float_t AliTRDefficiencyMC::fgThtCut = 50.; //[deg] | |
55 | //_____________________________________________________________________________ | |
56 | AliTRDefficiencyMC::AliTRDefficiencyMC() | |
57 | :AliTRDrecoTask() | |
58 | { | |
59 | // | |
60 | // Default constructor | |
61 | // | |
62 | } | |
63 | ||
64 | AliTRDefficiencyMC::AliTRDefficiencyMC(char* name) | |
65 | :AliTRDrecoTask(name, "Combined Tracking Efficiency") | |
66 | { | |
67 | // | |
68 | // Default constructor | |
69 | // | |
70 | } | |
71 | ||
72 | //_____________________________________________________________________________ | |
73 | void AliTRDefficiencyMC::UserCreateOutputObjects(){ | |
74 | // | |
75 | // Create output objects | |
76 | // | |
77 | ||
78 | fContainer = Histos(); | |
79 | } | |
80 | ||
81 | //_____________________________________________________________________________ | |
82 | void AliTRDefficiencyMC::UserExec(Option_t *){ | |
83 | // | |
84 | // Execute the task: | |
85 | // | |
86 | // Loop over TrackInfos | |
87 | // 1st: check if there is a trackTRD | |
88 | // 2nd: put conditions on the track: | |
89 | // - check if we did not register it before | |
90 | // - check if we have Track References for the track | |
91 | // 3rd: Register track: | |
92 | // - accepted if both conditions are fulfilled | |
93 | // - contamination if at least one condition is not fulfilled | |
94 | // 4th: check Monte-Carlo Track wheter findable or not if there is no TRD track in this track info | |
95 | // 5th: register MC track as rejected if findable and not jet registered | |
96 | // Match the registers accepted and rejected and clean register rejected | |
97 | // Fill the histograms | |
98 | // | |
99 | const Int_t kArraySize = 10000; // Store indices of track references in an array | |
100 | Int_t indexAccept[kArraySize], | |
101 | indexReject[kArraySize], | |
102 | indexContam[kArraySize]; | |
103 | memset(indexAccept, 0, sizeof(Int_t) * kArraySize); | |
104 | memset(indexReject, 0, sizeof(Int_t) * kArraySize); | |
105 | memset(indexContam, 0, sizeof(Int_t) * kArraySize); | |
106 | Int_t naccept(0), | |
107 | nreject(0), | |
108 | nfindnt(0), | |
109 | nkink(0), | |
110 | ncontam(0); | |
111 | Bool_t isContamination = kFALSE; | |
112 | ||
113 | fTracks = dynamic_cast<TObjArray *>(GetInputData(1)); | |
114 | if(!fTracks) return; | |
115 | Int_t nTrackInfos = fTracks->GetEntriesFast(); | |
116 | AliDebug(2, Form(" CANDIDATE TRACKS[%d]", nTrackInfos)); | |
117 | ||
118 | AliTRDtrackV1 *trackTRD(NULL); | |
119 | AliTRDtrackInfo *trkInf(NULL); | |
120 | for(Int_t itinf = 0; itinf < nTrackInfos; itinf++){ | |
121 | trkInf = dynamic_cast<AliTRDtrackInfo *>(fTracks->UncheckedAt(itinf)); | |
122 | if(!trkInf) continue; | |
123 | ||
124 | if(trkInf->GetTrack() || trkInf->GetNumberOfClustersRefit()){ | |
125 | isContamination = (IsRegistered(trkInf,indexAccept,naccept)>=0); | |
126 | if(!trkInf->GetNTrackRefs()){ | |
127 | // We reject the track since the Monte Carlo Information is missing | |
128 | AliDebug(2, Form("MC(Track Reference) missing @ label[%d]", trkInf->GetLabel())); | |
129 | isContamination = kTRUE; | |
130 | // Debugging | |
131 | if(trackTRD && DebugLevel()>5) FillStreamTrackWOMC(trkInf); | |
132 | } | |
133 | if(isContamination){ | |
134 | // reject kink (we count these only once) | |
135 | if(trkInf->GetKinkIndex()){ | |
136 | AliDebug(4, Form(" track @ idx[%d] MC[%d] is kink.", itinf, trkInf->GetLabel())); | |
137 | nkink++; | |
138 | continue; | |
139 | } | |
140 | // Register track as contamination | |
141 | AliDebug(4, Form(" track @ idx[%d] MC[%d] is contamination.", itinf, trkInf->GetLabel())); | |
142 | indexContam[ncontam++]=itinf; | |
143 | continue; | |
144 | } | |
145 | // Accept track | |
146 | AliDebug(4, Form(" track @ idx[%d] is ACCEPTED.", itinf)); | |
147 | // Register track as accepted | |
148 | indexAccept[naccept++] = itinf; | |
149 | }else{ | |
150 | Int_t code(0); | |
151 | if((code=IsFindableNot(trkInf))){ | |
152 | AliDebug(4, Form(" track @ idx[%d] MC[%d] not findable [%d].", itinf, trkInf->GetLabel(), code)); | |
153 | nfindnt++; | |
154 | } else { | |
155 | // register track as rejected if not already registered there | |
156 | // Attention: | |
157 | // These track infos are not!!! registered as contamination | |
158 | if(IsRegistered(trkInf, indexReject, nreject)<0){ | |
159 | AliDebug(4, Form(" track @ idx[%d] MC[%d] is missed.", itinf, trkInf->GetLabel())); | |
160 | indexReject[nreject++] = itinf; | |
161 | } | |
162 | } | |
163 | } | |
164 | } | |
165 | AliDebug(2, Form("TRACKS STATISTICS naccept[%d] ncontam[%d] nkink[%d] nmissed[%d] nfindnt[%d] ALL[%d] LOST[%d]", naccept, ncontam, nkink, nreject, nfindnt, naccept+ncontam+nkink, nTrackInfos-(naccept+nreject+ncontam+nkink+nfindnt))); | |
166 | ||
167 | // we have to check if the rejected tracks are registered as found | |
168 | // a possible source for this: | |
169 | // first the track is not found by the barrel tracking but it is later found | |
170 | // by the stand alone tracking, then two track info objects with the same | |
171 | // label would be created | |
172 | // Attention: | |
173 | // these tracks are not! registered as contamination | |
174 | Int_t tmprejected[kArraySize]; Int_t nrej = nreject; | |
175 | memcpy(tmprejected, indexReject, sizeof(Int_t) * nreject); | |
176 | nreject = 0; | |
177 | for(Int_t irej = 0; irej < nrej; irej++){ | |
178 | trkInf = dynamic_cast<AliTRDtrackInfo *>(fTracks->At(tmprejected[irej])); | |
179 | Int_t idx(-1); | |
180 | if((idx = IsRegistered(trkInf,indexAccept,naccept))<0){ | |
181 | indexReject[nreject++] = tmprejected[irej]; | |
182 | }else{ | |
183 | //printf("tracks @ accept[%d] missed[%d] are the same.\n", indexAccept[idx], tmprejected[irej]); | |
184 | } | |
185 | } | |
186 | ||
187 | // Fill Histograms | |
188 | FillHistograms(naccept, &indexAccept[0], kAccept); | |
189 | FillHistograms(nreject, &indexReject[0], kMiss); | |
190 | FillHistograms(ncontam, &indexContam[0], kFake); | |
191 | ||
192 | Int_t nall(naccept + nreject); | |
193 | AliInfo(Form("%3d Tracks: MC[%3d] TRD[%3d | %5.2f%%] Fake[%3d | %5.2f%%]", | |
194 | (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), | |
195 | nall, | |
196 | naccept, | |
197 | (nall ? 1.E2*Float_t(naccept)/Float_t(nall) : 0.), | |
198 | ncontam, | |
199 | (nall ? 1.E2*Float_t(ncontam)/Float_t(nall) : 0.))); | |
200 | ||
201 | PostData(1, fContainer); | |
202 | } | |
203 | ||
204 | ||
205 | //_____________________________________________________________________________ | |
206 | Bool_t AliTRDefficiencyMC::PostProcess() | |
207 | { | |
208 | // | |
209 | // Post Process | |
210 | // | |
211 | // Change the histogram style | |
212 | // For species histograms apply the colors connected with the given particle species | |
213 | fNRefFigures = 8; | |
214 | return kTRUE; | |
215 | } | |
216 | ||
217 | //_____________________________________________________________________________ | |
218 | Bool_t AliTRDefficiencyMC::GetRefFigure(Int_t ifig){ | |
219 | // | |
220 | // Plot the histograms | |
221 | // | |
222 | if(ifig >= fNRefFigures) return kFALSE; | |
223 | if(!gPad) return kFALSE; | |
224 | gPad->SetLogx(kTRUE); | |
225 | if(ifig < 2){ | |
226 | (dynamic_cast<TH1 *>(fContainer->At(ifig)))->Draw("e1"); | |
227 | return kTRUE; | |
228 | } | |
229 | TH1 *h(NULL); | |
230 | TLegend *leg=new TLegend(.65, .12, .85, .3); | |
231 | leg->SetHeader("Charge"); | |
232 | leg->SetBorderSize(1);leg->SetFillColor(kWhite); | |
233 | switch(ifig){ | |
234 | case 2: | |
235 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram)); | |
236 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
237 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+1)); | |
238 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
239 | leg->Draw(); | |
240 | break; | |
241 | case 3: | |
242 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+2)); | |
243 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
244 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+3)); | |
245 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
246 | leg->Draw(); | |
247 | break; | |
248 | case 4: | |
249 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+4)); | |
250 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
251 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+5)); | |
252 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
253 | leg->Draw(); | |
254 | break; | |
255 | case 5: | |
256 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+6)); | |
257 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
258 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+7)); | |
259 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
260 | leg->Draw(); | |
261 | break; | |
262 | case 6: | |
263 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+8)); | |
264 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
265 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+9)); | |
266 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
267 | leg->Draw(); | |
268 | break; | |
269 | case 7: | |
270 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+10)); | |
271 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
272 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+11)); | |
273 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
274 | leg->Draw(); | |
275 | break; | |
276 | case 8: | |
277 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+12)); | |
278 | h->Draw("e1"); leg->AddEntry(h, " -", "pl"); | |
279 | h=dynamic_cast<TH1 *>(fContainer->At(kEfficiencySpeciesHistogram+13)); | |
280 | h->Draw("e1same"); leg->AddEntry(h, " +", "pl"); | |
281 | leg->Draw(); | |
282 | break; | |
283 | } | |
284 | return kTRUE; | |
285 | } | |
286 | ||
287 | //_____________________________________________________________________________ | |
288 | TObjArray *AliTRDefficiencyMC::Histos(){ | |
289 | // | |
290 | // Create the histograms | |
291 | // | |
292 | ||
293 | if(fContainer) return fContainer; | |
294 | const Int_t nbins = AliTRDCalPID::kNMom; | |
295 | Float_t xbins[nbins+1] = {fgPCut, .7, .9, 1.3, 1.7, 2.4, 3.5, 4.5, 5.5, 7., 9., 11.}; | |
296 | const Int_t marker[2][AliPID::kSPECIES+1] = { | |
297 | {20, 21, 22, 23, 29, 2}, | |
298 | {24, 25, 26, 27, 30, 5} | |
299 | }; | |
300 | ||
301 | fContainer = new TObjArray();fContainer->Expand(14); | |
302 | ||
303 | TH1 *h(NULL); | |
304 | fContainer->AddAt(h=new TProfile("hEff", "Tracking Efficiency ALL", nbins, xbins), kEfficiencyHistogram); | |
305 | h->SetMarkerStyle(22); | |
306 | h->SetMarkerColor(kBlue); | |
307 | h->GetXaxis()->SetTitle("p [GeV/c]"); | |
308 | h->GetXaxis()->SetMoreLogLabels(); | |
309 | h->GetYaxis()->SetTitle("Efficiency"); | |
310 | h->GetYaxis()->SetRangeUser(0.2, 1.1); | |
311 | fContainer->AddAt(h=new TProfile("hFake", "Fake Tracks", nbins, xbins), kContaminationHistogram); | |
312 | h->SetMarkerStyle(22); | |
313 | h->SetMarkerColor(kBlue); | |
314 | h->GetXaxis()->SetTitle("p [GeV/c]"); | |
315 | h->GetXaxis()->SetMoreLogLabels(); | |
316 | h->GetYaxis()->SetTitle("Contamination"); | |
317 | ||
318 | Char_t sign[]={'+', '-'}; | |
319 | for(Int_t isign = 0; isign < 2; isign++){ | |
320 | for(Int_t ispec = 0; ispec < AliPID::kSPECIES; ispec++){ | |
321 | fContainer->AddAt(h=new TProfile( | |
322 | Form("hEff_%s%c", AliPID::ParticleShortName(ispec), sign[isign]), | |
323 | Form("Tracking Efficiency for %s", AliPID::ParticleName(ispec)), nbins, xbins), | |
324 | kEfficiencySpeciesHistogram+ispec*2+isign); | |
325 | h->SetMarkerStyle(marker[isign][ispec]); | |
326 | h->SetLineColor(AliTRDCalPID::GetPartColor(ispec)); | |
327 | h->SetMarkerColor(kBlack); | |
328 | h->GetXaxis()->SetTitle("p [GeV/c]"); | |
329 | h->GetXaxis()->SetMoreLogLabels(); | |
330 | h->GetYaxis()->SetTitle("Efficiency"); | |
331 | h->GetYaxis()->SetRangeUser(0.2, 1.1); | |
332 | } | |
333 | ||
334 | fContainer->AddAt(h=new TProfile(Form("hEff_PID%c", sign[isign]), "Tracking Efficiency no PID", nbins, xbins), kEfficiencySpeciesHistogram+AliPID::kSPECIES*2+isign); | |
335 | h->SetMarkerStyle(marker[isign][AliPID::kSPECIES]); | |
336 | h->SetMarkerColor(kBlack);h->SetLineColor(kBlack); | |
337 | h->GetXaxis()->SetTitle("p [GeV/c]"); | |
338 | h->GetXaxis()->SetMoreLogLabels(); | |
339 | h->GetYaxis()->SetTitle("Efficiency"); | |
340 | h->GetYaxis()->SetRangeUser(0.2, 1.1); | |
341 | } | |
342 | return fContainer; | |
343 | } | |
344 | ||
345 | //_____________________________________________________________________________ | |
346 | Int_t AliTRDefficiencyMC::IsFindableNot(AliTRDtrackInfo * const trkInf){ | |
347 | // | |
348 | // Apply Cuts on the Monte Carlo track references | |
349 | // return whether track is findable or not | |
350 | // | |
351 | ||
352 | ||
353 | const Float_t chmbHght = AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght(); | |
354 | const Float_t eps(1.E-3); | |
355 | Int_t ntr(trkInf->GetNTrackRefs()); | |
356 | ||
357 | AliDebug(10, Form(" CANDIDATE TrackRefs[%d]", ntr)); | |
358 | // Check if track is findable | |
359 | Double_t mom(0.), phi(0.), tht(0.); | |
360 | Float_t xmin = 10000.0, xmax = 0.0; | |
361 | Float_t ymin = 0.0, ymax = 0.0; | |
362 | Float_t zmin = 0.0, zmax = 0.0; | |
363 | Float_t lastx = 0.0, x = 0.0; | |
364 | Int_t nLayers(0), ntrTRD(0); | |
365 | Int_t sector[20]; | |
366 | AliTrackReference *trackRef(NULL); | |
367 | for(Int_t itr(0); itr<ntr; itr++){ | |
368 | if(!(trackRef = trkInf->GetTrackRef(itr))) continue; | |
369 | x = trackRef->LocalX(); | |
370 | // Be Sure that we are inside TRD | |
371 | if(x < AliTRDinfoGen::GetEndTPC() || x > AliTRDinfoGen::GetEndTRD()) continue; | |
372 | sector[ntrTRD] = Int_t(trackRef->Alpha()/AliTRDgeometry::GetAlpha()); | |
373 | AliDebug(10, Form(" [%2d] x[%7.2f] y[%7.2f] z[%7.2f] Sec[%2d]", itr, trackRef->LocalX(), trackRef->LocalY(), trackRef->Z(), sector[ntrTRD])); | |
374 | if(x < xmin){ | |
375 | xmin = trackRef->LocalX(); | |
376 | ymin = trackRef->LocalY(); | |
377 | zmin = trackRef->Z(); | |
378 | mom = trackRef->P(); | |
379 | } else if(x > xmax){ | |
380 | xmax = trackRef->LocalX(); | |
381 | ymax = trackRef->LocalY(); | |
382 | zmax = trackRef->Z(); | |
383 | } | |
384 | if(itr > 0){ | |
385 | Float_t dist = TMath::Abs(x - lastx); | |
386 | if(TMath::Abs(dist - chmbHght) < eps && sector[ntrTRD]==sector[0]){ | |
387 | AliDebug(10, Form(" dx = %7.2f", dist)); | |
388 | nLayers++; | |
389 | } | |
390 | } | |
391 | lastx = x; | |
392 | ntrTRD++; if(ntrTRD>=20) break; | |
393 | } | |
394 | Double_t dx(xmax - xmin); | |
395 | if(TMath::Abs(dx)<eps) return kNoChmb; | |
396 | ||
397 | phi = (ymax -ymin)/dx; | |
398 | tht = (zmax -zmin)/dx; | |
399 | phi=TMath::ATan(phi)*TMath::RadToDeg(); | |
400 | tht=TMath::ATan(tht)*TMath::RadToDeg(); | |
401 | Bool_t primary = trkInf->IsPrimary(); | |
402 | const AliTRDtrackInfo::AliESDinfo *esd(trkInf->GetESDinfo()); | |
403 | AliDebug(10, Form(" p=%6.3f[GeV/c] phi=%6.2f[deg] theta=%6.2f[deg] nLy[%d]", | |
404 | mom, phi, tht, nLayers)); | |
405 | if(DebugLevel()){ | |
406 | (*DebugStream()) << "IsFindable" | |
407 | << "P=" << mom | |
408 | << "Phi=" << phi | |
409 | << "Tht=" << tht | |
410 | << "Ntr=" << ntrTRD | |
411 | << "NLy=" << nLayers | |
412 | << "Primary=" << primary | |
413 | << "\n"; | |
414 | } | |
415 | ||
416 | // Apply cuts | |
417 | if(!nLayers) return kNoChmb; | |
418 | if(xmax < xmin) return kCurved; | |
419 | if(mom < fgPCut) return kPCut; | |
420 | ||
421 | ||
422 | if(TMath::Abs(phi) > fgPhiCut) return kPhiCut; | |
423 | if(TMath::Abs(tht) > fgThtCut) return kThtCut; | |
424 | ||
425 | if(nLayers < 4){ | |
426 | if(!esd)return kLayer; | |
427 | if(!(esd->GetStatus() & AliESDtrack::kTPCout)) return kLayer; | |
428 | } | |
429 | ||
430 | //if(!trkInf->IsPrimary()) {failCode=kPrimary; return kFALSE;} | |
431 | ||
432 | return kFindable; | |
433 | } | |
434 | ||
435 | //_____________________________________________________________________________ | |
436 | void AliTRDefficiencyMC::FillHistograms(Int_t nTracks, Int_t *indices, ETRDefficiencyMCstatus mode){ | |
437 | // | |
438 | // Fill Histograms in three different modes: | |
439 | // 1st tracks which are found and accepted | |
440 | // 2nd tracks which are not found and not already counted | |
441 | // 3rd contaminating tracks: either double counts (kinks!) or tracks with no MC hit inside TRD | |
442 | // | |
443 | ||
444 | TDatabasePDG *dbPDG(TDatabasePDG::Instance()); | |
445 | Double_t trkmom(0.); // the track momentum | |
446 | Int_t trkpdg(-1); // particle PDG code | |
447 | AliTRDtrackInfo *trkInf(NULL); | |
448 | for(Int_t itk = 0; itk < nTracks; itk++){ | |
449 | trkInf = dynamic_cast<AliTRDtrackInfo *>(fTracks->At(indices[itk])); | |
450 | if(trkInf->GetNTrackRefs()){ | |
451 | // use Monte-Carlo Information for Momentum and PID | |
452 | trkmom = trkInf->GetTrackRef(0)->P(); | |
453 | trkpdg = trkInf->GetPDG(); | |
454 | }else{ | |
455 | // Use TPC Momentum | |
456 | trkmom = trkInf->GetTrack()->P(); | |
457 | } | |
458 | ||
459 | const Char_t *cmode(NULL); | |
460 | switch(mode){ | |
461 | case kAccept: | |
462 | (dynamic_cast<TProfile *>(fContainer->At(kEfficiencyHistogram)))->Fill(trkmom, 1); | |
463 | (dynamic_cast<TProfile *>(fContainer->At(kContaminationHistogram)))->Fill(trkmom, 0); | |
464 | cmode="ACCEPT"; | |
465 | break; | |
466 | case kMiss: | |
467 | (dynamic_cast<TProfile *>(fContainer->At(kEfficiencyHistogram)))->Fill(trkmom, 0); | |
468 | (dynamic_cast<TProfile *>(fContainer->At(kContaminationHistogram)))->Fill(trkmom, 0); | |
469 | cmode="MISS"; | |
470 | break; | |
471 | case kFake: | |
472 | (dynamic_cast<TProfile *>(fContainer->At(kContaminationHistogram)))->Fill(trkmom, 1); | |
473 | cmode="FAKE"; | |
474 | break; | |
475 | } | |
476 | AliDebug(3, Form(" track[%d] MC[%d] Mode[%s]", indices[itk], trkInf->GetLabel(), cmode)); | |
477 | ||
478 | // Fill species histogram | |
479 | Int_t idxSpec = AliTRDpidUtil::Pdg2Pid(TMath::Abs(trkpdg)); | |
480 | Int_t sign = dbPDG->GetParticle(trkpdg)->Charge() > 0. ? 1 : 0; | |
481 | //printf("[%d]%s pdg[%d] sign[%d]\n", idxSpec, AliPID::ParticleName(idxSpec), trkpdg, sign); | |
482 | if(idxSpec < 0) idxSpec = AliPID::kSPECIES; | |
483 | (dynamic_cast<TProfile *>(fContainer->At(kEfficiencySpeciesHistogram + idxSpec*2+sign)))->Fill(trkmom, mode==kAccept?1:0); | |
484 | } | |
485 | } | |
486 | ||
487 | //_____________________________________________________________________________ | |
488 | void AliTRDefficiencyMC::FillStreamTrackWOMC(AliTRDtrackInfo * const trkInf){ | |
489 | // fill debug stream | |
490 | // we want to know: | |
491 | // 1. The event number | |
492 | // 2. The track label | |
493 | // 3. The TRD track label | |
494 | // 4. The frequency of the TRD Label | |
495 | // 5. Momentum from TPC (NO STAND ALONE TRACK) | |
496 | // 6. TPC Phi angle | |
497 | // 7. the TRD track | |
498 | // 8. Monte Carlo PID | |
499 | // 9. We check the Labels of the TRD track according to them we search the maching Monte-Carlo track. | |
500 | // From the matching Monte-Carlo track we store trackRefs, phi and momentum | |
501 | // 10. We may also want to keep the kink index | |
502 | Double_t mom = trkInf->GetESDinfo()->GetOuterParam()->P(); | |
503 | Int_t event = (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(); | |
504 | Int_t label = trkInf->GetLabel(); | |
505 | Int_t kinkIndex = trkInf->GetKinkIndex(); | |
506 | Int_t pdg = trkInf->GetPDG(); | |
507 | Double_t phiTPC = trkInf->GetESDinfo()->GetOuterParam()->Phi(); | |
508 | Int_t labelsTRD[180]; // Container for the cluster labels | |
509 | Int_t sortlabels[360]; // Cluster Labels sorted according their occurancy | |
510 | AliTRDseedV1 *tracklet(NULL); | |
511 | AliTRDcluster *c(NULL); | |
512 | Int_t nclusters(0); | |
513 | AliTRDtrackV1 *trackTRD = trkInf->GetTrack(); | |
514 | for(Int_t il = 0; il < AliTRDgeometry::kNlayer; il++){ | |
515 | tracklet = trackTRD->GetTracklet(il); | |
516 | if(!tracklet) continue; | |
517 | tracklet->ResetClusterIter(); | |
518 | c = NULL; | |
519 | while((c = tracklet->NextCluster())) labelsTRD[nclusters++] = c->GetLabel(0); | |
520 | } | |
521 | // Determine Label and Frequency | |
522 | AliMathBase::Freq(nclusters, const_cast<const Int_t *>(&labelsTRD[0]), &sortlabels[0], kTRUE); | |
523 | Int_t labelTRD = sortlabels[0]; | |
524 | Int_t freqTRD = sortlabels[1]; | |
525 | // find the track info object matching to the TRD track | |
526 | AliTRDtrackInfo *realtrack = 0; | |
527 | TObjArrayIter rtiter(fTracks); | |
528 | while((realtrack = (AliTRDtrackInfo *)rtiter())){ | |
529 | if(realtrack->GetLabel() != labelTRD) continue; | |
530 | break; | |
531 | } | |
532 | TClonesArray trackRefs("AliTrackReference"); | |
533 | Int_t realPdg = -1; | |
534 | Double_t realP = 0.; | |
535 | Double_t realPhi = 0.; | |
536 | if(realtrack){ | |
537 | // pack the track references into the trackRefsContainer | |
538 | for(Int_t iref = 0; iref < realtrack->GetNTrackRefs(); iref++){ | |
539 | new(trackRefs[iref])AliTrackReference(*(realtrack->GetTrackRef(iref))); | |
540 | } | |
541 | realPdg = realtrack->GetPDG(); | |
542 | if(realtrack->GetNTrackRefs()){ | |
543 | realP = realtrack->GetTrackRef(0)->P(); | |
544 | realPhi = realtrack->GetTrackRef(0)->Phi(); | |
545 | } | |
546 | } | |
547 | (*DebugStream()) << "EffMCfake" | |
548 | << "Event=" << event | |
549 | << "Label=" << label | |
550 | << "labelTRD=" << labelTRD | |
551 | << "FreqTRDlabel=" << freqTRD | |
552 | << "TPCp=" << mom | |
553 | << "phiTPC=" << phiTPC | |
554 | << "trackTRD=" << trackTRD | |
555 | << "PDG=" << pdg | |
556 | << "TrackRefs=" << &trackRefs | |
557 | << "RealPDG=" << realPdg | |
558 | << "RealP=" << realP | |
559 | << "RealPhi" << realPhi | |
560 | << "KinkIndex=" << kinkIndex | |
561 | << "\n"; | |
562 | } | |
563 | ||
564 | //_____________________________________________________________________________ | |
565 | Int_t AliTRDefficiencyMC::IsRegistered(AliTRDtrackInfo * const trkInf, Int_t *indices, Int_t nTracks){ | |
566 | // | |
567 | // Checks if track is registered in a given mode | |
568 | // | |
569 | ||
570 | Int_t label(trkInf->GetLabel()); | |
571 | for(Int_t il(nTracks); il--;){ | |
572 | if((dynamic_cast<AliTRDtrackInfo *>(fTracks->At(indices[il])))->GetLabel() == label) return il; | |
573 | } | |
574 | return -1; | |
575 | } | |
576 |