]>
Commit | Line | Data |
---|---|---|
70da6c5a | 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 | // Task for Efficiency studies | |
17 | // Used for testing classes AliHFEcontainer and AliHFEfilter | |
18 | // Creates Efficiency Histograms | |
19 | // | |
20 | // Author: | |
21 | // Markus Fasel <M.Fasel@gsi.de> | |
22 | // | |
23 | #include <TAxis.h> | |
24 | #include <TCanvas.h> | |
25 | #include <TFile.h> | |
26 | #include <TGraphErrors.h> | |
27 | #include <TH1F.h> | |
28 | #include <TH2.h> | |
29 | #include <TIterator.h> | |
30 | #include <TLegend.h> | |
31 | #include <TObjArray.h> | |
32 | #include <TPad.h> | |
33 | ||
91c7e1ec | 34 | #include "AliAnalysisManager.h" |
70da6c5a | 35 | #include "AliCFAcceptanceCuts.h" |
91c7e1ec | 36 | #include "AliCFTrackIsPrimaryCuts.h" |
70da6c5a | 37 | #include "AliCFContainer.h" |
38 | #include "AliCFEffGrid.h" | |
39 | #include "AliESDEvent.h" | |
40 | #include "AliHFEcollection.h" | |
41 | #include "AliHFEcontainer.h" | |
42 | #include "AliHFEcutStep.h" | |
43 | #include "AliHFEefficiency.h" | |
44 | #include "AliHFEextraCuts.h" | |
45 | #include "AliHFEtrackFilter.h" | |
46 | #include "AliHFEtools.h" | |
47 | #include "AliMCEvent.h" | |
91c7e1ec | 48 | #include "AliMCEventHandler.h" |
70da6c5a | 49 | |
50 | ClassImp(AliHFEefficiency) | |
51 | ||
52 | AliHFEefficiency::AliHFEefficiency(): | |
53 | AliAnalysisTaskSE(), | |
54 | fFilter(NULL), | |
55 | fMCcut(NULL), | |
56 | fAcceptanceCuts(NULL), | |
57 | fEfficiency(NULL), | |
58 | fOutput(NULL), | |
59 | fCutTRD(kFALSE) | |
60 | { | |
61 | // | |
62 | // Default constructor | |
63 | // | |
64 | } | |
65 | ||
66 | AliHFEefficiency::AliHFEefficiency(const Char_t *name): | |
67 | AliAnalysisTaskSE(name), | |
68 | fFilter(NULL), | |
69 | fMCcut(NULL), | |
70 | fAcceptanceCuts(NULL), | |
71 | fEfficiency(NULL), | |
72 | fOutput(NULL), | |
73 | fCutTRD(kFALSE) | |
74 | { | |
75 | // | |
76 | // Main Constructor | |
77 | // | |
78 | DefineOutput(1, AliHFEcontainer::Class()); | |
79 | DefineOutput(2, TList::Class()); | |
faee3b18 | 80 | |
81 | SetRunTerminate(); | |
70da6c5a | 82 | } |
83 | ||
84 | AliHFEefficiency::~AliHFEefficiency(){ | |
85 | // | |
86 | // Destructor | |
87 | // | |
88 | if(fFilter) delete fFilter; | |
89 | if(fEfficiency) delete fEfficiency; | |
90 | if(fAcceptanceCuts) delete fAcceptanceCuts; | |
91 | if(fOutput) delete fOutput; | |
92 | } | |
93 | ||
94 | void AliHFEefficiency::UserCreateOutputObjects(){ | |
95 | // | |
96 | // Create the output objects | |
97 | // | |
98 | fEfficiency = new AliHFEcontainer("Efficiency"); | |
99 | fEfficiency->SetNumberOfVariables(4); | |
100 | // InitializeVariables | |
101 | fEfficiency->SetVariableName(0, "pt"); | |
102 | fEfficiency->MakeLogarithmicBinning(0, 40, 0.1, 10); | |
103 | fEfficiency->SetVariableName(1, "Eta"); | |
104 | fEfficiency->MakeLinearBinning(1, 20, -0.9, 0.9); | |
105 | fEfficiency->SetVariableName(2, "phi"); | |
106 | fEfficiency->MakeLinearBinning(2 , 18, 0, 2 * TMath::Pi()); | |
107 | fEfficiency->SetVariableName(3, "Charge"); | |
108 | fEfficiency->MakeLinearBinning(3, 2, -1.1, 1.1); | |
109 | // New container | |
110 | fEfficiency->CreateContainer("MCFilter", "Efficiency for MC Studies", 2); | |
111 | AliCFContainer *cont = fEfficiency->GetCFContainer("MCFilter"); | |
112 | cont->SetStepTitle(0, "MC Signal"); | |
113 | cont->SetStepTitle(1, "MC Signal in acceptance"); | |
114 | ||
115 | fFilter = new AliHFEtrackFilter("testfilter"); | |
116 | fFilter->GenerateCutSteps(); | |
117 | fFilter->InitCF(fEfficiency); | |
118 | fMCcut = fFilter->GetMCSignalCuts(); | |
119 | if(fCutTRD){ | |
120 | AliHFEcutStep *hfeTRD = fFilter->GetCutStep("HFETRD"); | |
121 | AliHFEextraCuts *hfecuts = dynamic_cast<AliHFEextraCuts *>(hfeTRD->GetCut("HFETRDCuts")); | |
bf892a6a | 122 | if(hfecuts) hfecuts->SetMinTrackletsTRD(4); |
70da6c5a | 123 | } |
124 | AliHFEcutStep *hfeITS = fFilter->GetCutStep("HFEITS"); | |
125 | AliHFEextraCuts *hfeitscuts = dynamic_cast<AliHFEextraCuts *>(hfeITS->GetCut("HFEPixelsCuts")); | |
bf892a6a | 126 | if(hfeitscuts)hfeitscuts->SetRequireITSpixel(AliHFEextraCuts::kFirst); |
70da6c5a | 127 | |
91c7e1ec | 128 | AliHFEcutStep *primary = fFilter->GetCutStep("Primary"); |
129 | AliCFTrackIsPrimaryCuts *primcuts = dynamic_cast<AliCFTrackIsPrimaryCuts *>(primary->GetCut("PrimaryCuts")); | |
bf892a6a | 130 | if(primcuts){ |
131 | primcuts->SetMaxDCAToVertexXY(3); | |
132 | primcuts->SetMaxDCAToVertexZ(5); | |
133 | } | |
91c7e1ec | 134 | |
70da6c5a | 135 | fAcceptanceCuts = new AliCFAcceptanceCuts("Acceptance", "MC Acceptance Cuts"); |
136 | fAcceptanceCuts->SetMinNHitITS(3); | |
137 | fAcceptanceCuts->SetMinNHitTPC(2); | |
138 | fAcceptanceCuts->SetMinNHitTRD(0); | |
139 | ||
140 | // create additional histos for testing purpose | |
141 | fOutput = new AliHFEcollection("histos", "QA histograms"); | |
142 | fOutput->CreateTH1F("hNtracks","Number of Tracks per Event", 100, 0, 100); | |
faee3b18 | 143 | fOutput->CreateTH1F("hPt","Pt of the tracks", 40, 0.1, 10, 0); |
70da6c5a | 144 | fOutput->CreateTH1F("kinkIndex", "Kink Index", 400, -200, 200); |
145 | fOutput->CreateTH1F("itspixel", "ITS PIXEL", 2, 0, 2); | |
146 | fOutput->CreateTH1F("mcmother", "Mother PDG", 1000, -500, 500); | |
faee3b18 | 147 | fOutput->CreateTH2F("ptres", "Pt Resolution", 40, 0.1, 10, 200, -1.5, 1.5, 0); |
70da6c5a | 148 | } |
149 | ||
150 | void AliHFEefficiency::UserExec(Option_t *){ | |
151 | // | |
152 | // Event Loop | |
153 | // Filter track, fill Efficiency container | |
154 | // | |
e3ae862b | 155 | AliESDEvent *esdevent = dynamic_cast<AliESDEvent *>(fInputEvent); |
156 | if(!esdevent){ | |
157 | AliError("ESD Event required"); | |
158 | return; | |
159 | } | |
70da6c5a | 160 | fEfficiency->NewEvent(); |
161 | fFilter->SetRecEvent(fInputEvent); | |
162 | if(fMCEvent){ | |
91c7e1ec | 163 | AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()); |
215ffe88 | 164 | if ( ! mcH ) { |
165 | AliError("Cannot get MC truth event handler"); | |
166 | return; | |
167 | } | |
e3ae862b | 168 | if(mcH &&(!mcH->InitOk())) return; |
169 | if(mcH &&(!mcH->TreeK())) return; | |
170 | if(mcH &&(!mcH->TreeTR())) return; | |
70da6c5a | 171 | fFilter->SetMC(fMCEvent); |
172 | FilterMC(); | |
173 | } | |
e3ae862b | 174 | fFilter->FilterTracks(esdevent); |
70da6c5a | 175 | TObjArray *tracks = fFilter->GetFilteredTracks(); |
176 | TIterator *iter = tracks->MakeIterator(); | |
177 | fOutput->Fill("hNtracks", tracks->GetEntriesFast()); | |
178 | AliESDtrack *track = NULL; | |
179 | while((track = dynamic_cast<AliESDtrack *>(iter->Next()))){ | |
180 | fOutput->Fill("hPt", track->Pt()); | |
181 | fOutput->Fill("kinkIndex", track->GetKinkIndex(0)); | |
182 | if(TESTBIT(track->GetITSClusterMap(), 0)) | |
183 | fOutput->Fill("itspixel",1); | |
184 | else | |
185 | fOutput->Fill("itspixel",0); | |
faee3b18 | 186 | if(fMCEvent){ |
187 | AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(track->GetLabel()))); | |
188 | if(mctrack){ | |
189 | Int_t motherLabel = mctrack->Particle()->GetFirstMother(); | |
190 | if(motherLabel){ | |
191 | AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(motherLabel)); | |
192 | if(mother)fOutput->Fill("mcmother", mother->Particle()->GetPdgCode()); | |
193 | } | |
194 | fOutput->Fill("ptres", mctrack->Pt(), (track->Pt() - mctrack->Pt())/mctrack->Pt()); | |
91c7e1ec | 195 | } |
91c7e1ec | 196 | } |
70da6c5a | 197 | } |
198 | delete iter; | |
199 | fFilter->Flush(); | |
200 | PostData(1, fEfficiency); | |
201 | PostData(2, fOutput->GetList()); | |
202 | } | |
203 | ||
204 | void AliHFEefficiency::Terminate(Option_t *){ | |
205 | // | |
206 | // Evaluate Results | |
207 | // | |
208 | fEfficiency = dynamic_cast<AliHFEcontainer *>(GetOutputData(1)); | |
faee3b18 | 209 | if(!fEfficiency){ |
210 | AliError("No Output data available"); | |
211 | return; | |
212 | } | |
213 | ||
214 | if(!IsRunningTerminate()) return; | |
70da6c5a | 215 | PostProcess(); |
216 | ||
217 | TList *l = dynamic_cast<TList *>(GetOutputData(2)); | |
bf892a6a | 218 | if(l) DrawPtResolution(l); |
70da6c5a | 219 | } |
220 | ||
221 | void AliHFEefficiency::FilterMC(){ | |
222 | // | |
223 | // Cut MC tracks | |
224 | // | |
225 | Double_t cont[4] = {0., 0., 0., 0.}; | |
226 | AliMCParticle *track = NULL; | |
227 | for(Int_t itrack = 0; itrack < fMCEvent->GetNumberOfTracks(); itrack++){ | |
228 | track = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(itrack)); | |
e3ae862b | 229 | if(!track) continue; |
70da6c5a | 230 | if(!fMCcut->IsSelected(track)) continue; |
231 | cont[0] = track->Pt(); | |
232 | cont[1] = track->Eta(); | |
233 | cont[2] = track->Phi(); | |
234 | cont[3] = track->Charge()/3; | |
235 | //AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(track->Particle()->GetFirstMother())); | |
236 | //if(TMath::Abs(mother->Particle()->GetPdgCode()) != 443) continue; | |
237 | fEfficiency->FillCFContainer("MCFilter", 0, cont); | |
238 | if(!fAcceptanceCuts->IsSelected(track)) continue; | |
239 | fEfficiency->FillCFContainer("MCFilter", 1, cont); | |
240 | } | |
241 | } | |
242 | ||
243 | void AliHFEefficiency::Load(const char* filename){ | |
3a72645a | 244 | // |
245 | // Load results for post processing | |
246 | // | |
70da6c5a | 247 | TFile *input = TFile::Open(filename); |
3ccf8f4c | 248 | AliHFEcontainer *ctmp = dynamic_cast<AliHFEcontainer *>(input->Get("Efficiency")); |
249 | if(ctmp) fEfficiency = dynamic_cast<AliHFEcontainer *>(ctmp->Clone()); | |
70da6c5a | 250 | input->Close(); |
251 | delete input; | |
252 | } | |
253 | ||
254 | void AliHFEefficiency::PostProcess(){ | |
255 | // | |
256 | // Calculate Efficiencies | |
257 | // | |
258 | fEfficiency->Print(); | |
259 | ||
260 | // Prepare containers | |
261 | AliCFContainer *crec = fEfficiency->MakeMergedCFContainer("merged", "Merged Container", "MCFilter:testfilter_container_signal"); | |
262 | AliCFContainer *cmc = fEfficiency->MakeMergedCFContainer("merged", "Merged Container", "MCFilter:testfilter_container_signalMC"); | |
263 | AliCFEffGrid *gridRec = new AliCFEffGrid("effrec", "Efficiency Calculation using rec vars", *crec); | |
264 | AliCFEffGrid *gridMC = new AliCFEffGrid("effmc", "Efficiency Calculation using mc vars", *cmc); | |
265 | ||
266 | TCanvas *ctmp = NULL; | |
267 | for(Int_t ivar = 0; ivar < crec->GetNVar(); ivar++){ | |
268 | ctmp = new TCanvas(Form("cEffSignal%s", crec->GetVarTitle(ivar)), Form("Signal Efficiency (%s)", crec->GetVarTitle(ivar))); | |
269 | ctmp->cd(); | |
270 | DrawSignalEfficiency(gridRec, crec, ivar); | |
271 | ||
272 | ctmp = new TCanvas(Form("cCutEff%s", crec->GetVarTitle(ivar)), Form("Cut Step Efficiency (%s)", crec->GetVarTitle(ivar))); | |
273 | ctmp->cd(); | |
274 | DrawCutEfficiency(gridRec, crec, ivar); | |
275 | ||
276 | ctmp = new TCanvas(Form("cEffSignalMC%s", crec->GetVarTitle(ivar)), Form("Signal Efficiency (%s, MC)", crec->GetVarTitle(ivar))); | |
277 | ctmp->cd(); | |
278 | DrawSignalEfficiency(gridMC, cmc, ivar); | |
279 | ||
280 | ctmp = new TCanvas(Form("cCutEffMC%s", crec->GetVarTitle(ivar)), Form("Cut Step Efficiency (%s, MC)", crec->GetVarTitle(ivar))); | |
281 | ctmp->cd(); | |
282 | DrawCutEfficiency(gridMC, cmc, ivar); | |
283 | } | |
284 | ||
285 | CalculatePTsmearing(); | |
286 | ||
287 | delete gridRec; | |
288 | delete gridMC; | |
289 | delete crec; | |
290 | delete cmc; | |
291 | } | |
292 | ||
293 | void AliHFEefficiency::DrawSignalEfficiency(AliCFEffGrid *eff, AliCFContainer *cont, Int_t var){ | |
294 | // | |
295 | // Efficiency of a cut step with respect to the signal | |
296 | // | |
297 | TH1 *hEff = NULL; | |
298 | Bool_t kFirst = kTRUE; | |
299 | TLegend *leg = new TLegend(0.5, 0.7, 0.89, 0.89); | |
300 | for(Int_t istep = 1; istep < cont->GetNStep(); istep++){ | |
301 | eff->CalculateEfficiency(istep, 0); | |
302 | hEff = eff->Project(var); | |
303 | hEff->SetTitle(Form("Signal Efficiency (%s)", cont->GetVarTitle(var))); | |
304 | hEff->GetYaxis()->SetRangeUser(0., 1.6); | |
305 | hEff->GetYaxis()->SetTitle("Efficiency"); | |
306 | hEff->SetStats(kFALSE); | |
307 | hEff->SetMarkerStyle(20+istep); | |
308 | hEff->SetMarkerColor(kBlue - 5 + istep); | |
309 | hEff->Draw(kFirst ? "ep" : "epsame"); | |
310 | kFirst = kFALSE; | |
311 | leg->AddEntry(hEff, cont->GetStepTitle(istep), "p"); | |
312 | } | |
313 | leg->SetBorderSize(0); | |
314 | leg->SetFillStyle(0); | |
315 | leg->Draw(); | |
316 | gPad->Update(); | |
317 | } | |
318 | ||
319 | void AliHFEefficiency::DrawCutEfficiency(AliCFEffGrid *eff, AliCFContainer *cont, Int_t var){ | |
320 | // | |
321 | // Efficiency of a cut step with respect to the previous one | |
322 | // | |
323 | Bool_t kFirst = kTRUE; | |
324 | TH1 *hEff = NULL; | |
325 | TLegend *leg = new TLegend(0.5, 0.7, 0.89, 0.89); | |
326 | TString effTitle; | |
327 | Int_t start = 0, length = 0; | |
328 | for(Int_t istep = 1; istep < cont->GetNStep(); istep++){ | |
329 | eff->CalculateEfficiency(istep, istep -1); | |
330 | hEff = eff->Project(var); | |
331 | effTitle = hEff->GetTitle(); | |
332 | start = effTitle.Index("projected"); length = effTitle.Length() - start; | |
333 | effTitle.Remove(start - 1, length + 1); | |
334 | effTitle.ReplaceAll("Efficiency: ", ""); | |
335 | hEff->SetTitle(Form("Cut Efficiency (%s)", cont->GetVarTitle(var))); | |
336 | hEff->GetYaxis()->SetRangeUser(0., 1.6); | |
337 | hEff->GetYaxis()->SetTitle("Efficiency"); | |
338 | hEff->SetStats(kFALSE); | |
339 | hEff->SetMarkerStyle(20+istep); | |
340 | hEff->SetMarkerColor(kBlue - 5 + istep); | |
341 | hEff->Draw(kFirst ? "ep" : "epsame"); | |
342 | kFirst = kFALSE; | |
d44e3c84 | 343 | leg->AddEntry(hEff, cont->GetStepTitle(istep), "p"); |
70da6c5a | 344 | } |
345 | leg->SetBorderSize(0); | |
346 | leg->SetFillStyle(0); | |
347 | leg->Draw(); | |
348 | gPad->Update(); | |
349 | } | |
350 | ||
351 | void AliHFEefficiency::CalculatePTsmearing(){ | |
352 | // | |
353 | // Efficiency of a cut step with respect to the same cut step filled with MC | |
354 | // Information | |
355 | // | |
356 | AliCFContainer *cont = fEfficiency->MakeMergedCFContainer("merged", "Merged Container", "testfilter_container_signal:testfilter_container_signalMC"); | |
357 | AliCFEffGrid *grid = new AliCFEffGrid("effc", "Efficiency Calculation", *cont); | |
358 | Bool_t kFirst = kTRUE; | |
359 | TH1 *hEff = NULL; | |
360 | TCanvas *c2 = new TCanvas("cSmear", "pT modification"); | |
361 | c2->cd(); | |
362 | TLegend *leg1 = new TLegend(0.5, 0.7, 0.89, 0.89); | |
363 | Int_t nStep = cont->GetNStep()/2; | |
364 | for(Int_t istep = 0; istep < nStep; istep++){ | |
365 | grid->CalculateEfficiency(istep, istep + nStep); | |
366 | hEff = grid->Project(0); | |
367 | hEff->SetStats(kFALSE); | |
368 | hEff->SetMarkerStyle(20+istep); | |
369 | hEff->SetMarkerColor(kBlue - 5 + istep); | |
370 | hEff->Draw(kFirst ? "ep" : "epsame"); | |
371 | kFirst = kFALSE; | |
372 | leg1->AddEntry(hEff, cont->GetStepTitle(istep), "p"); | |
373 | } | |
374 | leg1->Draw(); | |
375 | gPad->Update(); | |
376 | ||
377 | delete cont; | |
378 | delete grid; | |
379 | } | |
380 | ||
3a72645a | 381 | void AliHFEefficiency::DrawPtResolution(const TList * const l){ |
70da6c5a | 382 | // |
383 | // Draw pt resolution | |
384 | // | |
385 | TH2 *h2 = dynamic_cast<TH2 *>(l->FindObject("ptres")); | |
bf892a6a | 386 | if(!h2) return; |
70da6c5a | 387 | TGraphErrors *mean = new TGraphErrors(h2->GetNbinsX()); |
388 | TGraphErrors *rms = new TGraphErrors(h2->GetNbinsX()); | |
389 | ||
390 | Double_t pt = 0., pterr = 0., mymean = 0., myrms = 0., mymeanerr = 0., myrmserr = 0.; | |
391 | TH1 *htmp = NULL; | |
392 | for(Int_t imom = 1; imom <= h2->GetXaxis()->GetLast(); imom++){ | |
393 | pt = h2->GetXaxis()->GetBinCenter(imom); | |
394 | pterr = h2->GetXaxis()->GetBinWidth(imom)/2.; | |
395 | htmp = h2->ProjectionY("py", imom, imom); | |
396 | mymean = htmp->GetMean(); | |
397 | myrms = htmp->GetRMS(); | |
398 | mymeanerr = htmp->GetMeanError(); | |
399 | myrmserr = htmp->GetRMSError(); | |
400 | delete htmp; | |
401 | mean->SetPoint(imom -1, pt, mymean); | |
402 | rms->SetPoint(imom -1, pt, myrms); | |
403 | mean->SetPointError(imom-1,pterr,mymeanerr); | |
404 | rms->SetPointError(imom-1,pterr,myrmserr); | |
405 | } | |
406 | ||
407 | TCanvas *c1 = new TCanvas("cPtShift", "pT Shift"); | |
408 | c1->cd(); | |
409 | mean->SetMarkerStyle(22); | |
410 | mean->SetMarkerColor(kBlue); | |
411 | mean->Draw("ape"); | |
412 | rms->SetMarkerStyle(22); | |
413 | rms->SetMarkerColor(kBlack); | |
414 | rms->Draw("pesame"); | |
415 | TLegend *leg = new TLegend(0.5, 0.7, 0.89, 0.89); | |
416 | leg->AddEntry(mean, "Mean", "lp"); | |
417 | leg->AddEntry(rms, "RMS", "lp"); | |
418 | leg->Draw(); | |
419 | gPad->Update(); | |
420 | } | |
421 |