1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //------------------------------
24 //------------------------------
34 #include "TParticle.h"
35 #include "AliESDEvent.h"
37 #include "AliESDcascade.h"
38 #include "AliESDMuonTrack.h"
39 #include "AliESDCaloCluster.h"
41 #include "AliMCEvent.h"
42 #include "AliGenEventHeader.h"
45 #include "AliAnalysisTaskCheckESD.h"
47 ClassImp(AliAnalysisTaskCheckESD)
49 AliAnalysisTaskCheckESD::AliAnalysisTaskCheckESD():
50 AliAnalysisTaskSE("AliAnalysisTaskCheckESD"),
72 // Default constructor
73 // Define input and output slots here
74 // Input slot #0 works with a TChain
75 DefineInput(0, TChain::Class());
76 // Output slot #1 TList
77 DefineOutput(1, TList::Class());
80 AliAnalysisTaskCheckESD::AliAnalysisTaskCheckESD(const char* name):
81 AliAnalysisTaskSE(name),
104 AliInfo("Constructor AliAnalysisTaskCheckESD");
105 // Define input and output slots here
106 // Input slot #0 works with a TChain
107 DefineInput(0, TChain::Class());
108 // Output slot #1 TList
109 DefineOutput(1, TList::Class());
112 TH1F * AliAnalysisTaskCheckESD::CreateHisto(const char* name, const char* title,Int_t nBins,
113 Double_t xMin, Double_t xMax,
114 const char* xLabel, const char* yLabel)
116 // create a histogram
117 TH1F* result = new TH1F(name, title, nBins, xMin, xMax);
118 result->SetOption("E");
119 if (xLabel) result->GetXaxis()->SetTitle(xLabel);
120 if (yLabel) result->GetYaxis()->SetTitle(yLabel);
121 result->SetMarkerStyle(kFullCircle);
125 TH1F *AliAnalysisTaskCheckESD::CreateEffHisto(const TH1F* hGen, const TH1F* hRec)
127 // create an efficiency histogram
128 Int_t nBins = hGen->GetNbinsX();
129 TH1F* hEff = (TH1F*) hGen->Clone("hEff");
131 hEff->SetStats(kFALSE);
132 hEff->SetMinimum(0.);
133 hEff->SetMaximum(110.);
134 hEff->GetYaxis()->SetTitle("#epsilon [%]");
136 for (Int_t iBin = 0; iBin <= nBins; iBin++) {
137 Double_t nGenEff = hGen->GetBinContent(iBin);
138 Double_t nRecEff = hRec->GetBinContent(iBin);
140 Double_t eff = nRecEff/nGenEff;
141 hEff->SetBinContent(iBin, 100. * eff);
142 Double_t error = sqrt(eff*(1.-eff) / nGenEff);
143 if (TMath::Abs(error) < 1.e-33) error = 0.0001;
144 hEff->SetBinError(iBin, 100. * error);
147 hEff->SetBinContent(iBin, -100.);
148 hEff->SetBinError(iBin, 0);
155 Bool_t AliAnalysisTaskCheckESD::FitHisto(TH1* histo, Double_t& res, Double_t& resError)
157 // fit a gaussian to a histogram
158 static TF1* fitFunc = new TF1("fitFunc", "gaus");
159 fitFunc->SetLineWidth(2);
160 fitFunc->SetFillStyle(0);
161 Double_t maxFitRange = 2;
163 if (histo->Integral() > 50) {
164 Float_t mean = histo->GetMean();
165 Float_t rms = histo->GetRMS();
166 fitFunc->SetRange(mean - maxFitRange*rms, mean + maxFitRange*rms);
167 fitFunc->SetParameters(mean, rms);
168 histo->Fit(fitFunc, "QRI0");
169 histo->GetFunction("fitFunc")->ResetBit(1<<9);
170 res = TMath::Abs(fitFunc->GetParameter(2));
171 resError = TMath::Abs(fitFunc->GetParError(2));
177 void AliAnalysisTaskCheckESD::UserCreateOutputObjects()
180 AliInfo("AliAnalysisTaskCheckESD::UserCreateOutputObjects");
181 // Create output container
182 fListOfHistos = new TList();
184 // efficiency and resolution histog
189 fGen = CreateHisto("hGen", "generated tracks", nBinsPt, minPt, maxPt, "p_{t} [GeV/c]", "N");
190 fRec = CreateHisto("hRec", "reconstructed tracks", nBinsPt, minPt, maxPt, "p_{t} [GeV/c]", "N");
191 fResPtInv = CreateHisto("hResPtInv", "", 100, -10, 10, "(p_{t,rec}^{-1}-p_{t,sim}^{-1}) / p_{t,sim}^{-1} [%]", "N");
192 fResPhi = CreateHisto("hResPhi", "", 100, -20, 20, "#phi_{rec}-#phi_{sim} [mrad]", "N");
193 fResTheta = CreateHisto("hResTheta", "", 100, -20, 20, "#theta_{rec}-#theta_{sim} [mrad]", "N");
196 fDEdxRight = new TH2F("hDEdxRight", "", 300, 0, 3, 100, 0, 400);
197 fDEdxRight->SetStats(kFALSE);
198 fDEdxRight->GetXaxis()->SetTitle("p [GeV/c]");
199 fDEdxRight->GetYaxis()->SetTitle("dE/dx_{TPC}");
200 fDEdxRight->SetMarkerStyle(kFullCircle);
201 fDEdxRight->SetMarkerSize(0.4);
202 fDEdxWrong = new TH2F("hDEdxWrong", "", 300, 0, 3, 100, 0, 400);
203 fDEdxWrong->SetStats(kFALSE);
204 fDEdxWrong->GetXaxis()->SetTitle("p [GeV/c]");
205 fDEdxWrong->GetYaxis()->SetTitle("dE/dx_{TPC}");
206 fDEdxWrong->SetMarkerStyle(kFullCircle);
207 fDEdxWrong->SetMarkerSize(0.4);
208 fDEdxWrong->SetMarkerColor(kRed);
210 fResTOFRight = CreateHisto("hResTOFRight", "", 100, -1000, 1000, "t_{TOF}-t_{track} [ps]", "N");
211 fResTOFWrong = CreateHisto("hResTOFWrong", "", 100, -1000, 1000, "t_{TOF}-t_{track} [ps]", "N");
212 fResTOFWrong->SetLineColor(kRed);
215 fEPHOS = CreateHisto("hEPHOS", "PHOS", 100, 0, 50, "E [GeV]", "N");
216 fEEMCAL = CreateHisto("hEEMCAL", "EMCAL", 100, 0, 50, "E [GeV]", "N");
219 fPtMUON = CreateHisto("hPtMUON", "MUON", 100, 0, 20, "p_{t} [GeV/c]", "N");
222 fMassK0 = CreateHisto("hMassK0", "K^{0}", 100, 0.4, 0.6, "M(#pi^{+}#pi^{-}) [GeV/c^{2}]", "N");
223 fMassLambda = CreateHisto("hMassLambda", "#Lambda", 100, 1.0, 1.2, "M(p#pi^{-}) [GeV/c^{2}]", "N");
225 fMassLambdaBar = CreateHisto("hMassLambdaBar", "#bar{#Lambda}", 100, 1.0, 1.2, "M(#bar{p}#pi^{+}) [GeV/c^{2}]", "N");
226 fMassXi = CreateHisto("hMassXi", "#Xi", 100, 1.2, 1.5, "M(#Lambda#pi) [GeV/c^{2}]", "N");
227 fMassOmega = CreateHisto("hMassOmega", "#Omega", 100, 1.5, 1.8, "M(#LambdaK) [GeV/c^{2}]", "N");
228 fScalars = new TH1F("hScalars","Container of scalars",8,0,8);
229 fArrayHist = new TH1F("hArrayHist","Container for Array",
230 (AliPID::kSPECIES+1)*AliPID::kSPECIES,0,(AliPID::kSPECIES+1)*AliPID::kSPECIES);
232 fListOfHistos->Add(fGen);
233 fListOfHistos->Add(fRec);
234 fListOfHistos->Add(fResPtInv);
235 fListOfHistos->Add(fResPhi);
236 fListOfHistos->Add(fResTheta);
237 fListOfHistos->Add(fDEdxRight);
238 fListOfHistos->Add(fDEdxWrong);
239 fListOfHistos->Add(fResTOFRight);
240 fListOfHistos->Add(fResTOFWrong);
241 fListOfHistos->Add(fEPHOS);
242 fListOfHistos->Add(fEEMCAL);
243 fListOfHistos->Add(fPtMUON);
244 fListOfHistos->Add(fMassK0);
245 fListOfHistos->Add(fMassLambda);
246 fListOfHistos->Add(fMassLambdaBar);
247 fListOfHistos->Add(fMassXi);
248 fListOfHistos->Add(fMassOmega);
249 fListOfHistos->Add(fScalars);
250 fListOfHistos->Add(fArrayHist);
253 void AliAnalysisTaskCheckESD::UserExec(Option_t */*option*/)
255 // check the content of the ESD
256 Double_t cutPtV0 = 0.3;
257 Double_t cutPtCascade = 0.5;
261 Int_t partCode[AliPID::kSPECIES] = {kElectron, kMuonMinus, kPiPlus, kKPlus, kProton};
262 Double_t partFrac[AliPID::kSPECIES] = {0.01, 0.01, 0.85, 0.10, 0.05};
264 AliMCEvent* mcEvent = MCEvent();
266 Printf("ERROR: Could not retrieve MC event");
270 //Primary vertex needed
272 mcEvent->GenEventHeader()->PrimaryVertex(vertex);
274 TObjArray selParticles; //Use TClonesArray?
276 TObjArray selCascades;
277 for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfTracks(); iTracks++) {//Check this loop again
278 AliMCParticle* track = (AliMCParticle*)mcEvent->GetTrack(iTracks);
279 TParticle* particle = track->Particle();
280 if (!particle) continue;
281 if (particle->Pt() < 0.001) continue;
282 if (TMath::Abs(particle->Eta()) > 0.9) continue;
283 TVector3 dVertex(particle->Vx() - vertex[0], particle->Vy() - vertex[1], particle->Vz() - vertex[2]);
284 if (dVertex.Mag() > 0.0001) continue;
286 switch (TMath::Abs(particle->GetPdgCode())) {
293 if (particle->Pt() > minPt) {
294 selParticles.Add(particle);
296 fGen->Fill(particle->Pt());
303 if (particle->Pt() > cutPtV0) {
305 selV0s.Add(particle);
312 if (particle->Pt() > cutPtCascade) {
314 selCascades.Add(particle);
323 AliVEvent* event = InputEvent();
325 Printf("ERROR: Could not retrieve event");
329 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
331 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfTracks(); iTrack++) {
332 AliESDtrack* track = esd->GetTrack(iTrack);
334 // select tracks of selected particles
335 Int_t label = TMath::Abs(track->GetLabel());
336 if (label > mcEvent->GetNumberOfTracks()) continue; // background
337 AliMCParticle* mctrack = (AliMCParticle*)mcEvent->GetTrack(label);
338 TParticle* particle = mctrack->Particle();
339 if (!selParticles.Contains(particle)) continue;
340 if ((track->GetStatus() & AliESDtrack::kITSrefit) == 0) continue;
341 if (track->GetConstrainedChi2() > 1e9) continue;
342 selParticles.Remove(particle); // don't count multiple tracks
345 fRec->Fill(particle->Pt());
346 if (track->GetLabel() < 0) {
351 fResPtInv->Fill(100. * (TMath::Abs(track->GetSigned1Pt()) - 1./particle->Pt()) *particle->Pt());
352 fResPhi->Fill(1000. * (track->Phi() - particle->Phi()));
353 fResTheta->Fill(1000. * (track->Theta() - particle->Theta()));
356 if ((track->GetStatus() & AliESDtrack::kESDpid) == 0) continue;
359 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
360 if (TMath::Abs(particle->GetPdgCode()) == partCode[i]) iGen = i;
362 Double_t probability[AliPID::kSPECIES];
363 track->GetESDpid(probability);
366 for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
367 probability[i] *= partFrac[i];
368 if (probability[i] > pMax) {
369 pMax = probability[i];
374 fArrayHist->Fill(AliPID::kSPECIES*iGen + iRec);
380 Double_t time[AliPID::kSPECIESC];
381 track->GetIntegratedTimes(time,AliPID::kSPECIESC);
383 fDEdxRight->Fill(particle->P(), track->GetTPCsignal());
384 if ((track->GetStatus() & AliESDtrack::kTOFpid) != 0) {
385 fResTOFRight->Fill(track->GetTOFsignal() - time[iRec]);
389 fDEdxWrong->Fill(particle->P(), track->GetTPCsignal());
390 if ((track->GetStatus() & AliESDtrack::kTOFpid) != 0) {
391 fResTOFWrong->Fill(track->GetTOFsignal() - time[iRec]);
396 // loop over muon tracks
397 for (Int_t iTrack = 0; iTrack < esd->GetNumberOfMuonTracks(); iTrack++) {
398 AliESDMuonTrack* muonTrack = esd->GetMuonTrack(iTrack);
399 Double_t ptInv = TMath::Abs(muonTrack->GetInverseBendingMomentum());
401 fPtMUON->Fill(1./ptInv);
406 for (Int_t iV0 = 0; iV0 < esd->GetNumberOfV0s(); iV0++) {
407 AliESDv0* v0 = esd->GetV0(iV0);
408 if (v0->GetOnFlyStatus()) continue;
409 v0->ChangeMassHypothesis(kK0Short);
410 fMassK0->Fill(v0->GetEffMass());
411 v0->ChangeMassHypothesis(kLambda0);
412 fMassLambda->Fill(v0->GetEffMass());
413 v0->ChangeMassHypothesis(kLambda0Bar);
414 fMassLambdaBar->Fill(v0->GetEffMass());
416 Int_t negLabel = TMath::Abs(esd->GetTrack(v0->GetNindex())->GetLabel());
417 if (negLabel > mcEvent->GetNumberOfTracks()) continue; // background
418 AliMCParticle* track = (AliMCParticle*)mcEvent->GetTrack(negLabel);
419 TParticle* particle = track->Particle();
420 Int_t negMother = particle->GetMother(0);
421 if (negMother < 0) continue;
422 Int_t posLabel = TMath::Abs(esd->GetTrack(v0->GetPindex())->GetLabel());
423 if (posLabel > mcEvent->GetNumberOfTracks()) continue; // background
424 track = (AliMCParticle*)mcEvent->GetTrack(posLabel);
425 particle = track->Particle();
426 Int_t posMother = particle->GetMother(0);
427 if (negMother != posMother) continue;
428 track = (AliMCParticle*)mcEvent->GetTrack(negMother);
429 particle = track->Particle();
430 if (!selV0s.Contains(particle)) continue;
431 selV0s.Remove(particle);
435 // loop over Cascades
436 for (Int_t iCascade = 0; iCascade < esd->GetNumberOfCascades(); iCascade++) {
437 AliESDcascade* cascade = esd->GetCascade(iCascade);
439 cascade->ChangeMassHypothesis(v0q,kXiMinus);
440 fMassXi->Fill(cascade->GetEffMassXi());
441 cascade->ChangeMassHypothesis(v0q,kOmegaMinus);
442 fMassOmega->Fill(cascade->GetEffMassXi());
444 Int_t negLabel = TMath::Abs(esd->GetTrack(cascade->GetNindex())->GetLabel());
445 if (negLabel > mcEvent->GetNumberOfTracks()) continue; // background
446 AliMCParticle* track = (AliMCParticle*)mcEvent->GetTrack(negLabel);
447 TParticle* particle = track->Particle();
448 Int_t negMother = particle->GetMother(0);
449 if (negMother < 0) continue;
450 Int_t posLabel = TMath::Abs(esd->GetTrack(cascade->GetPindex())->GetLabel());
451 if (posLabel > mcEvent->GetNumberOfTracks()) continue; // background
452 track = (AliMCParticle*)mcEvent->GetTrack(posLabel);
453 particle = track->Particle();
454 Int_t posMother = particle->GetMother(0);
455 if (negMother != posMother) continue;
456 track = (AliMCParticle*)mcEvent->GetTrack(negMother);
457 particle = track->Particle();
458 Int_t v0Mother = particle->GetMother(0);
459 if (v0Mother < 0) continue;
460 Int_t bacLabel = TMath::Abs(esd->GetTrack(cascade->GetBindex())->GetLabel());
461 if (bacLabel > mcEvent->GetNumberOfTracks()) continue; // background
462 track = (AliMCParticle*)mcEvent->GetTrack(bacLabel);
463 particle = track->Particle();
464 Int_t bacMother = particle->GetMother(0);
465 if (v0Mother != bacMother) continue;
466 track = (AliMCParticle*)mcEvent->GetTrack(v0Mother);
467 particle = track->Particle();
468 if (!selCascades.Contains(particle)) continue;
469 selCascades.Remove(particle);
473 // loop over the clusters
474 for (Int_t iCluster=0; iCluster<esd->GetNumberOfCaloClusters(); iCluster++) {
475 AliESDCaloCluster * clust = esd->GetCaloCluster(iCluster);
476 if (clust->IsPHOS()) fEPHOS->Fill(clust->E());
477 if (clust->IsEMCAL()) fEEMCAL->Fill(clust->E());
481 PostData(1, fListOfHistos);
484 void AliAnalysisTaskCheckESD::Terminate(Option_t *)
487 Int_t checkNGenLow = 1;
488 Double_t checkEffLow = 0.5;
489 Double_t checkEffSigma = 3;
490 Double_t checkFakeHigh = 0.5;
491 Double_t checkFakeSigma = 3;
493 Double_t checkResPtInvHigh = 5;
494 Double_t checkResPtInvSigma = 3;
495 Double_t checkResPhiHigh = 10;
496 Double_t checkResPhiSigma = 3;
497 Double_t checkResThetaHigh = 10;
498 Double_t checkResThetaSigma = 3;
500 Double_t checkPIDEffLow = 0.5;
501 Double_t checkPIDEffSigma = 3;
502 Double_t checkResTOFHigh = 500;
503 Double_t checkResTOFSigma = 3;
505 Double_t checkPHOSNLow = 5;
506 Double_t checkPHOSEnergyLow = 0.3;
507 Double_t checkPHOSEnergyHigh = 1.0;
508 Double_t checkEMCALNLow = 50;
509 Double_t checkEMCALEnergyLow = 0.05;
510 Double_t checkEMCALEnergyHigh = 1.0;
512 Double_t checkMUONNLow = 1;
513 Double_t checkMUONPtLow = 0.5;
514 Double_t checkMUONPtHigh = 10.;
516 Double_t checkV0EffLow = 0.02;
517 Double_t checkV0EffSigma = 3;
518 Double_t checkCascadeEffLow = 0.01;
519 Double_t checkCascadeEffSigma = 3;
521 const char* partName[AliPID::kSPECIES+1] = {"electron", "muon", "pion", "kaon", "proton", "other"};
522 fListOfHistos = dynamic_cast<TList*>(GetOutputData(1));
523 if (!fListOfHistos) {
524 Printf("ERROR: fListOfHistos not available");
528 fGen = dynamic_cast<TH1F*>(fListOfHistos->At(0));
529 fRec = dynamic_cast<TH1F*>(fListOfHistos->At(1));
530 fResPtInv = dynamic_cast<TH1F*>(fListOfHistos->At(2));
531 fResPhi = dynamic_cast<TH1F*>(fListOfHistos->At(3));
532 fResTheta = dynamic_cast<TH1F*>(fListOfHistos->At(4));
533 fDEdxRight = dynamic_cast<TH2F*>(fListOfHistos->At(5));
534 fDEdxWrong = dynamic_cast<TH2F*>(fListOfHistos->At(6));
535 fResTOFRight = dynamic_cast<TH1F*>(fListOfHistos->At(7));
536 fResTOFWrong = dynamic_cast<TH1F*>(fListOfHistos->At(8));
537 fEPHOS = dynamic_cast<TH1F*>(fListOfHistos->At(9));
538 fEEMCAL = dynamic_cast<TH1F*>(fListOfHistos->At(10));
539 fPtMUON = dynamic_cast<TH1F*>(fListOfHistos->At(11));
540 fMassK0 = dynamic_cast<TH1F*>(fListOfHistos->At(12));
541 fMassLambda = dynamic_cast<TH1F*>(fListOfHistos->At(13));
542 fMassLambdaBar = dynamic_cast<TH1F*>(fListOfHistos->At(14));
543 fMassXi = dynamic_cast<TH1F*>(fListOfHistos->At(15));
544 fMassOmega = dynamic_cast<TH1F*>(fListOfHistos->At(16));
545 fScalars = dynamic_cast<TH1F*>(fListOfHistos->At(17));
546 fArrayHist = dynamic_cast<TH1F*>(fListOfHistos->At(18));
548 Int_t nGen = Int_t(fScalars->GetBinContent(1));
549 Int_t nRec = Int_t(fScalars->GetBinContent(2));
550 Int_t nFake = Int_t(fScalars->GetBinContent(3));
551 Int_t nGenV0s = Int_t(fScalars->GetBinContent(4));
552 Int_t nRecV0s = Int_t(fScalars->GetBinContent(5));
553 Int_t nIdentified = Int_t(fScalars->GetBinContent(6));
554 Int_t nGenCascades = Int_t(fScalars->GetBinContent(7));
555 Int_t nRecCascades = Int_t(fScalars->GetBinContent(8));
559 Int_t identified[AliPID::kSPECIES+1][AliPID::kSPECIES];
560 for(Int_t i = 0; i < (AliPID::kSPECIES+1); i++)
561 for(Int_t j = 0; j < AliPID::kSPECIES; j++) {
562 identified[i][j] = Int_t(fArrayHist->GetBinContent(k));
568 if (nGen < checkNGenLow) {
569 Warning("CheckESD", "low number of generated particles: %d", Int_t(nGen));
572 TH1F* hEff = CreateEffHisto(fGen, fRec);
574 Info("CheckESD", "%d out of %d tracks reconstructed including %d "
575 "fake tracks", nRec, nGen, nFake);
578 Double_t eff = nRec*1./nGen;
579 Double_t effError = TMath::Sqrt(eff*(1.-eff) / nGen);
580 Double_t fake = nFake*1./nGen;
581 Double_t fakeError = TMath::Sqrt(fake*(1.-fake) / nGen);
582 Info("CheckESD", "eff = (%.1f +- %.1f) %% fake = (%.1f +- %.1f) %%",
583 100.*eff, 100.*effError, 100.*fake, 100.*fakeError);
584 if (eff < checkEffLow - checkEffSigma*effError) {
585 Warning("CheckESD", "low efficiency: (%.1f +- %.1f) %%",100.*eff, 100.*effError);
587 if (fake > checkFakeHigh + checkFakeSigma*fakeError) {
588 Warning("CheckESD", "high fake: (%.1f +- %.1f) %%",100.*fake, 100.*fakeError);
591 Double_t res, resError;
592 if (FitHisto(fResPtInv, res, resError)) {
593 Info("CheckESD", "relative inverse pt resolution = (%.1f +- %.1f) %%",res, resError);
594 if (res > checkResPtInvHigh + checkResPtInvSigma*resError) {
595 Warning("CheckESD", "bad pt resolution: (%.1f +- %.1f) %%",res, resError);
598 if (FitHisto(fResPhi, res, resError)) {
599 Info("CheckESD", "phi resolution = (%.1f +- %.1f) mrad", res, resError);
600 if (res > checkResPhiHigh + checkResPhiSigma*resError) {
601 Warning("CheckESD", "bad phi resolution: (%.1f +- %.1f) mrad",
606 if (FitHisto(fResTheta, res, resError)) {
607 Info("CheckESD", "theta resolution = (%.1f +- %.1f) mrad",
609 if (res > checkResThetaHigh + checkResThetaSigma*resError) {
610 Warning("CheckESD", "bad theta resolution: (%.1f +- %.1f) mrad",
616 Double_t eff = nIdentified*1./nRec;
617 Double_t effError = TMath::Sqrt(eff*(1.-eff) / nRec);
618 Info("CheckESD", "PID eff = (%.1f +- %.1f) %%",
619 100.*eff, 100.*effError);
620 if (eff < checkPIDEffLow - checkPIDEffSigma*effError) {
621 Warning("CheckESD", "low PID efficiency: (%.1f +- %.1f) %%",
622 100.*eff, 100.*effError);
626 printf("%9s:", "gen\\rec");
627 for (Int_t iRec = 0; iRec < AliPID::kSPECIES; iRec++) {
628 printf("%9s", partName[iRec]);
631 for (Int_t iGen = 0; iGen < AliPID::kSPECIES+1; iGen++) {
632 printf("%9s:", partName[iGen]);
633 for (Int_t iRec = 0; iRec < AliPID::kSPECIES; iRec++) {
634 printf("%9d", identified[iGen][iRec]);
639 if (FitHisto(fResTOFRight, res, resError)) {
640 Info("CheckESD", "TOF resolution = (%.1f +- %.1f) ps", res, resError);
641 if (res > checkResTOFHigh + checkResTOFSigma*resError) {
642 Warning("CheckESD", "bad TOF resolution: (%.1f +- %.1f) ps",
648 if (fEPHOS->Integral() < checkPHOSNLow) {
649 Warning("CheckESD", "low number of PHOS particles: %d",
650 Int_t(fEPHOS->Integral()));
653 Double_t mean = fEPHOS->GetMean();
654 if (mean < checkPHOSEnergyLow) {
655 Warning("CheckESD", "low mean PHOS energy: %.1f GeV", mean);
656 } else if (mean > checkPHOSEnergyHigh) {
657 Warning("CheckESD", "high mean PHOS energy: %.1f GeV", mean);
661 if (fEEMCAL->Integral() < checkEMCALNLow) {
662 Warning("CheckESD", "low number of EMCAL particles: %d",
663 Int_t(fEEMCAL->Integral()));
666 Double_t mean = fEEMCAL->GetMean();
667 if (mean < checkEMCALEnergyLow) {
668 Warning("CheckESD", "low mean EMCAL energy: %.1f GeV", mean);
670 else if (mean > checkEMCALEnergyHigh) {
671 Warning("CheckESD", "high mean EMCAL energy: %.1f GeV", mean);
676 if (fPtMUON->Integral() < checkMUONNLow) {
677 Warning("CheckESD", "low number of MUON particles: %d",
678 Int_t(fPtMUON->Integral()));
681 Double_t mean = fPtMUON->GetMean();
682 if (mean < checkMUONPtLow) {
683 Warning("CheckESD", "low mean MUON pt: %.1f GeV/c", mean);
685 else if (mean > checkMUONPtHigh) {
686 Warning("CheckESD", "high mean MUON pt: %.1f GeV/c", mean);
692 Double_t eff = nRecV0s*1./nGenV0s;
693 Double_t effError = TMath::Sqrt(eff*(1.-eff) / nGenV0s);
694 if (TMath::Abs(effError) < 1.e-33) effError = checkV0EffLow / TMath::Sqrt(1.*nGenV0s);
695 Info("CheckESD", "V0 eff = (%.1f +- %.1f) %%",
696 100.*eff, 100.*effError);
697 if (eff < checkV0EffLow - checkV0EffSigma*effError) {
698 Warning("CheckESD", "low V0 efficiency: (%.1f +- %.1f) %%",
699 100.*eff, 100.*effError);
704 if (nGenCascades > 0) {
705 Double_t eff = nRecCascades*1./nGenCascades;
706 Double_t effError = TMath::Sqrt(eff*(1.-eff) / nGenCascades);
707 if (effError == 0) effError = checkV0EffLow /
708 TMath::Sqrt(1.*nGenCascades);
709 Info("CheckESD", "Cascade eff = (%.1f +- %.1f) %%",
710 100.*eff, 100.*effError);
711 if (eff < checkCascadeEffLow - checkCascadeEffSigma*effError) {
712 Warning("CheckESD", "low Cascade efficiency: (%.1f +- %.1f) %%",
713 100.*eff, 100.*effError);
719 // draw the histograms if not in batch mode
720 if (!gROOT->IsBatch()) {
724 fResPtInv->DrawCopy("E");
726 fResPhi->DrawCopy("E");
728 fResTheta->DrawCopy("E");
730 fDEdxRight->DrawCopy();
731 fDEdxWrong->DrawCopy("SAME");
733 fResTOFRight->DrawCopy("E");
734 fResTOFWrong->DrawCopy("SAME");
736 fEPHOS->DrawCopy("E");
738 fEEMCAL->DrawCopy("E");
740 fPtMUON->DrawCopy("E");
742 fMassK0->DrawCopy("E");
744 fMassLambda->DrawCopy("E");
746 fMassLambdaBar->DrawCopy("E");
748 fMassXi->DrawCopy("E");
750 fMassOmega->DrawCopy("E");
753 // write the output histograms to a file
754 TFile* outputFile = TFile::Open("check.root", "recreate");
755 if (!outputFile || !outputFile->IsOpen())
757 Error("CheckESD", "opening output file check.root failed");
766 fResTOFRight->Write();
767 fResTOFWrong->Write();
772 fMassLambda->Write();
773 fMassLambdaBar->Write();
795 delete hMassLambdaBar;
799 delete hArrayHist; */
803 Info("CheckESD", "check of ESD was successfull");