]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TRD/qaRec/AliTRDpidChecker.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
improved tooltips
[u/mrichter/AliRoot.git] / TRD / qaRec / AliTRDpidChecker.cxx
CommitLineData
773d3f8c 1#include "TPDGCode.h"
a391a274 2#include "TCanvas.h"
773d3f8c 3#include "TF1.h"
4#include "TH1F.h"
422a2dc0 5#include "TH1D.h"
773d3f8c 6#include "TH2F.h"
c7cf2032 7#include "TProfile.h"
422a2dc0 8#include "TProfile2D.h"
773d3f8c 9#include "TGraph.h"
10#include "TGraphErrors.h"
11
12#include <TClonesArray.h>
13#include <TObjArray.h>
14#include <TList.h>
15
c7cf2032 16// #include "AliPID.h"
773d3f8c 17#include "AliESDEvent.h"
18#include "AliESDInputHandler.h"
19#include "AliTrackReference.h"
20
c7cf2032 21#include "AliAnalysisTask.h"
22// #include "AliAnalysisManager.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"
c7cf2032 29// #include "../Cal/AliTRDCalPID.h"
422a2dc0 30#include "AliTRDpidUtil.h"
773d3f8c 31
32
33#include "AliTRDpidChecker.h"
34#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
35
36// calculate pion efficiency at 90% electron efficiency for 11 momentum bins
37// this task should be used with simulated data only
38
39ClassImp(AliTRDpidChecker)
40
41//________________________________________________________________________
3d86166d 42AliTRDpidChecker::AliTRDpidChecker()
43 :AliTRDrecoTask("PID", "PID Checker")
773d3f8c 44 ,fReconstructor(0x0)
773d3f8c 45{
46 //
47 // Default constructor
48 //
49
50 fReconstructor = new AliTRDReconstructor();
51 fReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
773d3f8c 52}
53
54
55//________________________________________________________________________
56AliTRDpidChecker::~AliTRDpidChecker()
57{
3d86166d 58 if(fReconstructor) delete fReconstructor;
773d3f8c 59}
60
61
62//________________________________________________________________________
63void AliTRDpidChecker::CreateOutputObjects()
64{
65 // Create histograms
66 // Called once
67
68 OpenFile(0, "RECREATE");
422a2dc0 69 Int_t xBins = AliPID::kSPECIES*AliTRDCalPID::kNMom;
3d86166d 70 fContainer = new TObjArray();
422a2dc0 71 fContainer -> Expand(kGraphNN + 1);
773d3f8c 72
422a2dc0 73 const Float_t epsilon = 1./(2*(AliTRDpidUtil::kBins-1)); // get nice histos with bin center at 0 and 1
773d3f8c 74
75 // histos of the electron probability of all 5 particle species and 11 momenta for the 2-dim LQ method
422a2dc0 76 fContainer->AddAt(
77 new TH2F("PID_LQ", "",
78 xBins, -0.5, xBins - 0.5,
79 AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon)
80 ,kLQlikelihood);
81
773d3f8c 82
1a6d7c9a 83 // histos of the electron probability of all 5 particle species and 11 momenta for the neural network method
422a2dc0 84 fContainer->AddAt(
85 new TH2F("PID_NN", "",
86 xBins, -0.5, xBins - 0.5,
87 AliTRDpidUtil::kBins, 0.-epsilon, 1.+epsilon)
88 ,kNNlikelihood);
1a6d7c9a 89
90 // histos of the dE/dx distribution for all 5 particle species and 11 momenta
422a2dc0 91 fContainer->AddAt(
92 new TH2F("dEdx", "",
93 xBins, -0.5, xBins - 0.5,
94 200, 0, 10000)
95 ,kdEdx);
773d3f8c 96
b718144c 97 // histos of the dE/dx distribution for all 5 particle species and 11 momenta
98 fContainer->AddAt(
99 new TH2F("dEdxSlice", "",
100 xBins*AliTRDReconstructor::kLQslices, -0.5, xBins*AliTRDReconstructor::kLQslices - 0.5,
101 200, 0, 5000)
102 ,kdEdxSlice);
103
c7cf2032 104 // histos of the pulse height distribution for all 5 particle species and 11 momenta
422a2dc0 105 fContainer->AddAt(
106 new TProfile2D("PH", "",
107 xBins, -0.5, xBins - 0.5,
108 AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5)
109 ,kPH);
110
3afc1da3 111 // histos of the number of clusters distribution for all 5 particle species and 11 momenta
112 fContainer->AddAt(
113 new TH2F("NClus", "",
114 xBins, -0.5, xBins - 0.5,
115 AliTRDtrackerV1::GetNTimeBins(), -0.5, AliTRDtrackerV1::GetNTimeBins() - 0.5)
116 ,kNClus);
117
c7cf2032 118
119 // momentum distributions - absolute and in momentum bins
120 fContainer->AddAt(new TH1F("hMom", "momentum distribution", 100, 0., 12.),kMomentum);
121 fContainer->AddAt(new TH1F("hMomBin", "momentum distribution in momentum bins", AliTRDCalPID::kNMom, 0.5, 11.5),kMomentumBin);
122
c4c5bbfb 123
124 // TGraph of the pion efficiencies
125
422a2dc0 126 TGraphErrors *gEffisLQ = 0x0;
127 TGraphErrors *gEffisNN = 0x0;
c4c5bbfb 128
422a2dc0 129 fContainer->AddAt(gEffisLQ = new TGraphErrors(), kGraphLQ);
130 gEffisLQ->SetLineColor(kBlue);
131 gEffisLQ->SetMarkerColor(kBlue);
132 gEffisLQ->SetMarkerStyle(29);
c4c5bbfb 133
422a2dc0 134 fContainer -> AddAt(gEffisNN = new TGraphErrors(),kGraphNN);
135 gEffisNN->SetLineColor(kRed);
136 gEffisNN->SetMarkerColor(kRed);
137 gEffisNN->SetMarkerStyle(29);
c4c5bbfb 138
422a2dc0 139 gEffisLQ -> SetNameTitle("gEffisLQErr", "Efficiencies and Errors of the 2-dim LQ method");
140 gEffisNN -> SetNameTitle("gEffisNNErr", "Efficiencies and Errors of the NN method");
c4c5bbfb 141
773d3f8c 142}
143
144//________________________________________________________________________
145void AliTRDpidChecker::Exec(Option_t *)
146{
147 // Main loop
148 // Called for each event
149
150
151// if(!AliTracker::GetFieldMap()){
152// AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
153// AliTracker::SetFieldMap(field, kTRUE);
154// }
155
422a2dc0 156 TH2F *hPIDLQ;
157 TH2F *hPIDNN;
158 TH2F *hdEdx;
b718144c 159 TH2F *hdEdxSlice;
422a2dc0 160 TProfile2D *hPH;
3afc1da3 161 TH2F *hNClus;
422a2dc0 162
163 hPIDLQ = (TH2F*)fContainer->At(kLQlikelihood);
164 hPIDNN = (TH2F*)fContainer->At(kNNlikelihood);
165 hdEdx = (TH2F*)fContainer->At(kdEdx);
b718144c 166 hdEdxSlice = (TH2F*)fContainer->At(kdEdxSlice);
422a2dc0 167 hPH = (TProfile2D*)fContainer->At(kPH);
3afc1da3 168 hNClus = (TH2F*)fContainer->At(kNClus);
422a2dc0 169
c7cf2032 170 TH1F *hMom = (TH1F*)fContainer->At(kMomentum);
171 TH1F *hMomBin = (TH1F*)fContainer->At(kMomentumBin);
1a6d7c9a 172
773d3f8c 173 Int_t labelsacc[10000];
174 memset(labelsacc, 0, sizeof(Int_t) * 10000);
175
176 Float_t mom;
177 ULong_t status;
178 Int_t nTRD = 0;
1a6d7c9a 179 Float_t *fdEdx;
773d3f8c 180
181 AliTRDtrackInfo *track = 0x0;
c7cf2032 182 AliTRDtrackV1 *TRDtrack = 0x0;
773d3f8c 183 AliTrackReference *ref = 0x0;
184 AliExternalTrackParam *esd = 0x0;
1a6d7c9a 185
5d6dc395 186 AliTRDseedV1 *TRDtracklet[AliTRDgeometry::kNlayer];
187 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++)
1a6d7c9a 188 TRDtracklet[iChamb] = 0x0;
189
c7cf2032 190 AliTRDcluster *TRDcluster = 0x0;
191
422a2dc0 192 AliTRDpidUtil *util = new AliTRDpidUtil();
773d3f8c 193 for(Int_t itrk=0; itrk<fTracks->GetEntriesFast(); itrk++){
194 track = (AliTRDtrackInfo*)fTracks->UncheckedAt(itrk);
195 if(!track->HasESDtrack()) continue;
196 status = track->GetStatus();
197 if(!(status&AliESDtrack::kTPCout)) continue;
198
b718144c 199 if(!(TRDtrack = track->GetTrack())) continue;
773d3f8c 200 //&&(track->GetNumberOfClustersRefit()
1a6d7c9a 201
202 // use only tracks that hit 6 chambers
5d6dc395 203 if(!(TRDtrack->GetNumberOfTracklets() == AliTRDgeometry::kNlayer)) continue;
c4c5bbfb 204
773d3f8c 205 ref = track->GetTrackRef(0);
206 esd = track->GetOuterParam();
207 mom = ref ? ref->P(): esd->P();
208
209 labelsacc[nTRD] = track->GetLabel();
210 nTRD++;
211
773d3f8c 212
213 // set the 11 momentum bins
214 Int_t iMomBin = -1;
422a2dc0 215 iMomBin = util -> GetMomentumBin(mom);
216 if(fDebugLevel>=4) Printf("MomBin[%d] MomTot[%f]", iMomBin, mom);
217
773d3f8c 218
c7cf2032 219 // fill momentum histo to have the momentum distribution
220 hMom -> Fill(mom);
221 hMomBin -> Fill(iMomBin);
222
223
773d3f8c 224 // set the reconstructor
225 TRDtrack -> SetReconstructor(fReconstructor);
226
227
c7cf2032 228 // if no monte carlo data available -> use TRDpid
229 if(!HasMCdata()){
230 fReconstructor -> SetOption("nn");
231 TRDtrack -> CookPID();
232 if(TRDtrack -> GetPID(0) > TRDtrack -> GetPID(1) + TRDtrack -> GetPID(2) + TRDtrack -> GetPID(3) + TRDtrack -> GetPID(4)){
c4c5bbfb 233 track -> SetPDG(kElectron);
234 }
235 else if(TRDtrack -> GetPID(4) > TRDtrack -> GetPID(2) && TRDtrack -> GetPID(4) > TRDtrack -> GetPID(3) && TRDtrack -> GetPID(4) > TRDtrack -> GetPID(1)){
236 track -> SetPDG(kProton);
237 }
238 else if(TRDtrack -> GetPID(3) > TRDtrack -> GetPID(1) && TRDtrack -> GetPID(3) > TRDtrack -> GetPID(2)){
239 track -> SetPDG(kKPlus);
240 }
241 else if(TRDtrack -> GetPID(1) > TRDtrack -> GetPID(2)){
242 track -> SetPDG(kMuonPlus);
243 }
244 else{
245 track -> SetPDG(kPiPlus);
c7cf2032 246 }
247 }
248
249
250 // calculate the probabilities for electron probability using 2-dim LQ, the deposited charge per chamber and the pulse height spectra and fill histograms
773d3f8c 251 fReconstructor -> SetOption("!nn");
252 TRDtrack -> CookPID();
c4c5bbfb 253
254 if(fDebugLevel>=4) Printf("PIDmethod[%d] Slices[%d] PDG[%d] LQLike[%f]", fReconstructor->GetPIDMethod(), fReconstructor->GetNdEdxSlices(), track->GetPDG(), TRDtrack -> GetPID(0));
255
256
3afc1da3 257// Bool_t bChange = kFALSE;
258
5d6dc395 259 Float_t SumdEdx[AliTRDgeometry::kNlayer];
b718144c 260 Int_t iNClus[AliTRDgeometry::kNlayer];
261 Float_t dEdxSlice[AliTRDgeometry::kNlayer][AliTRDReconstructor::kLQslices];
3afc1da3 262
5d6dc395 263 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
1a6d7c9a 264 TRDtracklet[iChamb] = TRDtrack -> GetTracklet(iChamb);
3afc1da3 265
266// if(!(TRDtracklet[iChamb] -> GetNChange() == 0))
267// bChange = 1;
268
1a6d7c9a 269 SumdEdx[iChamb] = 0.;
3afc1da3 270 iNClus[iChamb] = 0;
271 iNClus[iChamb] = TRDtracklet[iChamb] -> GetN();
272// Printf("NClus[%d]", iNClus[iChamb]);
1a6d7c9a 273 fdEdx = TRDtracklet[iChamb] -> GetdEdx();
274 SumdEdx[iChamb] += fdEdx[0] + fdEdx[1] + fdEdx[2];
b718144c 275 for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
276 dEdxSlice[iChamb][iSlice] = fdEdx[iSlice];
277 }
1a6d7c9a 278 }
279
3afc1da3 280// if(bChange == kTRUE)
281// continue;
282
773d3f8c 283 switch(track->GetPDG()){
284 case kElectron:
285 case kPositron:
b718144c 286 hPIDLQ -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
5d6dc395 287 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
422a2dc0 288 hdEdx -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
b718144c 289 for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
290 hdEdxSlice -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
291 dEdxSlice[iChamb][iSlice]);
292 }
3afc1da3 293 hNClus -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
5d6dc395 294 for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
295 if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
296 continue;
422a2dc0 297 hPH -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
5d6dc395 298 }
1a6d7c9a 299 }
773d3f8c 300 break;
301 case kMuonPlus:
302 case kMuonMinus:
b718144c 303 hPIDLQ -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
5d6dc395 304 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
422a2dc0 305 hdEdx -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
b718144c 306 for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
307 hdEdxSlice -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
308 dEdxSlice[iChamb][iSlice]);
309 }
3afc1da3 310 hNClus -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
5d6dc395 311 for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
312 if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
313 continue;
422a2dc0 314 hPH -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
5d6dc395 315 }
1a6d7c9a 316 }
773d3f8c 317 break;
318 case kPiPlus:
319 case kPiMinus:
b718144c 320 hPIDLQ -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
5d6dc395 321 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
422a2dc0 322 hdEdx -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
b718144c 323 for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
324 hdEdxSlice -> Fill(AliPID::kPion * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
325 dEdxSlice[iChamb][iSlice]);
326 }
3afc1da3 327 hNClus -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
5d6dc395 328 for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
329 if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
330 continue;
422a2dc0 331 hPH -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
5d6dc395 332 }
1a6d7c9a 333 }
773d3f8c 334 break;
335 case kKPlus:
336 case kKMinus:
b718144c 337 hPIDLQ -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
5d6dc395 338 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
422a2dc0 339 hdEdx -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
b718144c 340 for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
341 hdEdxSlice -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
342 dEdxSlice[iChamb][iSlice]);
343 }
3afc1da3 344 hNClus -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
5d6dc395 345 for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
346 if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
347 continue;
422a2dc0 348 hPH -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
5d6dc395 349 }
1a6d7c9a 350 }
773d3f8c 351 break;
352 case kProton:
353 case kProtonBar:
b718144c 354 hPIDLQ -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
5d6dc395 355 for(Int_t iChamb = 0; iChamb < AliTRDgeometry::kNlayer; iChamb++){
422a2dc0 356 hdEdx -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, SumdEdx[iChamb]);
b718144c 357 for(Int_t iSlice = 0; iSlice < AliTRDReconstructor::kLQslices; iSlice++){
358 hdEdxSlice -> Fill(AliPID::kProton * AliTRDCalPID::kNMom * AliTRDReconstructor::kLQslices + iMomBin * AliTRDReconstructor::kLQslices + iSlice,
359 dEdxSlice[iChamb][iSlice]);
360 }
3afc1da3 361 hNClus -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, iNClus[iChamb]);
5d6dc395 362 for(Int_t iClus = 0; iClus < AliTRDtrackerV1::GetNTimeBins(); iClus++){
363 if(!(TRDcluster = (AliTRDcluster*)TRDtracklet[iChamb] -> GetClusters(iClus)))
364 continue;
422a2dc0 365 hPH -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, TRDcluster -> GetLocalTimeBin(), TRDtracklet[iChamb] -> GetdQdl(iClus));
5d6dc395 366 }
1a6d7c9a 367 }
773d3f8c 368 break;
369 }
370
371
372 // calculate the probabilities and fill histograms for electrons using NN
373 fReconstructor -> SetOption("nn");
374 TRDtrack->CookPID();
c4c5bbfb 375
376
377 if(fDebugLevel>=4) Printf("PIDmethod[%d] Slices[%d] PDG[%d] NNLike[%f]", fReconstructor->GetPIDMethod(), fReconstructor->GetNdEdxSlices(), track->GetPDG(), TRDtrack -> GetPID(0));
378
379
773d3f8c 380 switch(track->GetPDG()){
381 case kElectron:
382 case kPositron:
b718144c 383 hPIDNN -> Fill(AliPID::kElectron * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
773d3f8c 384 break;
385 case kMuonPlus:
386 case kMuonMinus:
b718144c 387 hPIDNN -> Fill(AliPID::kMuon * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
773d3f8c 388 break;
389 case kPiPlus:
390 case kPiMinus:
b718144c 391 hPIDNN -> Fill(AliPID::kPion * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
773d3f8c 392 break;
393 case kKPlus:
394 case kKMinus:
b718144c 395 hPIDNN -> Fill(AliPID::kKaon * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
773d3f8c 396 break;
397 case kProton:
398 case kProtonBar:
b718144c 399 hPIDNN -> Fill(AliPID::kProton * AliTRDCalPID::kNMom + iMomBin, TRDtrack -> GetPID(AliPID::kElectron));
773d3f8c 400 break;
401 }
422a2dc0 402
773d3f8c 403 }
404
422a2dc0 405 util -> Delete();
3d86166d 406 PostData(0, fContainer);
773d3f8c 407}
408
a391a274 409
c4c5bbfb 410//________________________________________________________
a391a274 411void AliTRDpidChecker::GetRefFigure(Int_t ifig)
c4c5bbfb 412{
a391a274 413 Bool_t FIRST = kTRUE;
c026e9cc 414 TGraphErrors *g = 0x0;
a391a274 415 TH1 *h1 = 0x0;
416 TH2 *h2 = 0x0;
c4c5bbfb 417 switch(ifig){
418 case 0:
c026e9cc 419 g = (TGraphErrors*)fContainer->At(kGraphStart);
420 g->Draw("apl");
421 g->GetHistogram()->GetXaxis()->SetTitle("p [GeV/c]");
422 g->GetHistogram()->GetXaxis()->SetRangeUser(.6, 10.5);
423 g->GetHistogram()->GetYaxis()->SetTitle("#epsilon_{#pi} [%]");
a391a274 424 ((TGraphErrors*)fContainer->At(kGraphStart+1))->Draw("pl");
425 gPad->SetLogy();
c026e9cc 426 gPad->SetLogx();
3afc1da3 427 gPad->SetGridy();
428 gPad->SetGridx();
c4c5bbfb 429 break;
a391a274 430 case 1:
431 // save 2.0 GeV projection as reference
432 FIRST = kTRUE;
433 h2 = (TH2F*)(fContainer->At(kdEdx));
434 for(Int_t is = AliPID::kSPECIES-1; is>=0; is--){
435 Int_t bin = is*AliTRDCalPID::kNMom+4;
436 h1 = h2->ProjectionY("px", bin, bin);
437 if(!h1->GetEntries()) continue;
438 h1->Scale(1./h1->Integral());
439 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
440 h1->DrawClone(FIRST ? "c" : "samec");
441 FIRST = kFALSE;
442 }
443 gPad->SetLogy();
c026e9cc 444 gPad->SetLogx(0);
3afc1da3 445 gPad->SetGridy();
446 gPad->SetGridx();
a391a274 447 break;
448 case 2:
449 // save 2.0 GeV projection as reference
450 FIRST = kTRUE;
451 h2 = (TH2F*)(fContainer->At(kPH));
452 for(Int_t is=0; is<AliPID::kSPECIES; is++){
453 Int_t bin = is*AliTRDCalPID::kNMom+4;
454 h1 = h2->ProjectionY("py", bin, bin);
455 if(!h1->GetEntries()) continue;
456 h1->SetMarkerStyle(24);
457 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
458 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
459 h1->DrawClone(FIRST ? "e2" : "same e2");
460 FIRST = kFALSE;
461 }
462 gPad->SetLogy(0);
c026e9cc 463 gPad->SetLogx(0);
3afc1da3 464 gPad->SetGridy();
465 gPad->SetGridx();
466 break;
467 case 3:
468 // save 2.0 GeV projection as reference
469 FIRST = kTRUE;
470 h2 = (TH2F*)(fContainer->At(kNClus));
471 for(Int_t is=0; is<AliPID::kSPECIES; is++){
472 Int_t bin = is*AliTRDCalPID::kNMom+4;
473 h1 = h2->ProjectionY("py", bin, bin);
474 if(!h1->GetEntries()) continue;
475 h1->SetMarkerStyle(24);
476 h1->SetMarkerColor(AliTRDCalPID::GetPartColor(is));
477 h1->SetLineColor(AliTRDCalPID::GetPartColor(is));
478 h1->DrawClone(FIRST ? "e2" : "same e2");
479 FIRST = kFALSE;
480 }
481 gPad->SetLogy(0);
482 gPad->SetLogx(0);
483 gPad->SetGridy();
484 gPad->SetGridx();
c4c5bbfb 485 break;
486 }
487}
488
773d3f8c 489//________________________________________________________________________
d85cd79c 490Bool_t AliTRDpidChecker::PostProcess()
773d3f8c 491{
c7cf2032 492 // Draw result to the screen
493 // Called once at the end of the query
494
3d86166d 495 if (!fContainer) {
773d3f8c 496 Printf("ERROR: list not available");
d85cd79c 497 return kFALSE;
773d3f8c 498 }
c7cf2032 499// return kTRUE; // testing protection
773d3f8c 500
1a6d7c9a 501
773d3f8c 502 // container for the pion efficiencies and the errors
503 Double_t PionEffiLQ[AliTRDCalPID::kNMom],
773d3f8c 504 PionEffiErrorLQ[AliTRDCalPID::kNMom],
422a2dc0 505 EleEffiLQ[AliTRDCalPID::kNMom],
506 ThresholdLQ[AliTRDCalPID::kNMom];
507
508 Double_t PionEffiNN[AliTRDCalPID::kNMom],
509 PionEffiErrorNN[AliTRDCalPID::kNMom],
510 EleEffiNN[AliTRDCalPID::kNMom],
511 ThresholdNN[AliTRDCalPID::kNMom];
773d3f8c 512
422a2dc0 513 Float_t mom = 0.;
514
515 TH1D *Histo1=0x0, *Histo2=0x0;
516
517 TH2F *hPIDLQ=0x0, *hPIDNN=0x0;
518 hPIDLQ = (TH2F*)fContainer->At(kLQlikelihood);
519 hPIDNN = (TH2F*)fContainer->At(kNNlikelihood);
773d3f8c 520
521 // calculate the pion efficiencies and the errors for 90% electron efficiency (2-dim LQ)
522 for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){
422a2dc0 523
524 AliTRDpidUtil *util = new AliTRDpidUtil();
525 mom = AliTRDCalPID::GetMomentum(iMom);
526
527 Histo1 = hPIDLQ -> ProjectionY("LQ_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1);
528 Histo2 = hPIDLQ -> ProjectionY("LQ_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1);
529
530 util -> CalculatePionEffi(Histo1, Histo2);
531
532 PionEffiLQ[iMom] = util -> GetPionEfficiency();
533 PionEffiErrorLQ[iMom] = util -> GetError();
534 EleEffiLQ[iMom] = util -> GetCalcElectronEfficiency();
535 ThresholdLQ[iMom] = util -> GetThreshold();
536
c7cf2032 537 if(fDebugLevel>=1) Printf("Pion Efficiency for 2-dim LQ is : %f +/- %f\n\n", PionEffiLQ[iMom], PionEffiErrorLQ[iMom]);
422a2dc0 538 util -> Delete();
773d3f8c 539 }
540
422a2dc0 541
542
773d3f8c 543 // calculate the pion efficiencies and the errors for 90% electron efficiency (NN)
544 for(Int_t iMom = 0; iMom < AliTRDCalPID::kNMom; iMom++){
422a2dc0 545
546 AliTRDpidUtil *util = new AliTRDpidUtil();
547 mom = AliTRDCalPID::GetMomentum(iMom);
548
549 Histo1 = hPIDNN -> ProjectionY("NN_ele",AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1,AliTRDCalPID::kNMom*AliPID::kElectron+iMom+1);
550 Histo2 = hPIDNN -> ProjectionY("NN_pio",AliTRDCalPID::kNMom*AliPID::kPion+iMom+1,AliTRDCalPID::kNMom*AliPID::kPion+iMom+1);
551
552 util -> CalculatePionEffi(Histo1, Histo2);
553
554 PionEffiNN[iMom] = util -> GetPionEfficiency();
555 PionEffiErrorNN[iMom] = util -> GetError();
556 EleEffiNN[iMom] = util -> GetCalcElectronEfficiency();
557 ThresholdNN[iMom] = util -> GetThreshold();
558
c7cf2032 559 if(fDebugLevel>=1) Printf("Pion Efficiency for NN is : %f +/- %f\n\n", PionEffiNN[iMom], PionEffiErrorNN[iMom]);
422a2dc0 560
561 util -> Delete();
773d3f8c 562 }
563
564
565 // create TGraph to plot the pion efficiencies
422a2dc0 566 TGraphErrors *gEffisLQ=0x0, *gEffisNN=0x0;
567 gEffisLQ = (TGraphErrors*)fContainer->At(kGraphLQ);
568 gEffisNN = (TGraphErrors*)fContainer->At(kGraphNN);
773d3f8c 569
773d3f8c 570
571 for(Int_t iBin = 0; iBin < AliTRDCalPID::kNMom; iBin++){
422a2dc0 572
773d3f8c 573 Float_t momentum = AliTRDCalPID::GetMomentum(iBin);
574 gEffisLQ->SetPoint(iBin, momentum, PionEffiLQ[iBin]);
422a2dc0 575 gEffisLQ->SetPointError(iBin, 0., PionEffiErrorLQ[iBin]);
773d3f8c 576
577 gEffisNN->SetPoint(iBin, momentum, PionEffiNN[iBin]);
422a2dc0 578 gEffisNN->SetPointError(iBin, 0., PionEffiErrorNN[iBin]);
773d3f8c 579 }
580
c7cf2032 581
3afc1da3 582 fNRefFigures = 4/*1*/;
c7cf2032 583 return kTRUE; // testing protection
d85cd79c 584}
585
c7cf2032 586
d85cd79c 587//________________________________________________________________________
588void AliTRDpidChecker::Terminate(Option_t *)
589{
590 // Draw result to the screen
591 // Called once at the end of the query
592
593 fContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
594 if (!fContainer) {
595 Printf("ERROR: list not available");
596 return;
597 }
773d3f8c 598}
599
600