]>
Commit | Line | Data |
---|---|---|
a589a812 | 1 | #include "TROOT.h" |
773d3f8c | 2 | #include "TPDGCode.h" |
a391a274 | 3 | #include "TCanvas.h" |
773d3f8c | 4 | #include "TF1.h" |
5 | #include "TH1F.h" | |
422a2dc0 | 6 | #include "TH1D.h" |
773d3f8c | 7 | #include "TH2F.h" |
c7cf2032 | 8 | #include "TProfile.h" |
422a2dc0 | 9 | #include "TProfile2D.h" |
773d3f8c | 10 | #include "TGraph.h" |
11 | #include "TGraphErrors.h" | |
c46a7947 | 12 | #include "TLegend.h" |
773d3f8c | 13 | |
14 | #include <TClonesArray.h> | |
15 | #include <TObjArray.h> | |
16 | #include <TList.h> | |
17 | ||
773d3f8c | 18 | #include "AliESDEvent.h" |
19 | #include "AliESDInputHandler.h" | |
20 | #include "AliTrackReference.h" | |
21 | ||
c7cf2032 | 22 | #include "AliAnalysisTask.h" |
773d3f8c | 23 | |
c7cf2032 | 24 | #include "AliTRDtrackerV1.h" |
773d3f8c | 25 | #include "AliTRDtrackV1.h" |
c7cf2032 | 26 | #include "AliTRDcluster.h" |
773d3f8c | 27 | #include "AliTRDReconstructor.h" |
28 | #include "AliCDBManager.h" | |
422a2dc0 | 29 | #include "AliTRDpidUtil.h" |
773d3f8c | 30 | |
773d3f8c | 31 | #include "AliTRDpidChecker.h" |
32 | #include "AliTRDtrackInfo/AliTRDtrackInfo.h" | |
33 | ||
34 | // calculate pion efficiency at 90% electron efficiency for 11 momentum bins | |
35 | // this task should be used with simulated data only | |
36 | ||
37 | ClassImp(AliTRDpidChecker) | |
38 | ||
39 | //________________________________________________________________________ | |
3d86166d | 40 | AliTRDpidChecker::AliTRDpidChecker() |
41 | :AliTRDrecoTask("PID", "PID Checker") | |
773d3f8c | 42 | ,fReconstructor(0x0) |
74b2e03d | 43 | ,fUtil(0x0) |
c46a7947 | 44 | ,fGraph(0x0) |
45 | ,fEfficiency(0x0) | |
773d3f8c | 46 | { |
47 | // | |
48 | // Default constructor | |
49 | // | |
50 | ||
51 | fReconstructor = new AliTRDReconstructor(); | |
52 | fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam()); | |
74b2e03d | 53 | |
54 | fUtil = new AliTRDpidUtil(); | |
55 | ||
56 | InitFunctorList(); | |
773d3f8c | 57 | } |
58 | ||
59 | ||
60 | //________________________________________________________________________ | |
61 | AliTRDpidChecker::~AliTRDpidChecker() | |
62 | { | |
c46a7947 | 63 | if(fGraph){fGraph->Delete(); delete fGraph;} |
64 | if(fReconstructor) delete fReconstructor; | |
74b2e03d | 65 | if(fUtil) delete fUtil; |
773d3f8c | 66 | } |
67 | ||
68 | ||
69 | //________________________________________________________________________ | |
70 | void AliTRDpidChecker::CreateOutputObjects() | |
71 | { | |
72 | // Create histograms | |
73 | // Called once | |
74 | ||
75 | OpenFile(0, "RECREATE"); | |
74b2e03d | 76 | fContainer = Histos(); |
77 | } | |
78 | ||
79 | ||
80 | //_______________________________________________________ | |
81 | TObjArray * AliTRDpidChecker::Histos(){ | |
82 | ||
83 | // | |
84 | // Create QA histograms | |
85 | // | |
86 | if(fContainer) return fContainer; | |
87 | ||
422a2dc0 | 88 | Int_t xBins = AliPID::kSPECIES*AliTRDCalPID::kNMom; |
c46a7947 | 89 | fContainer = new TObjArray(); fContainer->Expand(7); |
773d3f8c | 90 | |
422a2dc0 | 91 | const Float_t epsilon = 1./(2*(AliTRDpidUtil::kBins-1)); // get nice histos with bin center at 0 and 1 |
773d3f8c | 92 | |
93 | // histos of the electron probability of all 5 particle species and 11 momenta for the 2-dim LQ method | |
b1957d3c | 94 | fEfficiency = new TObjArray(3); |
95 | fEfficiency->SetOwner(); fEfficiency->SetName("Efficiency"); | |
c46a7947 | 96 | fContainer->AddAt(fEfficiency, kEfficiency); |
a589a812 | 97 | |
98 | TH1 *h = 0x0; | |
99 | if(!(h = (TH2F*)gROOT->FindObject("PID_LQ"))){ | |
100 | h = new TH2F("PID_LQ", "", xBins, -0.5, xBins - 0.5, | |
101 | AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon); | |
102 | } else h->Reset(); | |
103 | fEfficiency->AddAt(h, kLQ); | |
773d3f8c | 104 | |
1a6d7c9a | 105 | // histos of the electron probability of all 5 particle species and 11 momenta for the neural network method |
a589a812 | 106 | if(!(h = (TH2F*)gROOT->FindObject("PID_NN"))){ |
107 | h = new TH2F("PID_NN", "", | |
422a2dc0 | 108 | xBins, -0.5, xBins - 0.5, |
a589a812 | 109 | AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon); |
110 | } else h->Reset(); | |
111 | fEfficiency->AddAt(h, kNN); | |
c46a7947 | 112 | |
113 | // histos of the electron probability of all 5 particle species and 11 momenta for the ESD output | |
a589a812 | 114 | if(!(h = (TH2F*)gROOT->FindObject("PID_ESD"))){ |
115 | h = new TH2F("PID_ESD", "", | |
c46a7947 | 116 | xBins, -0.5, xBins - 0.5, |
a589a812 | 117 | AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon); |
118 | } else h->Reset(); | |
119 | fEfficiency->AddAt(h, kESD); | |
1a6d7c9a | 120 | |
121 | // histos of the dE/dx distribution for all 5 particle species and 11 momenta | |
a589a812 | 122 | if(!(h = (TH2F*)gROOT->FindObject("dEdx"))){ |
123 | h = new TH2F("dEdx", "", | |
422a2dc0 | 124 | xBins, -0.5, xBins - 0.5, |
a589a812 | 125 | 200, 0, 10000); |
126 | } else h->Reset(); | |
127 | fContainer->AddAt(h, kdEdx); | |
773d3f8c | 128 | |
c46a7947 | 129 | // histos of the dE/dx slices for all 5 particle species and 11 momenta |
a589a812 | 130 | if(!(h = (TH2F*)gROOT->FindObject("dEdxSlice"))){ |
131 | h = new TH2F("dEdxSlice", "", | |
b718144c | 132 | xBins*AliTRDReconstructor::kLQslices, -0.5, xBins*AliTRDReconstructor::kLQslices - 0.5, |
a589a812 | 133 | 200, 0, 5000); |
134 | } else h->Reset(); | |
135 | fContainer->AddAt(h, kdEdxSlice); | |
b718144c | 136 | |
c7cf2032 | 137 | // histos of the pulse height distribution for all 5 particle species and 11 momenta |
a589a812 | 138 | if(!(h = (TH2F*)gROOT->FindObject("PH"))){ |
139 | h = new TProfile2D("PH", "", | |
422a2dc0 | 140 | xBins, -0.5, xBins - 0.5, |
a589a812 | 141 | AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5); |
142 | } else h->Reset(); | |
143 | fContainer->AddAt(h, kPH); | |
422a2dc0 | 144 | |
3afc1da3 | 145 | // histos of the number of clusters distribution for all 5 particle species and 11 momenta |
a589a812 | 146 | if(!(h = (TH2F*)gROOT->FindObject("NClus"))){ |
147 | h = new TH2F("NClus", "", | |
3afc1da3 | 148 | xBins, -0.5, xBins - 0.5, |
a589a812 | 149 | AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5); |
150 | } else h->Reset(); | |
151 | fContainer->AddAt(h, kNClus); | |
3afc1da3 | 152 | |
c7cf2032 | 153 | |
154 | // momentum distributions - absolute and in momentum bins | |
a589a812 | 155 | if(!(h = (TH1F*)gROOT->FindObject("hMom"))){ |
156 | h = new TH1F("hMom", "momentum distribution", 100, 0., 12.); | |
157 | } else h->Reset(); | |
158 | fContainer->AddAt(h, kMomentum); | |
159 | ||
160 | if(!(h = (TH1F*)gROOT->FindObject("hMomBin"))){ | |
161 | h = new TH1F("hMomBin", "momentum distribution in momentum bins", AliTRDCalPID::kNMom, 0.5, 11.5); | |
162 | } else h->Reset(); | |
163 | fContainer->AddAt(h, kMomentumBin); | |
c7cf2032 | 164 | |
c4c5bbfb | 165 | |
74b2e03d | 166 | return fContainer; |
167 | } | |
168 | ||
169 | ||
170 | //________________________________________________________________________ | |
171 | Bool_t AliTRDpidChecker::CheckTrackQuality(const AliTRDtrackV1* track) | |
172 | { | |
173 | // | |
174 | // Check if the track is ok for PID | |
175 | // | |
176 | ||
177 | if(track->GetNumberOfTracklets() == AliTRDgeometry::kNlayer) return 1; | |
178 | // if(!fESD) | |
179 | // return 0; | |
180 | ||
181 | return 0; | |
773d3f8c | 182 | } |
183 | ||
74b2e03d | 184 | |
773d3f8c | 185 | //________________________________________________________________________ |
74b2e03d | 186 | Int_t AliTRDpidChecker::CalcPDG(AliTRDtrackV1* track) |
773d3f8c | 187 | { |
773d3f8c | 188 | |
74b2e03d | 189 | track -> SetReconstructor(fReconstructor); |
773d3f8c | 190 | |
74b2e03d | 191 | fReconstructor -> SetOption("nn"); |
192 | track -> CookPID(); | |
773d3f8c | 193 | |
74b2e03d | 194 | if(track -> GetPID(AliPID::kElectron) > track -> GetPID(AliPID::kMuon) + track -> GetPID(AliPID::kPion) + track -> GetPID(AliPID::kKaon) + track -> GetPID(AliPID::kProton)){ |
195 | return kElectron; | |
196 | } | |
197 | else if(track -> GetPID(kProton) > track -> GetPID(AliPID::kPion) && track -> GetPID(AliPID::kProton) > track -> GetPID(AliPID::kKaon) && track -> GetPID(AliPID::kProton) > track -> GetPID(AliPID::kMuon)){ | |
198 | return kProton; | |
199 | } | |
200 | else if(track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kMuon) && track -> GetPID(AliPID::kKaon) > track -> GetPID(AliPID::kPion)){ | |
201 | return kKPlus; | |
202 | } | |
203 | else if(track -> GetPID(AliPID::kMuon) > track -> GetPID(AliPID::kPion)){ | |
204 | return kMuonPlus; | |
205 | } | |
206 | else{ | |
207 | return kPiPlus; | |
208 | } | |
209 | } | |
210 | ||
211 | ||
212 | //_______________________________________________________ | |
213 | TH1 *AliTRDpidChecker::PlotLQ(const AliTRDtrackV1 *track) | |
214 | { | |
215 | // | |
216 | // Plot the probabilities for electrons using 2-dim LQ | |
217 | // | |
218 | ||
219 | if(track) fTrack = track; | |
220 | if(!fTrack){ | |
221 | AliWarning("No Track defined."); | |
222 | return 0x0; | |
223 | } | |
224 | ||
225 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
226 | ||
c46a7947 | 227 | if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){ |
228 | AliWarning("No Histogram defined."); | |
229 | return 0x0; | |
230 | } | |
231 | TH2F *hPIDLQ = 0x0; | |
232 | if(!(hPIDLQ = dynamic_cast<TH2F *>(fEfficiency->At(kLQ)))){ | |
74b2e03d | 233 | AliWarning("No Histogram defined."); |
234 | return 0x0; | |
235 | } | |
236 | ||
237 | AliTRDtrackV1 cTrack(*fTrack); | |
238 | cTrack.SetReconstructor(fReconstructor); | |
239 | ||
240 | Int_t pdg = 0; | |
241 | Float_t momentum = 0.; | |
242 | ||
243 | if(fMC){ | |
244 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
245 | pdg = fMC->GetPDG(); | |
246 | } else{ | |
247 | //AliWarning("No MC info available!"); | |
248 | momentum = cTrack.GetMomentum(0); | |
249 | pdg = CalcPDG(&cTrack); | |
250 | } | |
c46a7947 | 251 | if(momentum < 0.4) return 0x0;; |
252 | if(momentum > 12.) return 0x0;; | |
74b2e03d | 253 | |
254 | fReconstructor -> SetOption("!nn"); | |
255 | cTrack.CookPID(); | |
c46a7947 | 256 | Int_t iMomBin = fUtil->GetMomentumBin(momentum); |
74b2e03d | 257 | |
258 | switch(pdg){ | |
259 | case kElectron: | |
260 | case kPositron: | |
261 | hPIDLQ -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
262 | break; | |
263 | case kMuonPlus: | |
264 | case kMuonMinus: | |
265 | hPIDLQ -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, cTrack .GetPID(AliPID::kElectron)); | |
266 | break; | |
267 | case kPiPlus: | |
268 | case kPiMinus: | |
269 | hPIDLQ -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, cTrack .GetPID(AliPID::kElectron)); | |
270 | break; | |
271 | case kKPlus: | |
272 | case kKMinus: | |
273 | hPIDLQ -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, cTrack .GetPID(AliPID::kElectron)); | |
274 | break; | |
275 | case kProton: | |
276 | case kProtonBar: | |
277 | hPIDLQ -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
278 | break; | |
279 | } | |
280 | ||
281 | return hPIDLQ; | |
282 | } | |
283 | ||
422a2dc0 | 284 | |
74b2e03d | 285 | //_______________________________________________________ |
286 | TH1 *AliTRDpidChecker::PlotNN(const AliTRDtrackV1 *track) | |
287 | { | |
288 | // | |
289 | // Plot the probabilities for electrons using 2-dim LQ | |
290 | // | |
291 | ||
292 | if(track) fTrack = track; | |
293 | if(!fTrack){ | |
294 | AliWarning("No Track defined."); | |
295 | return 0x0; | |
296 | } | |
297 | ||
298 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
422a2dc0 | 299 | |
c46a7947 | 300 | if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){ |
301 | AliWarning("No Histogram defined."); | |
302 | return 0x0; | |
303 | } | |
74b2e03d | 304 | TH2F *hPIDNN; |
48b51864 | 305 | if(!(hPIDNN = dynamic_cast<TH2F *>(fEfficiency->At(kNN)))){ |
74b2e03d | 306 | AliWarning("No Histogram defined."); |
307 | return 0x0; | |
308 | } | |
309 | ||
310 | ||
311 | AliTRDtrackV1 cTrack(*fTrack); | |
312 | cTrack.SetReconstructor(fReconstructor); | |
313 | ||
314 | Int_t pdg = 0; | |
315 | Float_t momentum = 0.; | |
74b2e03d | 316 | if(fMC){ |
317 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
318 | pdg = fMC->GetPDG(); | |
319 | } else { | |
320 | //AliWarning("No MC info available!"); | |
321 | momentum = cTrack.GetMomentum(0); | |
322 | pdg = CalcPDG(&cTrack); | |
323 | } | |
c46a7947 | 324 | if(momentum < 0.4) return 0x0;; |
325 | if(momentum > 12.) return 0x0;; | |
74b2e03d | 326 | |
327 | fReconstructor -> SetOption("nn"); | |
328 | cTrack.CookPID(); | |
c46a7947 | 329 | Int_t iMomBin = fUtil -> GetMomentumBin(momentum); |
74b2e03d | 330 | |
331 | switch(pdg){ | |
332 | case kElectron: | |
333 | case kPositron: | |
334 | hPIDNN -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
335 | break; | |
336 | case kMuonPlus: | |
337 | case kMuonMinus: | |
338 | hPIDNN -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
339 | break; | |
340 | case kPiPlus: | |
341 | case kPiMinus: | |
342 | hPIDNN -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
343 | break; | |
344 | case kKPlus: | |
345 | case kKMinus: | |
346 | hPIDNN -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
347 | break; | |
348 | case kProton: | |
349 | case kProtonBar: | |
350 | hPIDNN -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, cTrack.GetPID(AliPID::kElectron)); | |
351 | break; | |
352 | } | |
353 | return hPIDNN; | |
354 | } | |
355 | ||
356 | ||
c46a7947 | 357 | //_______________________________________________________ |
358 | TH1 *AliTRDpidChecker::PlotESD(const AliTRDtrackV1 *track) | |
359 | { | |
360 | // | |
361 | // Plot the probabilities for electrons using 2-dim LQ | |
362 | // | |
363 | ||
364 | if(!fESD){ | |
365 | AliWarning("No ESD info available."); | |
366 | return 0x0; | |
367 | } | |
368 | ||
369 | if(track) fTrack = track; | |
370 | if(!fTrack){ | |
371 | AliWarning("No Track defined."); | |
372 | return 0x0; | |
373 | } | |
374 | ||
375 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
376 | ||
377 | if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){ | |
378 | AliWarning("No Histogram defined."); | |
379 | return 0x0; | |
380 | } | |
381 | TH2F *hPIDESD = 0x0; | |
48b51864 | 382 | if(!(hPIDESD = dynamic_cast<TH2F *>(fEfficiency->At(kESD)))){ |
c46a7947 | 383 | AliWarning("No Histogram defined."); |
384 | return 0x0; | |
385 | } | |
386 | ||
387 | ||
388 | Int_t pdg = 0; | |
389 | Float_t momentum = 0.; | |
390 | if(fMC){ | |
391 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
392 | pdg = fMC->GetPDG(); | |
393 | } else { | |
394 | //AliWarning("No MC info available!"); | |
395 | AliTRDtrackV1 cTrack(*fTrack); | |
396 | cTrack.SetReconstructor(fReconstructor); | |
397 | momentum = cTrack.GetMomentum(0); | |
398 | pdg = CalcPDG(&cTrack); | |
399 | } | |
400 | if(momentum < 0.4) return 0x0;; | |
401 | if(momentum > 12.) return 0x0;; | |
402 | ||
403 | ||
404 | Int_t iMomBin = fUtil->GetMomentumBin(momentum); | |
405 | ||
406 | ||
407 | // Double32_t pidESD[AliPID::kSPECIES]; | |
408 | const Double32_t *pidESD = fESD->GetResponseIter(); | |
409 | ||
410 | switch(pdg){ | |
411 | case kElectron: | |
412 | case kPositron: | |
413 | hPIDESD -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, pidESD[0]); | |
414 | break; | |
415 | case kMuonPlus: | |
416 | case kMuonMinus: | |
417 | hPIDESD -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, pidESD[0]); | |
418 | break; | |
419 | case kPiPlus: | |
420 | case kPiMinus: | |
421 | hPIDESD -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, pidESD[0]); | |
422 | break; | |
423 | case kKPlus: | |
424 | case kKMinus: | |
425 | hPIDESD -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, pidESD[0]); | |
426 | break; | |
427 | case kProton: | |
428 | case kProtonBar: | |
429 | hPIDESD -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, pidESD[0]); | |
430 | break; | |
431 | } | |
432 | return hPIDESD; | |
433 | } | |
434 | ||
435 | ||
74b2e03d | 436 | //_______________________________________________________ |
437 | TH1 *AliTRDpidChecker::PlotdEdx(const AliTRDtrackV1 *track) | |
438 | { | |
439 | // | |
440 | // Plot the probabilities for electrons using 2-dim LQ | |
441 | // | |
442 | ||
443 | if(track) fTrack = track; | |
444 | if(!fTrack){ | |
445 | AliWarning("No Track defined."); | |
446 | return 0x0; | |
447 | } | |
1a6d7c9a | 448 | |
74b2e03d | 449 | if(!CheckTrackQuality(fTrack)) return 0x0; |
773d3f8c | 450 | |
74b2e03d | 451 | TH2F *hdEdx; |
452 | if(!(hdEdx = dynamic_cast<TH2F *>(fContainer->At(kdEdx)))){ | |
453 | AliWarning("No Histogram defined."); | |
454 | return 0x0; | |
455 | } | |
456 | ||
773d3f8c | 457 | |
74b2e03d | 458 | Int_t pdg = 0; |
459 | Float_t momentum = 0.; | |
74b2e03d | 460 | if(fMC){ |
461 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
462 | pdg = fMC->GetPDG(); | |
463 | } else { | |
464 | //AliWarning("No MC info available!"); | |
c46a7947 | 465 | AliTRDtrackV1 cTrack(*fTrack); |
466 | cTrack.SetReconstructor(fReconstructor); | |
74b2e03d | 467 | momentum = cTrack.GetMomentum(0); |
468 | pdg = CalcPDG(&cTrack); | |
469 | } | |
c46a7947 | 470 | if(momentum < 0.4) return 0x0;; |
471 | if(momentum > 12.) return 0x0;; | |
472 | ||
473 | Int_t iMomBin = fUtil -> GetMomentumBin(momentum); | |
74b2e03d | 474 | |
74b2e03d | 475 | |
74b2e03d | 476 | |
477 | Float_t SumdEdx[AliTRDgeometry::kNlayer]; | |
5d6dc395 | 478 | AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer]; |
74b2e03d | 479 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0; |
480 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb); | |
481 | ||
482 | Float_t *fdEdx; | |
483 | Float_t dEdxSlice[AliTRDgeometry::kNlayer][AliTRDReconstructor::kLQslices]; | |
484 | ||
485 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
486 | SumdEdx[iChamb] = 0.; | |
487 | fdEdx = TRDtracklet[iChamb] -> GetdEdx(); | |
488 | SumdEdx[iChamb] += fdEdx[0] + fdEdx[1] + fdEdx[2]; | |
489 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){ | |
490 | dEdxSlice[iChamb][iSlice] = fdEdx[iSlice]; | |
491 | } | |
492 | } | |
1a6d7c9a | 493 | |
74b2e03d | 494 | switch(pdg){ |
495 | case kElectron: | |
496 | case kPositron: | |
497 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
498 | hdEdx -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]); | |
499 | break; | |
500 | case kMuonPlus: | |
501 | case kMuonMinus: | |
502 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
503 | hdEdx -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]); | |
504 | break; | |
505 | case kPiPlus: | |
506 | case kPiMinus: | |
507 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
508 | hdEdx -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]); | |
509 | break; | |
510 | case kKPlus: | |
511 | case kKMinus: | |
512 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
513 | hdEdx -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]); | |
514 | break; | |
515 | case kProton: | |
516 | case kProtonBar: | |
517 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
518 | hdEdx -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]); | |
519 | break; | |
520 | } | |
c7cf2032 | 521 | |
74b2e03d | 522 | return hdEdx; |
523 | } | |
773d3f8c | 524 | |
1a6d7c9a | 525 | |
74b2e03d | 526 | //_______________________________________________________ |
527 | TH1 *AliTRDpidChecker::PlotdEdxSlice(const AliTRDtrackV1 *track) | |
528 | { | |
529 | // | |
530 | // Plot the probabilities for electrons using 2-dim LQ | |
531 | // | |
773d3f8c | 532 | |
74b2e03d | 533 | if(track) fTrack = track; |
534 | if(!fTrack){ | |
535 | AliWarning("No Track defined."); | |
536 | return 0x0; | |
537 | } | |
538 | ||
539 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
540 | ||
541 | TH2F *hdEdxSlice; | |
542 | if(!(hdEdxSlice = dynamic_cast<TH2F *>(fContainer->At(kdEdxSlice)))){ | |
543 | AliWarning("No Histogram defined."); | |
544 | return 0x0; | |
545 | } | |
773d3f8c | 546 | |
74b2e03d | 547 | Int_t pdg = 0; |
548 | Float_t momentum = 0.; | |
74b2e03d | 549 | if(fMC){ |
550 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
551 | pdg = fMC->GetPDG(); | |
552 | } else { | |
553 | //AliWarning("No MC info available!"); | |
c46a7947 | 554 | AliTRDtrackV1 cTrack(*fTrack); |
555 | cTrack.SetReconstructor(fReconstructor); | |
74b2e03d | 556 | momentum = cTrack.GetMomentum(0); |
557 | pdg = CalcPDG(&cTrack); | |
558 | } | |
c46a7947 | 559 | if(momentum < 0.4) return 0x0; |
560 | if(momentum > 12.) return 0x0;; | |
561 | ||
562 | Int_t iMomBin = fUtil -> GetMomentumBin(momentum); | |
c7cf2032 | 563 | |
773d3f8c | 564 | |
74b2e03d | 565 | |
566 | AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer]; | |
567 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0; | |
568 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb); | |
569 | ||
570 | Float_t *fdEdx; | |
571 | Float_t dEdxSlice[AliTRDgeometry::kNlayer][AliTRDReconstructor::kLQslices]; | |
572 | ||
573 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
574 | fdEdx = TRDtracklet[iChamb] -> GetdEdx(); | |
575 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){ | |
576 | dEdxSlice[iChamb][iSlice] = fdEdx[iSlice]; | |
c7cf2032 | 577 | } |
74b2e03d | 578 | } |
579 | ||
580 | switch(pdg){ | |
581 | case kElectron: | |
582 | case kPositron: | |
c46a7947 | 583 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ |
584 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){ | |
585 | hdEdxSlice -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice, dEdxSlice[iChamb][iSlice]); | |
586 | } | |
587 | } | |
74b2e03d | 588 | break; |
589 | case kMuonPlus: | |
590 | case kMuonMinus: | |
c46a7947 | 591 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ |
592 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){ | |
593 | hdEdxSlice -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice, | |
594 | dEdxSlice[iChamb][iSlice]); | |
595 | } | |
596 | } | |
74b2e03d | 597 | break; |
598 | case kPiPlus: | |
599 | case kPiMinus: | |
600 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
601 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++) | |
602 | hdEdxSlice -> Fill(AliPID::kPion * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice, | |
603 | dEdxSlice[iChamb][iSlice]); | |
604 | break; | |
605 | case kKPlus: | |
606 | case kKMinus: | |
607 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
608 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++) | |
609 | hdEdxSlice -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice, | |
610 | dEdxSlice[iChamb][iSlice]); | |
611 | break; | |
612 | case kProton: | |
613 | case kProtonBar: | |
614 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
615 | for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++) | |
616 | hdEdxSlice -> Fill(AliPID::kProton * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice, | |
617 | dEdxSlice[iChamb][iSlice]); | |
618 | break; | |
619 | } | |
620 | ||
621 | return hdEdxSlice; | |
622 | ||
623 | } | |
624 | ||
625 | ||
626 | //_______________________________________________________ | |
627 | TH1 *AliTRDpidChecker::PlotPH(const AliTRDtrackV1 *track) | |
628 | { | |
629 | // | |
630 | // Plot the probabilities for electrons using 2-dim LQ | |
631 | // | |
632 | ||
633 | if(track) fTrack = track; | |
634 | if(!fTrack){ | |
635 | AliWarning("No Track defined."); | |
636 | return 0x0; | |
637 | } | |
638 | ||
639 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
640 | ||
641 | TProfile2D *hPH; | |
642 | if(!(hPH = dynamic_cast<TProfile2D *>(fContainer->At(kPH)))){ | |
643 | AliWarning("No Histogram defined."); | |
644 | return 0x0; | |
645 | } | |
646 | ||
74b2e03d | 647 | Int_t pdg = 0; |
648 | Float_t momentum = 0.; | |
74b2e03d | 649 | if(fMC){ |
650 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
651 | pdg = fMC->GetPDG(); | |
652 | } else { | |
653 | //AliWarning("No MC info available!"); | |
c46a7947 | 654 | AliTRDtrackV1 cTrack(*fTrack); |
655 | cTrack.SetReconstructor(fReconstructor); | |
74b2e03d | 656 | momentum = cTrack.GetMomentum(0); |
657 | pdg = CalcPDG(&cTrack); | |
658 | } | |
74b2e03d | 659 | if(momentum < 0.4) return 0x0;; |
c46a7947 | 660 | if(momentum > 12.) return 0x0;; |
661 | ||
662 | Int_t iMomBin = fUtil -> GetMomentumBin(momentum); | |
c4c5bbfb | 663 | |
74b2e03d | 664 | AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer]; |
665 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0; | |
666 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb); | |
3afc1da3 | 667 | |
74b2e03d | 668 | AliTRDcluster *TRDcluster = 0x0; |
3afc1da3 | 669 | |
74b2e03d | 670 | switch(pdg){ |
671 | case kElectron: | |
672 | case kPositron: | |
5d6dc395 | 673 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ |
74b2e03d | 674 | for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){ |
675 | if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus))) | |
676 | continue; | |
677 | hPH -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus)); | |
b718144c | 678 | } |
1a6d7c9a | 679 | } |
74b2e03d | 680 | break; |
681 | case kMuonPlus: | |
682 | case kMuonMinus: | |
683 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
684 | for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){ | |
685 | if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus))) | |
686 | continue; | |
687 | hPH -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus)); | |
1a6d7c9a | 688 | } |
74b2e03d | 689 | } |
690 | break; | |
691 | case kPiPlus: | |
692 | case kPiMinus: | |
693 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
694 | for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){ | |
695 | if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus))) | |
696 | continue; | |
697 | hPH -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus)); | |
1a6d7c9a | 698 | } |
74b2e03d | 699 | } |
700 | break; | |
701 | case kKPlus: | |
702 | case kKMinus: | |
703 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
704 | for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){ | |
705 | if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus))) | |
706 | continue; | |
707 | hPH -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus)); | |
1a6d7c9a | 708 | } |
74b2e03d | 709 | } |
710 | break; | |
711 | case kProton: | |
712 | case kProtonBar: | |
713 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
714 | for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){ | |
715 | if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus))) | |
716 | continue; | |
717 | hPH -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus)); | |
1a6d7c9a | 718 | } |
773d3f8c | 719 | } |
74b2e03d | 720 | break; |
721 | } | |
722 | ||
723 | return hPH; | |
724 | } | |
773d3f8c | 725 | |
726 | ||
74b2e03d | 727 | //_______________________________________________________ |
728 | TH1 *AliTRDpidChecker::PlotNClus(const AliTRDtrackV1 *track) | |
729 | { | |
730 | // | |
731 | // Plot the probabilities for electrons using 2-dim LQ | |
732 | // | |
733 | ||
734 | if(track) fTrack = track; | |
735 | if(!fTrack){ | |
736 | AliWarning("No Track defined."); | |
737 | return 0x0; | |
738 | } | |
739 | ||
740 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
741 | ||
742 | TH2F *hNClus; | |
743 | if(!(hNClus = dynamic_cast<TH2F *>(fContainer->At(kNClus)))){ | |
744 | AliWarning("No Histogram defined."); | |
745 | return 0x0; | |
746 | } | |
747 | ||
748 | ||
74b2e03d | 749 | Int_t pdg = 0; |
750 | Float_t momentum = 0.; | |
74b2e03d | 751 | if(fMC){ |
752 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
753 | pdg = fMC->GetPDG(); | |
754 | } else { | |
755 | //AliWarning("No MC info available!"); | |
c46a7947 | 756 | AliTRDtrackV1 cTrack(*fTrack); |
757 | cTrack.SetReconstructor(fReconstructor); | |
74b2e03d | 758 | momentum = cTrack.GetMomentum(0); |
759 | pdg = CalcPDG(&cTrack); | |
760 | } | |
c46a7947 | 761 | if(momentum < 0.4) return 0x0;; |
762 | if(momentum > 12.) return 0x0;; | |
74b2e03d | 763 | |
c46a7947 | 764 | Int_t iMomBin = fUtil -> GetMomentumBin(momentum); |
74b2e03d | 765 | |
74b2e03d | 766 | |
767 | Int_t iNClus[AliTRDgeometry::kNlayer]; | |
768 | memset(iNClus, 0, sizeof(Int_t) * AliTRDgeometry::kNlayer); | |
769 | ||
770 | AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer]; | |
771 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) TRDtracklet[iChamb] = 0x0; | |
772 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){ | |
773 | TRDtracklet[iChamb] = fTrack -> GetTracklet(iChamb); | |
774 | iNClus[iChamb] = TRDtracklet[iChamb] -> GetN(); | |
775 | } | |
776 | ||
777 | switch(pdg){ | |
778 | case kElectron: | |
779 | case kPositron: | |
780 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
781 | hNClus -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]); | |
782 | break; | |
783 | case kMuonPlus: | |
784 | case kMuonMinus: | |
785 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
786 | hNClus -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]); | |
787 | break; | |
788 | case kPiPlus: | |
789 | case kPiMinus: | |
790 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
791 | hNClus -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]); | |
792 | break; | |
793 | case kKPlus: | |
794 | case kKMinus: | |
795 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
796 | hNClus -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]); | |
797 | break; | |
798 | case kProton: | |
799 | case kProtonBar: | |
800 | for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++) | |
801 | hNClus -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]); | |
802 | break; | |
803 | } | |
804 | ||
805 | return hNClus; | |
806 | } | |
807 | ||
808 | ||
809 | //_______________________________________________________ | |
810 | TH1 *AliTRDpidChecker::PlotMom(const AliTRDtrackV1 *track) | |
811 | { | |
812 | // | |
813 | // Plot the probabilities for electrons using 2-dim LQ | |
814 | // | |
422a2dc0 | 815 | |
74b2e03d | 816 | if(track) fTrack = track; |
817 | if(!fTrack){ | |
818 | AliWarning("No Track defined."); | |
819 | return 0x0; | |
820 | } | |
821 | ||
822 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
823 | ||
824 | TH1F *hMom; | |
825 | if(!(hMom = dynamic_cast<TH1F *>(fContainer->At(kMomentum)))){ | |
826 | AliWarning("No Histogram defined."); | |
827 | return 0x0; | |
773d3f8c | 828 | } |
829 | ||
74b2e03d | 830 | |
74b2e03d | 831 | Int_t pdg = 0; |
832 | Float_t momentum = 0.; | |
74b2e03d | 833 | if(fMC){ |
834 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
835 | pdg = fMC->GetPDG(); | |
836 | } else { | |
837 | //AliWarning("No MC info available!"); | |
c46a7947 | 838 | AliTRDtrackV1 cTrack(*fTrack); |
839 | cTrack.SetReconstructor(fReconstructor); | |
74b2e03d | 840 | momentum = cTrack.GetMomentum(0); |
841 | pdg = CalcPDG(&cTrack); | |
842 | } | |
74b2e03d | 843 | if(momentum < 0.4) return 0x0; |
c46a7947 | 844 | if(momentum > 12.) return 0x0;; |
74b2e03d | 845 | |
846 | hMom -> Fill(momentum); | |
847 | return hMom; | |
848 | } | |
849 | ||
850 | ||
851 | //_______________________________________________________ | |
852 | TH1 *AliTRDpidChecker::PlotMomBin(const AliTRDtrackV1 *track) | |
853 | { | |
854 | // | |
855 | // Plot the probabilities for electrons using 2-dim LQ | |
856 | // | |
857 | ||
858 | if(track) fTrack = track; | |
859 | if(!fTrack){ | |
860 | AliWarning("No Track defined."); | |
861 | return 0x0; | |
862 | } | |
863 | ||
864 | if(!CheckTrackQuality(fTrack)) return 0x0; | |
865 | ||
866 | TH1F *hMomBin; | |
867 | if(!(hMomBin = dynamic_cast<TH1F *>(fContainer->At(kMomentumBin)))){ | |
868 | AliWarning("No Histogram defined."); | |
869 | return 0x0; | |
870 | } | |
871 | ||
872 | ||
74b2e03d | 873 | Int_t pdg = 0; |
874 | Float_t momentum = 0.; | |
875 | ||
876 | if(fMC){ | |
877 | if(fMC->GetTrackRef()) momentum = fMC->GetTrackRef()->P(); | |
878 | pdg = fMC->GetPDG(); | |
879 | } else { | |
880 | //AliWarning("No MC info available!"); | |
c46a7947 | 881 | AliTRDtrackV1 cTrack(*fTrack); |
882 | cTrack.SetReconstructor(fReconstructor); | |
74b2e03d | 883 | momentum = cTrack.GetMomentum(0); |
884 | pdg = CalcPDG(&cTrack); | |
885 | } | |
74b2e03d | 886 | if(momentum < 0.4) return 0x0; |
c46a7947 | 887 | if(momentum > 12.) return 0x0;; |
74b2e03d | 888 | |
c46a7947 | 889 | Int_t iMomBin = fUtil -> GetMomentumBin(momentum); |
74b2e03d | 890 | hMomBin -> Fill(iMomBin); |
891 | return hMomBin; | |
773d3f8c | 892 | } |
893 | ||
a391a274 | 894 | |
c4c5bbfb | 895 | //________________________________________________________ |
e15179be | 896 | Bool_t AliTRDpidChecker::GetRefFigure(Int_t ifig) |
c4c5bbfb | 897 | { |
a391a274 | 898 | Bool_t FIRST = kTRUE; |
c46a7947 | 899 | TLegend *leg = new TLegend(.7, .7, .98, .98); |
900 | leg->SetBorderSize(1); | |
c026e9cc | 901 | TGraphErrors *g = 0x0; |
c46a7947 | 902 | TAxis *ax = 0x0; |
903 | TH1 *h1 = 0x0, *h=0x0; | |
a391a274 | 904 | TH2 *h2 = 0x0; |
ef66c5d4 | 905 | TObjArray *arr = 0x0; |
c4c5bbfb | 906 | switch(ifig){ |
c46a7947 | 907 | case kEfficiency: |
ef66c5d4 | 908 | arr = (TObjArray*)fGraph->At(ifig); |
909 | if(!(g = (TGraphErrors*)arr->At(kLQ))) break; | |
c46a7947 | 910 | if(!g->GetN()) break; |
911 | leg->SetHeader("PID Method"); | |
c026e9cc | 912 | g->Draw("apl"); |
c46a7947 | 913 | ax = g->GetHistogram()->GetXaxis(); |
914 | ax->SetTitle("p [GeV/c]"); | |
915 | ax->SetRangeUser(.6, 10.5); | |
916 | ax->SetMoreLogLabels(); | |
afdc0e7a | 917 | ax = g->GetHistogram()->GetYaxis(); |
918 | ax->SetTitle("#epsilon_{#pi} [%]"); | |
919 | ax->SetRangeUser(1.e-3, 1.e-1); | |
c46a7947 | 920 | leg->AddEntry(g, "2D LQ", "pl"); |
ef66c5d4 | 921 | if(! (g = (TGraphErrors*)arr->At(kNN))) break; |
c46a7947 | 922 | g->Draw("pl"); |
923 | leg->AddEntry(g, "NN", "pl"); | |
ef66c5d4 | 924 | if(! (g = (TGraphErrors*)arr->At(kESD))) break; |
afdc0e7a | 925 | g->Draw("p"); |
c46a7947 | 926 | leg->AddEntry(g, "ESD", "pl"); |
927 | leg->Draw(); | |
a391a274 | 928 | gPad->SetLogy(); |
c026e9cc | 929 | gPad->SetLogx(); |
3afc1da3 | 930 | gPad->SetGridy(); |
931 | gPad->SetGridx(); | |
e15179be | 932 | return kTRUE; |
c46a7947 | 933 | case kdEdx: |
a391a274 | 934 | // save 2.0 GeV projection as reference |
935 | FIRST = kTRUE; | |
c46a7947 | 936 | if(!(h2 = (TH2F*)(fContainer->At(kdEdx)))) break; |
937 | leg->SetHeader("Particle Species"); | |
a391a274 | 938 | for(Int_t is = AliPID::kSPECIES-1; is>=0; is--){ |
939 | Int_t bin = is*AliTRDCalPID::kNMom+4; | |
940 | h1 = h2->ProjectionY("px", bin, bin); | |
941 | if(!h1->GetEntries()) continue; | |
942 | h1->Scale(1./h1->Integral()); | |
943 | h1->SetLineColor(AliTRDCalPID::GetPartColor(is)); | |
c46a7947 | 944 | h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec"); |
afdc0e7a | 945 | leg->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l"); |
a391a274 | 946 | FIRST = kFALSE; |
947 | } | |
c46a7947 | 948 | if(FIRST) break; |
949 | leg->Draw(); | |
a391a274 | 950 | gPad->SetLogy(); |
c026e9cc | 951 | gPad->SetLogx(0); |
3afc1da3 | 952 | gPad->SetGridy(); |
953 | gPad->SetGridx(); | |
e15179be | 954 | return kTRUE; |
c46a7947 | 955 | case kdEdxSlice: |
a391a274 | 956 | break; |
c46a7947 | 957 | case kPH: |
a391a274 | 958 | // save 2.0 GeV projection as reference |
959 | FIRST = kTRUE; | |
c46a7947 | 960 | if(!(h2 = (TH2F*)(fContainer->At(kPH)))) break;; |
961 | leg->SetHeader("Particle Species"); | |
a391a274 | 962 | for(Int_t is=0; is<AliPID::kSPECIES; is++){ |
963 | Int_t bin = is*AliTRDCalPID::kNMom+4; | |
964 | h1 = h2->ProjectionY("py", bin, bin); | |
965 | if(!h1->GetEntries()) continue; | |
966 | h1->SetMarkerStyle(24); | |
967 | h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is)); | |
968 | h1->SetLineColor(AliTRDCalPID::GetPartColor(is)); | |
afdc0e7a | 969 | if(FIRST){ |
970 | h1->GetXaxis()->SetTitle("tb[1/100 ns^{-1}]"); | |
971 | h1->GetYaxis()->SetTitle("<PH> [a.u.]"); | |
972 | } | |
973 | h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec"); | |
c46a7947 | 974 | leg->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "pl"); |
a391a274 | 975 | FIRST = kFALSE; |
976 | } | |
c46a7947 | 977 | if(FIRST) break; |
978 | leg->Draw(); | |
a391a274 | 979 | gPad->SetLogy(0); |
c026e9cc | 980 | gPad->SetLogx(0); |
3afc1da3 | 981 | gPad->SetGridy(); |
982 | gPad->SetGridx(); | |
e15179be | 983 | return kTRUE; |
c46a7947 | 984 | case kNClus: |
3afc1da3 | 985 | // save 2.0 GeV projection as reference |
986 | FIRST = kTRUE; | |
c46a7947 | 987 | if(!(h2 = (TH2F*)(fContainer->At(kNClus)))) break; |
988 | leg->SetHeader("Particle Species"); | |
3afc1da3 | 989 | for(Int_t is=0; is<AliPID::kSPECIES; is++){ |
990 | Int_t bin = is*AliTRDCalPID::kNMom+4; | |
991 | h1 = h2->ProjectionY("py", bin, bin); | |
992 | if(!h1->GetEntries()) continue; | |
afdc0e7a | 993 | //h1->SetMarkerStyle(24); |
994 | //h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is)); | |
3afc1da3 | 995 | h1->SetLineColor(AliTRDCalPID::GetPartColor(is)); |
afdc0e7a | 996 | if(FIRST) h1->GetXaxis()->SetTitle("N^{cl}/tracklet"); |
997 | h = (TH1F*)h1->DrawClone(FIRST ? "c" : "samec"); | |
998 | leg->AddEntry(h, Form("%s", AliTRDCalPID::GetPartName(is)), "l"); | |
3afc1da3 | 999 | FIRST = kFALSE; |
1000 | } | |
c46a7947 | 1001 | if(FIRST) break; |
1002 | leg->Draw(); | |
afdc0e7a | 1003 | gPad->SetLogy(); |
3afc1da3 | 1004 | gPad->SetLogx(0); |
1005 | gPad->SetGridy(); | |
1006 | gPad->SetGridx(); | |
e15179be | 1007 | return kTRUE; |
c46a7947 | 1008 | case kMomentum: |
1009 | case kMomentumBin: | |
1010 | break; | |
1011 | case kThresh: | |
ef66c5d4 | 1012 | arr = (TObjArray*)fGraph->FindObject("Thresholds"); |
1013 | if(!(g = (TGraphErrors*)arr->At(kLQ))) break; | |
c46a7947 | 1014 | if(!g->GetN()) break; |
1015 | leg->SetHeader("PID Method"); | |
1016 | g->Draw("apl"); | |
afdc0e7a | 1017 | ax = g->GetHistogram()->GetXaxis(); |
1018 | ax->SetTitle("p [GeV/c]"); | |
1019 | ax->SetRangeUser(.6, 10.5); | |
1020 | ax->SetMoreLogLabels(); | |
1021 | ax = g->GetHistogram()->GetYaxis(); | |
1022 | ax->SetTitle("threshold"); | |
1023 | ax->SetRangeUser(5.e-2, 1.); | |
c46a7947 | 1024 | leg->AddEntry(g, "2D LQ", "pl"); |
ef66c5d4 | 1025 | if(!(g = (TGraphErrors*)arr->At(kNN))) break; |
c46a7947 | 1026 | g->Draw("pl"); |
1027 | leg->AddEntry(g, "NN", "pl"); | |
ef66c5d4 | 1028 | if(!(g = (TGraphErrors*)arr->At(kESD))) break; |
afdc0e7a | 1029 | g->Draw("p"); |
c46a7947 | 1030 | leg->AddEntry(g, "ESD", "pl"); |
1031 | leg->Draw(); | |
1032 | gPad->SetLogx(); | |
1033 | gPad->SetGridy(); | |
1034 | gPad->SetGridx(); | |
e15179be | 1035 | return kTRUE; |
c4c5bbfb | 1036 | } |
c46a7947 | 1037 | AliInfo(Form("Reference plot [%d] missing result", ifig)); |
e15179be | 1038 | return kFALSE; |
c4c5bbfb | 1039 | } |
1040 | ||
773d3f8c | 1041 | //________________________________________________________________________ |
d85cd79c | 1042 | Bool_t AliTRDpidChecker::PostProcess() |
773d3f8c | 1043 | { |
c7cf2032 | 1044 | // Draw result to the screen |
1045 | // Called once at the end of the query | |
1046 | ||
3d86166d | 1047 | if (!fContainer) { |
773d3f8c | 1048 | Printf("ERROR: list not available"); |
d85cd79c | 1049 | return kFALSE; |
773d3f8c | 1050 | } |
c46a7947 | 1051 | if(!(fEfficiency = dynamic_cast<TObjArray *>(fContainer->At(kEfficiency)))){ |
1052 | AliError("Efficiency container missing."); | |
1053 | return 0x0; | |
1054 | } | |
c46a7947 | 1055 | if(!fGraph){ |
ef66c5d4 | 1056 | fGraph = new TObjArray(2); |
c46a7947 | 1057 | fGraph->SetOwner(); |
c46a7947 | 1058 | } |
ef66c5d4 | 1059 | EvaluatePionEfficiency(fEfficiency, fGraph, 0.9); |
1bd44ae8 | 1060 | fNRefFigures = 8; |
1061 | return kTRUE; | |
1062 | } | |
773d3f8c | 1063 | |
1bd44ae8 | 1064 | //________________________________________________________________________ |
ef66c5d4 | 1065 | void AliTRDpidChecker::EvaluatePionEfficiency(TObjArray *histoContainer, TObjArray *results, Float_t electron_efficiency) |
1066 | { | |
1067 | if(!histoContainer || !results) return; | |
1068 | ||
1bd44ae8 | 1069 | TGraphErrors *g = 0x0; |
1070 | fUtil->SetElectronEfficiency(electron_efficiency); | |
1071 | ||
1072 | // efficiency graphs | |
ef66c5d4 | 1073 | TObjArray *eff = new TObjArray(3); eff->SetOwner(); eff->SetName("Efficiencies"); |
1074 | results->AddAt(eff, 0); | |
1075 | eff->AddAt(g = new TGraphErrors(), kLQ); | |
1bd44ae8 | 1076 | g->SetLineColor(kBlue); |
1077 | g->SetMarkerColor(kBlue); | |
1078 | g->SetMarkerStyle(7); | |
ef66c5d4 | 1079 | eff->AddAt(g = new TGraphErrors(), kNN); |
1bd44ae8 | 1080 | g->SetLineColor(kGreen); |
1081 | g->SetMarkerColor(kGreen); | |
1082 | g->SetMarkerStyle(7); | |
ef66c5d4 | 1083 | eff -> AddAt(g = new TGraphErrors(), kESD); |
1bd44ae8 | 1084 | g->SetLineColor(kRed); |
1085 | g->SetMarkerColor(kRed); | |
1086 | g->SetMarkerStyle(24); | |
1087 | ||
ef66c5d4 | 1088 | // Threshold graphs |
1089 | TObjArray *thres = new TObjArray(3); thres->SetOwner(); thres->SetName("Thresholds"); | |
1090 | results->AddAt(thres, 1); | |
1091 | thres->AddAt(g = new TGraphErrors(), kLQ); | |
1bd44ae8 | 1092 | g->SetLineColor(kBlue); |
1093 | g->SetMarkerColor(kBlue); | |
1094 | g->SetMarkerStyle(7); | |
ef66c5d4 | 1095 | thres->AddAt(g = new TGraphErrors(), kNN); |
1bd44ae8 | 1096 | g->SetLineColor(kGreen); |
1097 | g->SetMarkerColor(kGreen); | |
1098 | g->SetMarkerStyle(7); | |
ef66c5d4 | 1099 | thres -> AddAt(g = new TGraphErrors(), kESD); |
1bd44ae8 | 1100 | g->SetLineColor(kRed); |
1101 | g->SetMarkerColor(kRed); | |
1102 | g->SetMarkerStyle(24); | |
1103 | ||
422a2dc0 | 1104 | Float_t mom = 0.; |
422a2dc0 | 1105 | TH1D *Histo1=0x0, *Histo2=0x0; |
1106 | ||
773d3f8c | 1107 | // calculate the pion efficiencies and the errors for 90% electron efficiency (2-dim LQ) |
1bd44ae8 | 1108 | TH2F *hPIDLQ = (TH2F*)histoContainer->At(kLQ); |
773d3f8c | 1109 | for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){ |
422a2dc0 | 1110 | mom = AliTRDCalPID::GetMomentum(iMom); |
1111 | ||
1112 | Histo1 = hPIDLQ -> ProjectionY("LQ_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1); | |
1113 | Histo2 = hPIDLQ -> ProjectionY("LQ_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1); | |
1114 | ||
c46a7947 | 1115 | if(!fUtil->CalculatePionEffi(Histo1, Histo2)) continue; |
422a2dc0 | 1116 | |
ef66c5d4 | 1117 | g = (TGraphErrors*)eff->At(kLQ); |
c46a7947 | 1118 | g->SetPoint(iMom, mom, fUtil->GetPionEfficiency()); |
1119 | g->SetPointError(iMom, 0., fUtil->GetError()); | |
ef66c5d4 | 1120 | g = (TGraphErrors*)thres->At(kLQ); |
c46a7947 | 1121 | g->SetPoint(iMom, mom, fUtil->GetThreshold()); |
1122 | g->SetPointError(iMom, 0., 0.); | |
422a2dc0 | 1123 | |
c46a7947 | 1124 | if(fDebugLevel>=2) Printf("Pion Efficiency for 2-dim LQ is : %f +/- %f\n\n", fUtil->GetPionEfficiency(), fUtil->GetError()); |
773d3f8c | 1125 | } |
1126 | ||
422a2dc0 | 1127 | |
773d3f8c | 1128 | // calculate the pion efficiencies and the errors for 90% electron efficiency (NN) |
1bd44ae8 | 1129 | TH2F *hPIDNN = (TH2F*)histoContainer->At(kNN); |
773d3f8c | 1130 | for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){ |
422a2dc0 | 1131 | mom = AliTRDCalPID::GetMomentum(iMom); |
1132 | ||
1133 | Histo1 = hPIDNN -> ProjectionY("NN_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1); | |
1134 | Histo2 = hPIDNN -> ProjectionY("NN_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1); | |
1135 | ||
c46a7947 | 1136 | if(!fUtil -> CalculatePionEffi(Histo1, Histo2)) continue; |
422a2dc0 | 1137 | |
ef66c5d4 | 1138 | g = (TGraphErrors*)eff->At(kNN); |
c46a7947 | 1139 | g->SetPoint(iMom, mom, fUtil->GetPionEfficiency()); |
1140 | g->SetPointError(iMom, 0., fUtil->GetError()); | |
014bede1 | 1141 | g = (TGraphErrors*)thres->At(kNN); |
c46a7947 | 1142 | g->SetPoint(iMom, mom, fUtil->GetThreshold()); |
1143 | g->SetPointError(iMom, 0., 0.); | |
422a2dc0 | 1144 | |
c46a7947 | 1145 | if(fDebugLevel>=2) Printf("Pion Efficiency for NN is : %f +/- %f\n\n", fUtil->GetPionEfficiency(), fUtil->GetError()); |
773d3f8c | 1146 | } |
773d3f8c | 1147 | |
773d3f8c | 1148 | |
c46a7947 | 1149 | // calculate the pion efficiencies and the errors for 90% electron efficiency (ESD) |
1bd44ae8 | 1150 | TH2F *hPIDESD = (TH2F*)histoContainer->At(kESD); |
c46a7947 | 1151 | for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){ |
1152 | mom = AliTRDCalPID::GetMomentum(iMom); | |
1153 | ||
1154 | Histo1 = hPIDESD -> ProjectionY("NN_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1); | |
1155 | Histo2 = hPIDESD -> ProjectionY("NN_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1); | |
773d3f8c | 1156 | |
c46a7947 | 1157 | if(!fUtil->CalculatePionEffi(Histo1, Histo2)) continue; |
422a2dc0 | 1158 | |
ef66c5d4 | 1159 | g = (TGraphErrors*)eff->At(kESD); |
c46a7947 | 1160 | g->SetPoint(iMom, mom, fUtil->GetPionEfficiency()); |
1161 | g->SetPointError(iMom, 0., fUtil->GetError()); | |
014bede1 | 1162 | g = (TGraphErrors*)thres->At(kESD); |
c46a7947 | 1163 | g->SetPoint(iMom, mom, fUtil->GetThreshold()); |
1164 | g->SetPointError(iMom, 0., 0.); | |
773d3f8c | 1165 | |
c46a7947 | 1166 | if(fDebugLevel>=2) Printf("Pion Efficiency for ESD is : %f +/- %f\n\n", fUtil->GetPionEfficiency(), fUtil->GetError()); |
773d3f8c | 1167 | } |
1168 | ||
d85cd79c | 1169 | } |
1170 | ||
1171 | //________________________________________________________________________ | |
1172 | void AliTRDpidChecker::Terminate(Option_t *) | |
1173 | { | |
1174 | // Draw result to the screen | |
1175 | // Called once at the end of the query | |
1176 | ||
1177 | fContainer = dynamic_cast<TObjArray*>(GetOutputData(0)); | |
1178 | if (!fContainer) { | |
1179 | Printf("ERROR: list not available"); | |
1180 | return; | |
1181 | } | |
773d3f8c | 1182 | } |
1183 | ||
1184 |