]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliHFEtpcPIDqa.cxx
Merge branch 'TPCdev' of https://git.cern.ch/reps/AliRoot into TPCdev
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliHFEtpcPIDqa.cxx
CommitLineData
3a72645a 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// Class AliHFEtpcPIDqa
17// Monitoring TPC PID in the HFE PID montioring framework. The following
18// quantities are monitored:
19// TPC dE/dx (Number of sigmas)
20// TPC dE/dx (Absolute values)
21// (Always as function of momentum, particle species and centrality
22// before and after cut)
23// More information about the PID monitoring framework can be found in
24// AliHFEpidQAmanager.cxx and AliHFEdetPIDqa.cxx
25//
26// Author:
27// Markus Fasel <M.Fasel@gsi.de>
28//
bf892a6a 29#include <TBrowser.h>
3a72645a 30#include <TClass.h>
31#include <TH2.h>
32#include <THnSparse.h>
33#include <TString.h>
34
35#include "AliAODTrack.h"
36#include "AliESDtrack.h"
3a72645a 37#include "AliLog.h"
38#include "AliPID.h"
8c1c76e9 39#include "AliPIDResponse.h"
3a72645a 40
41#include "AliHFEcollection.h"
42#include "AliHFEpidBase.h"
6555e2ad 43#include "AliHFEpidQAmanager.h"
44#include "AliHFEpidTPC.h"
3a72645a 45#include "AliHFEtools.h"
46#include "AliHFEtpcPIDqa.h"
47
6555e2ad 48ClassImp(AliHFEtpcPIDqa)
49
3a72645a 50//_________________________________________________________
51AliHFEtpcPIDqa::AliHFEtpcPIDqa():
52 AliHFEdetPIDqa()
53 , fHistos(NULL)
e156c3bb 54 , fBrowseCentrality(-1)
3a72645a 55{
56 //
57 // Dummy constructor
58 //
59}
60
61//_________________________________________________________
62AliHFEtpcPIDqa::AliHFEtpcPIDqa(const char* name):
63 AliHFEdetPIDqa(name, "QA for TPC")
64 , fHistos(NULL)
e156c3bb 65 , fBrowseCentrality(-1)
3a72645a 66{
67 //
68 // Default constructor
69 //
70}
71
72//_________________________________________________________
73AliHFEtpcPIDqa::AliHFEtpcPIDqa(const AliHFEtpcPIDqa &o):
74 AliHFEdetPIDqa(o)
71478e06 75 , fHistos(NULL)
e156c3bb 76 , fBrowseCentrality(o.fBrowseCentrality)
3a72645a 77{
78 //
79 // Copy constructor
80 //
81 o.Copy(*this);
82}
83
84//_________________________________________________________
85AliHFEtpcPIDqa &AliHFEtpcPIDqa::operator=(const AliHFEtpcPIDqa &o){
86 //
87 // Do assignment
88 //
89 AliHFEdetPIDqa::operator=(o);
90 if(&o != this) o.Copy(*this);
91 return *this;
92}
93
94//_________________________________________________________
95AliHFEtpcPIDqa::~AliHFEtpcPIDqa(){
96 //
97 // Destructor
98 //
99 if(fHistos) delete fHistos;
100}
101
102//_________________________________________________________
103void AliHFEtpcPIDqa::Copy(TObject &o) const {
104 //
105 // Make copy
106 //
107 AliHFEtpcPIDqa &target = dynamic_cast<AliHFEtpcPIDqa &>(o);
108 if(target.fHistos){
109 delete target.fHistos;
110 target.fHistos = NULL;
111 }
112 if(fHistos) target.fHistos = new AliHFEcollection(*fHistos);
e156c3bb 113 target.fBrowseCentrality = fBrowseCentrality;
3a72645a 114}
115
116//_________________________________________________________
117Long64_t AliHFEtpcPIDqa::Merge(TCollection *coll){
118 //
119 // Merge with other objects
120 //
121 if(!coll) return 0;
122 if(coll->IsEmpty()) return 1;
123
124 TIter it(coll);
125 AliHFEtpcPIDqa *refQA = NULL;
126 TObject *o = NULL;
127 Long64_t count = 0;
128 TList listHistos;
129 while((o = it())){
130 refQA = dynamic_cast<AliHFEtpcPIDqa *>(o);
131 if(!refQA) continue;
132
133 listHistos.Add(refQA->fHistos);
134 count++;
135 }
136 fHistos->Merge(&listHistos);
137 return count + 1;
138}
139
bf892a6a 140//_________________________________________________________
141void AliHFEtpcPIDqa::Browse(TBrowser *b){
142 //
143 // Browse the PID QA
144 //
145 if(b){
146 if(fHistos){
147 b->Add(fHistos, fHistos->GetName());
148
149 // Make Projections of the dE/dx Spectra and add them to a new Folder
c2690925 150 TString specnames[AliPID::kSPECIES+1] = {"All", "Electrons", "Muon", "Pions", "Kaon", "Protons"};
151 Int_t specind[AliPID::kSPECIES+1] = {-1, AliPID::kElectron, AliPID::kMuon, AliPID::kPion, AliPID::kKaon, AliPID::kProton};
bf892a6a 152 TList *listdEdx = new TList;
153 listdEdx->SetOwner();
154 TList *listNsigma = new TList;
155 listNsigma->SetOwner();
156
157 TH2 *hptr = NULL;
c2690925 158 for(Int_t ispec = 0; ispec < AliPID::kSPECIES+1; ispec++){
bf892a6a 159 for(Int_t istep = 0; istep < 2; istep++){
e156c3bb 160 hptr = MakeSpectrumdEdx(static_cast<AliHFEdetPIDqa::EStep_t>(istep), specind[ispec], fBrowseCentrality);
161 hptr->SetName(Form("hTPCdEdx%s%s%s", specnames[ispec].Data(), istep == 0 ? "Before" : "After" , fBrowseCentrality == -1 ? "MinBias" : Form("Cent%d", fBrowseCentrality)));
bf892a6a 162 listdEdx->Add(hptr);
e156c3bb 163 hptr = MakeSpectrumNSigma(static_cast<AliHFEdetPIDqa::EStep_t>(istep), specind[ispec], fBrowseCentrality);
164 hptr->SetName(Form("hTPCnsigma%s%s%s", specnames[ispec].Data(), istep == 0 ? "Before" : "After", fBrowseCentrality == -1 ? "MinBias" : Form("Cent%d", fBrowseCentrality)));
bf892a6a 165 listNsigma->Add(hptr);
166 }
167 }
168
169 b->Add(listdEdx, "Projections dE/dx");
170 b->Add(listNsigma, "Projections NSigma");
171 }
172 }
173}
174
3a72645a 175//_________________________________________________________
176void AliHFEtpcPIDqa::Initialize(){
177 //
178 // Define Histograms
179 //
180
181 fHistos = new AliHFEcollection("tpcqahistos", "Collection of TPC QA histograms");
182
183 // Make common binning
a8d09079 184 const Int_t kNdim = 6;
3a72645a 185 const Int_t kPIDbins = AliPID::kSPECIES + 1;
3a72645a 186 const Int_t kSteps = 2;
6555e2ad 187 const Int_t kCentralityBins = 11;
a8d09079 188 const Int_t kEtabins=18;
3a72645a 189 const Double_t kMinPID = -1;
190 const Double_t kMinP = 0.;
191 const Double_t kMaxPID = (Double_t)AliPID::kSPECIES;
192 const Double_t kMaxP = 20.;
a8d09079 193 const Double_t kMinEta = -0.9;
194 const Double_t kMaxEta = 0.9;
195
bf892a6a 196 // Quantities where one can switch between low and high resolution
197 Int_t kPbins = fQAmanager->HasHighResolutionHistos() ? 1000 : 100;
198 Int_t kDedxbins = fQAmanager->HasHighResolutionHistos() ? 400 : 200;
199 Int_t kSigmaBins = fQAmanager->HasHighResolutionHistos() ? 1400 : 240;
200
3a72645a 201 // 1st histogram: TPC dEdx: (species, p, dEdx, step)
a8d09079 202 Int_t nBinsdEdx[kNdim] = {kPIDbins, kPbins, kDedxbins, kSteps, kCentralityBins, kEtabins};
203 Double_t mindEdx[kNdim] = {kMinPID, kMinP, 0., 0., 0., kMinEta};
204 Double_t maxdEdx[kNdim] = {kMaxPID, kMaxP, 200, 2., 11., kMaxEta};
205 fHistos->CreateTHnSparse("tpcDedx", "TPC signal; species; p [GeV/c]; TPC signal [a.u.]; Selection Step; Centrality; Eta", kNdim, nBinsdEdx, mindEdx, maxdEdx);
3a72645a 206 // 2nd histogram: TPC sigmas: (species, p nsigma, step)
a8d09079 207 Int_t nBinsSigma[kNdim] = {kPIDbins, kPbins, kSigmaBins, kSteps, kCentralityBins, kEtabins};
208 Double_t minSigma[kNdim] = {kMinPID, kMinP, -12., 0., 0., kMinEta};
209 Double_t maxSigma[kNdim] = {kMaxPID, kMaxP, 12., 2., 11., kMaxEta};
210 fHistos->CreateTHnSparse("tpcnSigma", "TPC signal; species; p [GeV/c]; TPC signal [a.u.]; Selection Step; Centrality; Eta", kNdim, nBinsSigma, minSigma, maxSigma);
3a72645a 211
212 // General TPC QA
213}
214
215//_________________________________________________________
6555e2ad 216void AliHFEtpcPIDqa::ProcessTrack(const AliHFEpidObject *track, AliHFEdetPIDqa::EStep_t step){
3a72645a 217 //
218 // Fill TPC histograms
219 //
220 AliDebug(1, Form("QA started for TPC PID for step %d", (Int_t)step));
6555e2ad 221 AliHFEpidObject::AnalysisType_t anatype = track->IsESDanalysis() ? AliHFEpidObject::kESDanalysis : AliHFEpidObject::kAODanalysis;
3a72645a 222 Float_t centrality = track->GetCentrality();
6555e2ad 223 Int_t species = track->GetAbInitioPID();
3a72645a 224 if(species >= AliPID::kSPECIES) species = -1;
3a72645a 225
6555e2ad 226 AliHFEpidTPC *tpcpid = dynamic_cast<AliHFEpidTPC *>(fQAmanager->GetDetectorPID(AliHFEpid::kTPCpid));
8c1c76e9 227 const AliPIDResponse *pidResponse = tpcpid ? tpcpid->GetPIDResponse() : NULL;
a8d09079 228 Double_t contentSignal[6];
3a72645a 229 contentSignal[0] = species;
bf892a6a 230 contentSignal[1] = tpcpid ? tpcpid->GetP(track->GetRecTrack(), anatype) : 0.;
6555e2ad 231 contentSignal[2] = GetTPCsignal(track->GetRecTrack(), anatype);
3a72645a 232 contentSignal[3] = step;
233 contentSignal[4] = centrality;
a8d09079 234 contentSignal[5] = GetEta(track->GetRecTrack(), anatype);
bf892a6a 235 fHistos->Fill("tpcDedx", contentSignal);
3a72645a 236
f06b970d 237 contentSignal[2] = pidResponse ? pidResponse->NumberOfSigmasTPC(track->GetRecTrack(), AliPID::kElectron) : -100.;
bf892a6a 238 fHistos->Fill("tpcnSigma", contentSignal);
3a72645a 239}
240
241//_________________________________________________________
6555e2ad 242Double_t AliHFEtpcPIDqa::GetTPCsignal(const AliVParticle *track, AliHFEpidObject::AnalysisType_t anatype){
3a72645a 243 //
6555e2ad 244 // Get TPC signal for ESD and AOD track
3a72645a 245 //
6555e2ad 246 Double_t tpcSignal = 0.;
247 if(anatype == AliHFEpidObject::kESDanalysis){
215ffe88 248 const AliESDtrack *esdtrack = static_cast<const AliESDtrack *>(track);
6555e2ad 249 tpcSignal = esdtrack->GetTPCsignal();
250 } else {
215ffe88 251 const AliAODTrack *aodtrack = static_cast<const AliAODTrack *>(track);
6555e2ad 252 tpcSignal = aodtrack->GetDetPid() ? aodtrack->GetDetPid()->GetTPCsignal() : 0.;
253 }
254 return tpcSignal;
3a72645a 255}
256
a8d09079 257
258//_________________________________________________________
259Double_t AliHFEtpcPIDqa::GetEta(const AliVParticle *track, AliHFEpidObject::AnalysisType_t anatype){
260 //
261 // Get TPC signal for ESD and AOD track
262 //
263 Double_t eta = 0.;
264 if(anatype == AliHFEpidObject::kESDanalysis){
265 const AliESDtrack *esdtrack = static_cast<const AliESDtrack *>(track);
266 eta = esdtrack->Eta();
267 } else {
268 const AliAODTrack *aodtrack = static_cast<const AliAODTrack *>(track);
269 eta = aodtrack->Eta();
270 }
271 return eta;
272}
273
3a72645a 274//_________________________________________________________
e156c3bb 275TH2 *AliHFEtpcPIDqa::MakeSpectrumdEdx(AliHFEdetPIDqa::EStep_t istep, Int_t species, Int_t centralityClass){
3a72645a 276 //
277 // Plot the Spectrum
278 //
279 THnSparseF *hSignal = dynamic_cast<THnSparseF *>(fHistos->Get("tpcDedx"));
bf892a6a 280 if(!hSignal) return NULL;
3a72645a 281 hSignal->GetAxis(3)->SetRange(istep + 1, istep + 1);
bf892a6a 282 if(species >= 0 && species < AliPID::kSPECIES)
3a72645a 283 hSignal->GetAxis(0)->SetRange(2 + species, 2 + species);
bf892a6a 284 TString hname = Form("hTPCsignal%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after"),
285 htitle = Form("TPC dE/dx Spectrum %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
3a72645a 286 if(species > -1){
bf892a6a 287 hname += AliPID::ParticleName(species);
288 htitle += Form(" for %ss", AliPID::ParticleName(species));
3a72645a 289 }
e156c3bb 290 TString centname, centtitle;
291 Bool_t hasCentralityInfo = kTRUE;
292 if(centralityClass > -1){
293 if(hSignal->GetNdimensions() < 5){
294 AliError("Centrality Information not available");
295 centname = centtitle = "MinBias";
296 hasCentralityInfo= kFALSE;
297 } else {
298 // Project centrality classes
299 // -1 is Min. Bias
300 hSignal->GetAxis(4)->SetRange(centralityClass+1, centralityClass+1);
301 centname = Form("Cent%d", centralityClass);
302 centtitle = Form(" Centrality %d", centralityClass);
303 }
304 } else {
305 centname = centtitle = "MinBias";
306 hasCentralityInfo= kFALSE;
307 }
308 hname += centtitle;
309 htitle += centtitle;
310
311 TH2 *hTmp = hSignal->Projection(2,1);
bf892a6a 312 hTmp->SetName(hname.Data());
313 hTmp->SetTitle(htitle.Data());
3a72645a 314 hTmp->SetStats(kFALSE);
315 hTmp->GetXaxis()->SetTitle("p [GeV/c]");
316 hTmp->GetYaxis()->SetTitle("TPC signal [a.u.]");
317 hSignal->GetAxis(3)->SetRange(0, hSignal->GetAxis(3)->GetNbins());
318 hSignal->GetAxis(0)->SetRange(0, hSignal->GetAxis(0)->GetNbins());
e156c3bb 319 if(hasCentralityInfo) hSignal->GetAxis(4)->SetRange(0, hSignal->GetAxis(4)->GetNbins());
3a72645a 320 return hTmp;
321}
322
323//_________________________________________________________
e156c3bb 324TH2 *AliHFEtpcPIDqa::MakeSpectrumNSigma(AliHFEdetPIDqa::EStep_t istep, Int_t species, Int_t centralityClass){
3a72645a 325 //
326 // Plot the Spectrum
327 //
328 THnSparseF *hSignal = dynamic_cast<THnSparseF *>(fHistos->Get("tpcnSigma"));
bf892a6a 329 if(!hSignal) return NULL;
3a72645a 330 hSignal->GetAxis(3)->SetRange(istep + 1, istep + 1);
331 if(species >= 0 && species < AliPID::kSPECIES)
332 hSignal->GetAxis(0)->SetRange(2 + species, 2 + species);
bf892a6a 333 TString hname = Form("hTPCsigma%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after"),
334 htitle = Form("TPC dE/dx Spectrum[#sigma] %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
3a72645a 335 if(species > -1){
bf892a6a 336 hname += AliPID::ParticleName(species);
337 htitle += Form(" for %ss", AliPID::ParticleName(species));
3a72645a 338 }
e156c3bb 339 TString centname, centtitle;
340 Bool_t hasCentralityInfo = kTRUE;
341 if(centralityClass > -1){
342 if(hSignal->GetNdimensions() < 5){
343 AliError("Centrality Information not available");
344 centname = centtitle = "MinBias";
345 hasCentralityInfo= kFALSE;
346 } else {
347 // Project centrality classes
348 // -1 is Min. Bias
349 hSignal->GetAxis(4)->SetRange(centralityClass+1, centralityClass+1);
350 centname = Form("Cent%d", centralityClass);
351 centtitle = Form(" Centrality %d", centralityClass);
352 }
353 } else {
354 centname = centtitle = "MinBias";
355 hasCentralityInfo= kFALSE;
356 }
357 hname += centtitle;
358 htitle += centtitle;
359
360 TH2 *hTmp = hSignal->Projection(2,1);
bf892a6a 361 hTmp->SetName(hname.Data());
362 hTmp->SetTitle(htitle.Data());
3a72645a 363 hTmp->SetStats(kFALSE);
364 hTmp->GetXaxis()->SetTitle("p [GeV/c]");
365 hTmp->GetYaxis()->SetTitle("TPC dE/dx - <dE/dx>|_{el} [#sigma]");
366 hSignal->GetAxis(3)->SetRange(0, hSignal->GetAxis(3)->GetNbins());
367 hSignal->GetAxis(0)->SetRange(0, hSignal->GetAxis(0)->GetNbins());
e156c3bb 368 if(hasCentralityInfo) hSignal->GetAxis(4)->SetRange(0, hSignal->GetAxis(4)->GetNbins());
3a72645a 369 return hTmp;
370}
371