]>
Commit | Line | Data |
---|---|---|
684b93f0 | 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 | //------------------------------ | |
17 | // Analysis task for quality-assurance | |
18 | // of forward detectors ESD | |
19 | // | |
20 | // 12/06/2009 cvetan.cheshkov@cern.ch | |
21 | //------------------------------ | |
22 | ||
23 | #include "TChain.h" | |
24 | #include "TROOT.h" | |
25 | #include "TFile.h" | |
26 | #include "TH1F.h" | |
27 | #include "TH2F.h" | |
28 | #include "TF1.h" | |
29 | #include "TCanvas.h" | |
30 | #include "TVector3.h" | |
31 | #include "TParticle.h" | |
8b04dae1 | 32 | #include "AliVParticle.h" |
33 | #include "AliMCParticle.h" | |
684b93f0 | 34 | #include "AliESDEvent.h" |
35 | #include "AliESDv0.h" | |
36 | #include "AliESDcascade.h" | |
37 | #include "AliESDMuonTrack.h" | |
38 | #include "AliESDCaloCluster.h" | |
39 | #include "AliRun.h" | |
40 | #include "AliMCEvent.h" | |
41 | #include "AliStack.h" | |
42 | #include "AliGenEventHeader.h" | |
43 | #include "AliPID.h" | |
44 | #include "AliESDVZERO.h" | |
45 | ||
46 | #include "AliAnaFwdDetsQA.h" | |
47 | ||
48 | ClassImp(AliAnaFwdDetsQA) | |
49 | ||
50 | AliAnaFwdDetsQA::AliAnaFwdDetsQA(): | |
51 | AliAnalysisTaskSE("AliAnaFwdDetsQA"), | |
52 | fListOfHistos(0), | |
53 | fT0vtxRec(0), | |
54 | fT0vtxRecGen(0), | |
55 | fT0time(0), | |
c543bd05 | 56 | fT0time2(0), |
684b93f0 | 57 | fT0mult(0), |
58 | fT0vtxRes(0), | |
c543bd05 | 59 | fT0ampl(0), |
684b93f0 | 60 | fV0a(0), |
61 | fV0c(0), | |
62 | fV0multA(0), | |
63 | fV0multC(0), | |
64 | fV0multAcorr(0), | |
65 | fV0multCcorr(0), | |
66 | fV0Acorr(0), | |
c543bd05 | 67 | fV0Ccorr(0), |
68 | fV0ampl(0) | |
684b93f0 | 69 | { |
70 | // Default constructor | |
71 | // Define input and output slots here | |
72 | // Input slot #0 works with a TChain | |
73 | DefineInput(0, TChain::Class()); | |
74 | // Output slot #1 TList | |
75 | DefineOutput(1, TList::Class()); | |
76 | } | |
77 | ||
78 | AliAnaFwdDetsQA::AliAnaFwdDetsQA(const char* name): | |
79 | AliAnalysisTaskSE(name), | |
80 | fListOfHistos(0), | |
81 | fT0vtxRec(0), | |
82 | fT0vtxRecGen(0), | |
83 | fT0time(0), | |
c543bd05 | 84 | fT0time2(0), |
684b93f0 | 85 | fT0mult(0), |
86 | fT0vtxRes(0), | |
c543bd05 | 87 | fT0ampl(0), |
684b93f0 | 88 | fV0a(0), |
89 | fV0c(0), | |
90 | fV0multA(0), | |
91 | fV0multC(0), | |
92 | fV0multAcorr(0), | |
93 | fV0multCcorr(0), | |
94 | fV0Acorr(0), | |
c543bd05 | 95 | fV0Ccorr(0), |
96 | fV0ampl(0) | |
684b93f0 | 97 | { |
98 | // Constructor | |
684b93f0 | 99 | // Define input and output slots here |
100 | // Input slot #0 works with a TChain | |
fc6c8a4a | 101 | AliInfo("Constructor AliAnaFwdDetsQA"); |
684b93f0 | 102 | DefineInput(0, TChain::Class()); |
103 | // Output slot #1 TList | |
104 | DefineOutput(1, TList::Class()); | |
105 | } | |
106 | ||
107 | TH1F * AliAnaFwdDetsQA::CreateHisto(const char* name, const char* title,Int_t nBins, | |
108 | Double_t xMin, Double_t xMax, | |
109 | const char* xLabel, const char* yLabel) | |
110 | { | |
fc6c8a4a | 111 | // helper method which can be used |
112 | // in order to create a histogram | |
684b93f0 | 113 | TH1F* result = new TH1F(name, title, nBins, xMin, xMax); |
114 | result->SetOption("E"); | |
115 | if (xLabel) result->GetXaxis()->SetTitle(xLabel); | |
116 | if (yLabel) result->GetYaxis()->SetTitle(yLabel); | |
117 | result->SetMarkerStyle(kFullCircle); | |
118 | return result; | |
119 | } | |
120 | ||
121 | TH1F *AliAnaFwdDetsQA::CreateEffHisto(const TH1F* hGen, const TH1F* hRec) | |
122 | { | |
fc6c8a4a | 123 | // helper method which can be used |
124 | // in order create an efficiency histogram | |
684b93f0 | 125 | Int_t nBins = hGen->GetNbinsX(); |
126 | TH1F* hEff = (TH1F*) hGen->Clone("hEff"); | |
127 | hEff->SetTitle(""); | |
128 | hEff->SetStats(kFALSE); | |
129 | hEff->SetMinimum(0.); | |
130 | hEff->SetMaximum(110.); | |
131 | hEff->GetYaxis()->SetTitle("#epsilon [%]"); | |
132 | ||
133 | for (Int_t iBin = 0; iBin <= nBins; iBin++) { | |
134 | Double_t nGenEff = hGen->GetBinContent(iBin); | |
135 | Double_t nRecEff = hRec->GetBinContent(iBin); | |
136 | if (nGenEff > 0) { | |
137 | Double_t eff = nRecEff/nGenEff; | |
138 | hEff->SetBinContent(iBin, 100. * eff); | |
139 | Double_t error = sqrt(eff*(1.-eff) / nGenEff); | |
fc6c8a4a | 140 | if (error < 1e-12) error = 0.0001; |
684b93f0 | 141 | hEff->SetBinError(iBin, 100. * error); |
142 | } | |
143 | else { | |
144 | hEff->SetBinContent(iBin, -100.); | |
145 | hEff->SetBinError(iBin, 0); | |
146 | } | |
147 | } | |
148 | return hEff; | |
149 | } | |
150 | ||
151 | ||
152 | Bool_t AliAnaFwdDetsQA::FitHisto(TH1* histo, Double_t& res, Double_t& resError) | |
153 | { | |
fc6c8a4a | 154 | // fit a gaussian to |
155 | // a histogram | |
684b93f0 | 156 | static TF1* fitFunc = new TF1("fitFunc", "gaus"); |
157 | fitFunc->SetLineWidth(2); | |
158 | fitFunc->SetFillStyle(0); | |
159 | Double_t maxFitRange = 2; | |
160 | ||
161 | if (histo->Integral() > 50) { | |
162 | Float_t mean = histo->GetMean(); | |
163 | Float_t rms = histo->GetRMS(); | |
164 | fitFunc->SetRange(mean - maxFitRange*rms, mean + maxFitRange*rms); | |
165 | fitFunc->SetParameters(mean, rms); | |
166 | histo->Fit(fitFunc, "QRI0"); | |
167 | histo->GetFunction("fitFunc")->ResetBit(1<<9); | |
168 | res = TMath::Abs(fitFunc->GetParameter(2)); | |
169 | resError = TMath::Abs(fitFunc->GetParError(2)); | |
170 | return kTRUE; | |
171 | } | |
172 | return kFALSE; | |
173 | } | |
174 | ||
175 | void AliAnaFwdDetsQA::UserCreateOutputObjects() | |
176 | { | |
177 | // Create histograms | |
684b93f0 | 178 | // Create output container |
fc6c8a4a | 179 | AliInfo("AliAnaFwdDetsQA::UserCreateOutputObjects"); |
684b93f0 | 180 | fListOfHistos = new TList(); |
181 | ||
182 | fT0vtxRec = CreateHisto("hT0vtxRec", "Z vertex reconstructed with T0", 100, -25, 25, "Z_{vtx} [cm]", ""); | |
183 | fT0time = CreateHisto("hT0time", "Time0 reconstruction with T0", 5000, 10000, 20000, "t_{0} [ps]", ""); | |
c543bd05 | 184 | fT0time2 = CreateHisto("hT0time2", "Time0 reconstruction with T0 (mult > 10)", 5000, 10000, 20000, "t_{0} [ps]", ""); |
684b93f0 | 185 | fT0mult = CreateHisto("hT0mult", "Total reconstructed multiplicity in T0", 100, -0.5, 99.5, "Multiplicity", ""); |
186 | fT0vtxRes = CreateHisto("hT0vtxRes", "T0 Z vertex resolution", 100, -25, 25, "Delta(Z_{vtx}) [cm]", ""); | |
c543bd05 | 187 | fT0ampl = CreateHisto("hT0ampl","T0 multiplicity in single channel (all T0 channels)",400,-0.5,99.5); |
684b93f0 | 188 | |
189 | fT0vtxRecGen = new TH2F("hT0vtxRecGen", "T0 reconstructed vs generated Z vertex", 100, -25, 25, 100, -25, 25); | |
190 | fT0vtxRecGen->GetXaxis()->SetTitle("Generated Z vertex"); | |
191 | fT0vtxRecGen->GetYaxis()->SetTitle("Reconstructed Z vertex"); | |
192 | fT0vtxRecGen->SetMarkerStyle(kFullCircle); | |
193 | fT0vtxRecGen->SetMarkerSize(0.4); | |
194 | ||
195 | fV0a = CreateHisto("hV0a","Number of fired PMTs (V0A)",65,-0.5,64.5); | |
196 | fV0c = CreateHisto("hV0c","Number of fired PMTs (V0C)",65,-0.5,64.5); | |
197 | fV0multA = CreateHisto("hV0multA","Total reconstructed multiplicity (V0A)",100,0.,1000.); | |
198 | fV0multC = CreateHisto("hV0multC","Total reconstructed multiplicity (V0C)",100,0.,1000.); | |
c543bd05 | 199 | fV0ampl = CreateHisto("hV0ampl","V0 multiplicity in single channel (all V0 channels)",400,-0.5,99.5); |
684b93f0 | 200 | |
201 | fV0multAcorr = new TH2F("hV0multAcorr", "Reconstructed vs generated (primaries only) multiplicity (V0A)",100,0.,500.,100,0.,1000.); | |
202 | fV0multAcorr->GetXaxis()->SetTitle("# of primaries in V0A acceptance"); | |
203 | fV0multAcorr->GetYaxis()->SetTitle("Reconstructed mutiplicity in V0A"); | |
204 | fV0multCcorr = new TH2F("hV0multCcorr", "Reconstructed vs generated (primaries only) multiplicity (V0C)",100,0.,500.,100,0.,1000.); | |
205 | fV0multCcorr->GetXaxis()->SetTitle("# of primaries in V0C acceptance"); | |
206 | fV0multCcorr->GetYaxis()->SetTitle("Reconstructed mutiplicity in V0C"); | |
207 | ||
208 | fV0Acorr = new TH2F("hV0Acorr", "Number of fired PMTs vs generated (primaries only) multiplicity (V0A)",100,0.,500.,65,-0.5,64.5); | |
209 | fV0Acorr->GetXaxis()->SetTitle("# of primaries in V0A acceptance"); | |
210 | fV0Acorr->GetYaxis()->SetTitle("Number of fired PMTs in V0A"); | |
211 | fV0Ccorr = new TH2F("hV0Ccorr", "Number of fired PMTs vs generated (primaries only) multiplicity (V0C)",100,0.,500.,65,-0.5,64.5); | |
212 | fV0Ccorr->GetXaxis()->SetTitle("# of primaries in V0C acceptance"); | |
213 | fV0Ccorr->GetYaxis()->SetTitle("Number of fired PMTs in V0C"); | |
214 | ||
215 | fListOfHistos->Add(fT0vtxRec); | |
216 | fListOfHistos->Add(fT0time); | |
217 | fListOfHistos->Add(fT0mult); | |
218 | fListOfHistos->Add(fT0vtxRecGen); | |
219 | fListOfHistos->Add(fT0vtxRes); | |
220 | fListOfHistos->Add(fV0a); | |
221 | fListOfHistos->Add(fV0c); | |
222 | fListOfHistos->Add(fV0multA); | |
223 | fListOfHistos->Add(fV0multC); | |
224 | fListOfHistos->Add(fV0multAcorr); | |
225 | fListOfHistos->Add(fV0multCcorr); | |
226 | fListOfHistos->Add(fV0Acorr); | |
227 | fListOfHistos->Add(fV0Ccorr); | |
c543bd05 | 228 | fListOfHistos->Add(fT0ampl); |
229 | fListOfHistos->Add(fV0ampl); | |
230 | fListOfHistos->Add(fT0time2); | |
684b93f0 | 231 | } |
232 | ||
233 | void AliAnaFwdDetsQA::UserExec(Option_t */*option*/) | |
234 | { | |
fc6c8a4a | 235 | // The analysis code |
236 | // goes here | |
684b93f0 | 237 | AliMCEvent* mcEvent = MCEvent(); |
238 | if (!mcEvent) { | |
239 | Printf("ERROR: Could not retrieve MC event"); | |
240 | return; | |
241 | } | |
242 | ||
243 | //Primary vertex needed | |
244 | TArrayF mcvtx(3); | |
245 | mcEvent->GenEventHeader()->PrimaryVertex(mcvtx); | |
246 | ||
247 | AliStack *stack = mcEvent->Stack(); | |
248 | if (!stack) { | |
249 | Printf("ERROR: Could not retrieve MC stack"); | |
250 | return; | |
251 | } | |
252 | ||
253 | Int_t nV0A = 0; | |
254 | Int_t nV0C = 0; | |
255 | for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfTracks(); iTracks++) {//Check this loop again | |
256 | if (!stack->IsPhysicalPrimary(iTracks)) continue; | |
8b04dae1 | 257 | AliMCParticle* track = (AliMCParticle*)mcEvent->GetTrack(iTracks); |
684b93f0 | 258 | TParticle* particle = track->Particle(); |
259 | if (!particle) continue; | |
c543bd05 | 260 | if (track->Charge() == 0) continue; |
684b93f0 | 261 | Double_t eta = particle->Eta(); |
262 | if (eta > 2.8 && eta < 5.1) { | |
263 | nV0A++; | |
264 | } | |
265 | if (eta > -3.7 && eta < -1.7) { | |
266 | nV0C++; | |
267 | } | |
268 | } | |
269 | ||
270 | // ESD information | |
271 | AliVEvent* event = InputEvent(); | |
272 | if (!event) { | |
273 | Printf("ERROR: Could not retrieve event"); | |
274 | return; | |
275 | } | |
276 | ||
277 | AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event); | |
278 | const AliESDTZERO* esdT0 = esd->GetESDTZERO(); | |
279 | Double_t t0zvtx = esdT0->GetT0zVertex(); | |
280 | Double_t t0time = esdT0->GetT0(); | |
281 | ||
282 | fT0vtxRec->Fill(t0zvtx); | |
283 | fT0time->Fill(t0time); | |
284 | fT0vtxRecGen->Fill(mcvtx[2],t0zvtx); | |
285 | t0zvtx *= -1.0; | |
286 | fT0vtxRes->Fill(t0zvtx - mcvtx[2]); | |
287 | ||
288 | Double_t t0mult = 0; | |
289 | for(Int_t i = 0; i < 24; i++) { | |
290 | t0mult += esdT0->GetT0amplitude()[i]; | |
c543bd05 | 291 | fT0ampl->Fill(esdT0->GetT0amplitude()[i]); |
684b93f0 | 292 | } |
293 | ||
294 | fT0mult->Fill(t0mult); | |
c543bd05 | 295 | if (t0mult > 10) |
296 | fT0time2->Fill(t0time); | |
684b93f0 | 297 | |
298 | AliESDVZERO* esdV0 = esd->GetVZEROData(); | |
299 | fV0a->Fill(esdV0->GetNbPMV0A()); | |
300 | fV0c->Fill(esdV0->GetNbPMV0C()); | |
301 | fV0multA->Fill(esdV0->GetMTotV0A()); | |
302 | fV0multC->Fill(esdV0->GetMTotV0C()); | |
303 | ||
304 | fV0multAcorr->Fill((Float_t)nV0A,esdV0->GetMTotV0A()); | |
305 | fV0multCcorr->Fill((Float_t)nV0C,esdV0->GetMTotV0C()); | |
306 | fV0Acorr->Fill((Float_t)nV0A,(Float_t)esdV0->GetNbPMV0A()); | |
307 | fV0Ccorr->Fill((Float_t)nV0C,(Float_t)esdV0->GetNbPMV0C()); | |
308 | ||
c543bd05 | 309 | for(Int_t i = 0; i < 64; i++) { |
310 | fV0ampl->Fill(esdV0->GetMultiplicity(i)); | |
311 | } | |
684b93f0 | 312 | // Post output data. |
313 | PostData(1, fListOfHistos); | |
314 | } | |
315 | ||
316 | void AliAnaFwdDetsQA::Terminate(Option_t *) | |
317 | { | |
fc6c8a4a | 318 | // Terminate |
319 | // Store the output histos | |
684b93f0 | 320 | fListOfHistos = dynamic_cast<TList*>(GetOutputData(1)); |
321 | if (!fListOfHistos) { | |
322 | Printf("ERROR: fListOfHistos not available"); | |
323 | return; | |
324 | } | |
325 | ||
326 | fT0vtxRec = dynamic_cast<TH1F*>(fListOfHistos->At(0)); | |
327 | fT0time = dynamic_cast<TH1F*>(fListOfHistos->At(1)); | |
328 | fT0mult = dynamic_cast<TH1F*>(fListOfHistos->At(2)); | |
329 | fT0vtxRecGen = dynamic_cast<TH2F*>(fListOfHistos->At(3)); | |
330 | fT0vtxRes = dynamic_cast<TH1F*>(fListOfHistos->At(4)); | |
331 | ||
332 | fV0a = dynamic_cast<TH1F*>(fListOfHistos->At(5)); | |
333 | fV0c = dynamic_cast<TH1F*>(fListOfHistos->At(6)); | |
334 | fV0multA = dynamic_cast<TH1F*>(fListOfHistos->At(7)); | |
335 | fV0multC = dynamic_cast<TH1F*>(fListOfHistos->At(8)); | |
336 | fV0multAcorr = dynamic_cast<TH2F*>(fListOfHistos->At(9)); | |
337 | fV0multCcorr = dynamic_cast<TH2F*>(fListOfHistos->At(10)); | |
338 | fV0Acorr = dynamic_cast<TH2F*>(fListOfHistos->At(11)); | |
339 | fV0Ccorr = dynamic_cast<TH2F*>(fListOfHistos->At(12)); | |
340 | ||
c543bd05 | 341 | fT0ampl = dynamic_cast<TH1F*>(fListOfHistos->At(13)); |
342 | fV0ampl = dynamic_cast<TH1F*>(fListOfHistos->At(14)); | |
343 | fT0time2 = dynamic_cast<TH1F*>(fListOfHistos->At(15)); | |
344 | ||
345 | ||
684b93f0 | 346 | // draw the histograms if not in batch mode |
347 | if (!gROOT->IsBatch()) { | |
348 | new TCanvas; | |
349 | fT0vtxRec->DrawCopy("E"); | |
350 | new TCanvas; | |
351 | fT0time->DrawCopy("E"); | |
352 | new TCanvas; | |
353 | fT0mult->DrawCopy("E"); | |
354 | new TCanvas; | |
355 | fT0vtxRes->DrawCopy("E"); | |
356 | new TCanvas; | |
357 | fT0vtxRecGen->DrawCopy(); | |
358 | new TCanvas; | |
359 | fV0a->DrawCopy("E"); | |
360 | new TCanvas; | |
361 | fV0c->DrawCopy("E"); | |
362 | new TCanvas; | |
363 | fV0multA->DrawCopy("E"); | |
364 | new TCanvas; | |
365 | fV0multC->DrawCopy("E"); | |
366 | new TCanvas; | |
367 | fV0multAcorr->DrawCopy(); | |
368 | new TCanvas; | |
369 | fV0multCcorr->DrawCopy(); | |
370 | new TCanvas; | |
371 | fV0Acorr->DrawCopy(); | |
372 | new TCanvas; | |
373 | fV0Ccorr->DrawCopy(); | |
c543bd05 | 374 | new TCanvas; |
375 | fT0ampl->DrawCopy("E"); | |
376 | new TCanvas; | |
377 | fV0ampl->DrawCopy("E"); | |
378 | new TCanvas; | |
379 | fT0time2->DrawCopy("E"); | |
684b93f0 | 380 | } |
381 | ||
382 | // write the output histograms to a file | |
383 | TFile* outputFile = TFile::Open("FwdDetsQA.root", "recreate"); | |
384 | if (!outputFile || !outputFile->IsOpen()) | |
385 | { | |
386 | Error("AliAnaFwdDetsQA", "opening output file FwdDetsQA.root failed"); | |
387 | return; | |
388 | } | |
389 | fT0vtxRec->Write(); | |
390 | fT0time->Write(); | |
391 | fT0mult->Write(); | |
392 | fT0vtxRecGen->Write(); | |
393 | fT0vtxRes->Write(); | |
c543bd05 | 394 | fT0ampl->Write(); |
395 | fT0time2->Write(); | |
684b93f0 | 396 | |
397 | fV0a->Write(); | |
398 | fV0c->Write(); | |
399 | fV0multA->Write(); | |
400 | fV0multC->Write(); | |
401 | fV0multAcorr->Write(); | |
402 | fV0multCcorr->Write(); | |
403 | fV0Acorr->Write(); | |
404 | fV0Ccorr->Write(); | |
c543bd05 | 405 | fV0ampl->Write(); |
684b93f0 | 406 | |
407 | outputFile->Close(); | |
408 | delete outputFile; | |
409 | ||
410 | //delete esd; | |
411 | Info("AliAnaFwdDetsQA", "Successfully finished"); | |
412 | } | |
413 | ||
414 | ||
415 | ||
416 | ||
417 | ||
418 | ||
419 | ||
420 | ||
421 | ||
422 | ||
423 | ||
424 | ||
425 | ||
426 | ||
427 | ||
428 | ||
429 | ||
430 | ||
431 | ||
432 | ||
433 | ||
434 | ||
435 | ||
436 | ||
437 | ||
438 | ||
439 | ||
440 | ||
441 | ||
442 | ||
443 | ||
444 | ||
445 | ||
446 | ||
447 | ||
448 | ||
449 | ||
450 | ||
451 | ||
452 | ||
453 | ||
454 | ||
455 | ||
456 | ||
457 | ||
458 | ||
459 | ||
460 | ||
461 | ||
462 | ||
463 | ||
464 | ||
465 | ||
466 | ||
467 | ||
468 | ||
469 | ||
470 | ||
471 | ||
472 | ||
473 | ||
474 |