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 **************************************************************************/
17 // Interface to the user task
18 // Several strategies for Electron identification implemented.
19 // In addition users can combine different detectors or use
20 // single detector PID
23 // Markus Fasel <M.Fasel@gsi.de>
28 #include <TIterator.h>
30 #include <TObjArray.h>
31 #include <TObjString.h>
34 #include "AliAODpidUtil.h"
35 #include "AliESDpid.h"
38 #include "AliVParticle.h"
40 #include "AliHFEcontainer.h"
41 #include "AliHFEpid.h"
42 #include "AliHFEpidQAmanager.h"
43 #include "AliHFEpidITS.h"
44 #include "AliHFEpidTPC.h"
45 #include "AliHFEpidTRD.h"
46 #include "AliHFEpidTOF.h"
47 #include "AliHFEpidEMCAL.h"
48 #include "AliHFEpidMC.h"
49 #include "AliHFEvarManager.h"
53 const Char_t* AliHFEpid::fgkDetectorName[AliHFEpid::kNdetectorPID + 1] = {
64 //____________________________________________________________
65 AliHFEpid::AliHFEpid():
73 // Default constructor
75 memset(fDetectorPID, 0, sizeof(AliHFEpidBase *) * kNdetectorPID);
76 memset(fDetectorOrder, kUndefined, sizeof(UInt_t) * kNdetectorPID);
77 memset(fDetectorOrder, 0, sizeof(UInt_t) * kNdetectorPID);
80 //____________________________________________________________
81 AliHFEpid::AliHFEpid(const Char_t *name):
89 // Default constructor
90 // Create PID objects for all detectors
92 memset(fDetectorPID, 0, sizeof(AliHFEpidBase *) * kNdetectorPID);
93 memset(fDetectorOrder, kUndefined, sizeof(UInt_t) * kNdetectorPID);
94 memset(fSortedOrder, 0, sizeof(UInt_t) * kNdetectorPID);
96 fDetectorPID[kMCpid] = new AliHFEpidMC("MCPID");
97 fDetectorPID[kTPCpid] = new AliHFEpidTPC("TPCPID");
98 fDetectorPID[kTRDpid] = new AliHFEpidTRD("TRDPID");
99 fDetectorPID[kTOFpid] = new AliHFEpidTOF("TOFPID");
100 fDetectorPID[kEMCALpid] = new AliHFEpidEMCAL("EMCALPID");
104 //____________________________________________________________
105 AliHFEpid::AliHFEpid(const AliHFEpid &c):
107 fEnabledDetectors(c.fEnabledDetectors),
108 fNPIDdetectors(c.fNPIDdetectors),
109 fVarManager(c.fVarManager),
115 memset(fDetectorPID, 0, sizeof(AliHFEpidBase *) * kNdetectorPID);
119 //____________________________________________________________
120 AliHFEpid& AliHFEpid::operator=(const AliHFEpid &c){
122 // Assignment operator
124 if(&c != this) c.Copy(*this);
128 //____________________________________________________________
129 AliHFEpid::~AliHFEpid(){
133 for(Int_t idet = 0; idet < kNdetectorPID; idet++)
134 if(fDetectorPID[idet]) delete fDetectorPID[idet];
135 ClearCommonObjects();
138 //____________________________________________________________
139 void AliHFEpid::Copy(TObject &o) const{
145 AliHFEpid &target = dynamic_cast<AliHFEpid &>(o);
146 target.ClearCommonObjects();
148 target.fEnabledDetectors = fEnabledDetectors;
149 target.fNPIDdetectors = fNPIDdetectors;
150 target.fVarManager = fVarManager;
152 // Copy detector PIDs
153 for(Int_t idet = 0; idet < kNdetectorPID; idet++){
154 //Cleanup pointers in case of assignment
155 if(target.fDetectorPID[idet])
156 delete target.fDetectorPID[idet];
157 if(fDetectorPID[idet])
158 target.fDetectorPID[idet] = dynamic_cast<AliHFEpidBase *>(fDetectorPID[idet]->Clone());
160 memcpy(target.fDetectorOrder, fDetectorOrder, sizeof(UInt_t) * kNdetectorPID);
161 memcpy(target.fSortedOrder, fSortedOrder, sizeof(UInt_t) * kNdetectorPID);
164 //____________________________________________________________
165 void AliHFEpid::AddCommonObject(TObject * const o){
167 // Add common object to the garbage collection
169 if(!fCommonObjects) fCommonObjects = new TObjArray;
170 fCommonObjects->Add(o);
173 //____________________________________________________________
174 void AliHFEpid::ClearCommonObjects(){
176 // empty garbage collection
179 fCommonObjects->Delete();
180 delete fCommonObjects;
181 fCommonObjects = NULL;
185 //____________________________________________________________
186 void AliHFEpid::AddDetector(TString detector, UInt_t position){
188 // Add Detector in position
190 UInt_t detectorID = kUndefined;
192 if(!detector.CompareTo("MC")) detectorID = kMCpid;
193 else if(!detector.CompareTo("TPC")) detectorID = kTPCpid;
194 else if(!detector.CompareTo("TRD")) detectorID = kTRDpid;
195 else if(!detector.CompareTo("TOF")) detectorID = kTOFpid;
196 else if(!detector.CompareTo("EMCAL")) detectorID = kEMCALpid;
197 else AliError("Detector not available");
199 if(detectorID == kUndefined) return;
200 if(IsDetectorOn(detectorID)) return;
201 SwitchOnDetector(detectorID);
202 fDetectorOrder[detectorID] = position;
206 //____________________________________________________________
207 Bool_t AliHFEpid::InitializePID(){
209 // Initializes the PID object
212 TMath::Sort(static_cast<UInt_t>(kNdetectorPID), fDetectorOrder, fSortedOrder, kFALSE);
213 // Initialize PID Objects
214 Bool_t status = kTRUE;
215 for(Int_t idet = 0; idet < kNdetectorPID; idet++){
216 if(!IsDetectorOn(idet)) continue;
217 if(fDetectorPID[idet]){
218 status &= fDetectorPID[idet]->InitializePID();
219 if(HasMCData() && status) fDetectorPID[idet]->SetHasMCData();
226 //____________________________________________________________
227 Bool_t AliHFEpid::IsSelected(const AliHFEpidObject * const track, AliHFEcontainer *cont, const Char_t *contname, AliHFEpidQAmanager *pidqa){
231 Bool_t isSelected = kTRUE;
232 AliDebug(1, Form("Particle used for PID, QA available: %s", pidqa ? "Yes" : "No"));
233 for(UInt_t idet = 0; idet < fNPIDdetectors; idet++){
234 AliDebug(2, Form("Using Detector %s\n", SortedDetectorName(idet)));
235 if(TMath::Abs(fDetectorPID[fSortedOrder[idet]]->IsSelected(track, pidqa)) != 11){
239 AliDebug(2, "Particlae selected by detector");
240 if(fVarManager && cont){
241 TString reccontname = contname; reccontname += "Reco";
242 AliDebug(2, Form("Filling container %s", reccontname.Data()));
243 if(fVarManager->IsSignalTrack())
244 fVarManager->FillContainerStepname(cont, reccontname.Data(), SortedDetectorName(idet));
246 TString mccontname = contname; mccontname += "MC";
247 AliDebug(2, Form("MC Information available, Filling container %s", mccontname.Data()));
248 if(fVarManager->IsSignalTrack()) {
249 fVarManager->FillContainerStepname(cont, mccontname.Data(), SortedDetectorName(idet), kTRUE);
250 if(cont->GetCorrelationMatrix("correlationstepafterTOF")){
251 TString tstept("TOFPID");
252 if(!tstept.CompareTo(SortedDetectorName(idet))) {
253 fVarManager->FillCorrelationMatrix(cont->GetCorrelationMatrix("correlationstepafterTOF"));
254 //printf("Step %s\n",(const char*) SortedDetectorName(idet));
259 // The PID will NOT fill the double counting information
265 //____________________________________________________________
266 void AliHFEpid::SetESDpid(AliESDpid *pid){
268 // Set ESD PID to the Detector PID objects
270 for(Int_t idet = 0; idet < kNdetectorPID; idet++){
271 if(fDetectorPID[idet]) fDetectorPID[idet]->SetESDpid(pid);
275 //____________________________________________________________
276 void AliHFEpid::SetAODpid(AliAODpidUtil *pid){
278 // Set ESD PID to the Detector PID objects
280 for(Int_t idet = 0; idet < kNdetectorPID; idet++){
281 if(fDetectorPID[idet]) fDetectorPID[idet]->SetAODpid(pid);
285 //____________________________________________________________
286 void AliHFEpid::ConfigureTPCasymmetric(Double_t pmin, Double_t pmax, Double_t sigmamin, Double_t sigmamax){
288 // TPC alone, symmetric 3 sigma cut and asymmetric sigma cut in the momentum region between 2GeV/c and 10 GeV/c and sigma between -1 and 100
290 AliHFEpidTPC *pid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
292 pid->SetTPCnSigma(3);
293 pid->SetAsymmetricTPCsigmaCut(pmin, pmax, sigmamin, sigmamax);
297 //____________________________________________________________
298 void AliHFEpid::ConfigureTPCrejectionSimple(){
300 // TPC alone, symmetric 3 sigma cut and 2 - -100 sigma pion rejection
302 AliHFEpidTPC *pid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
304 pid->SetTPCnSigma(3);
305 pid->SetRejectParticle(AliPID::kPion, 0., -100., 10., 1.);
309 //____________________________________________________________
310 void AliHFEpid::ConfigureTOF(Float_t TOFCut){
312 // Set Number of sigmas for TOF PID
314 AliHFEpidTOF *tofpid = dynamic_cast<AliHFEpidTOF *>(fDetectorPID[kTOFpid]);
315 if(tofpid) tofpid->SetTOFnSigma(TOFCut);
318 //____________________________________________________________
319 void AliHFEpid::ConfigureTPCcentralityCut(Int_t centralityBin, const char *lowerCutParam, const Double_t * const params, Float_t upperTPCCut){
321 // Cofigure centrality dependent cut function for TPC PID
323 ConfigureTPCcut(centralityBin, lowerCutParam, params, upperTPCCut);
326 //____________________________________________________________
327 void AliHFEpid::ConfigureTPCdefaultCut(const char *lowerCutParam, const Double_t * const params, Float_t upperTPCCut){
329 // Cofigure default cut function for TPC PID
331 ConfigureTPCcut(-1, lowerCutParam, params, upperTPCCut);
334 //____________________________________________________________
335 void AliHFEpid::ConfigureTPCcut(Int_t centralityBin, const char *lowerCutParam, const Double_t * const params, Float_t upperTPCCut){
337 // Cofigure cut function for TPC PID
338 // if no function parameterizaion is given, then the default one (exponential) is chosen
341 if(HasMCData()) AliInfo("Configuring TPC for MC\n");
342 AliHFEpidTPC *tpcpid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
343 //TF1 *upperCut = new TF1("upperCut", "[0] * TMath::Exp([1]*x)", 0, 20);
344 TF1 *upperCut = new TF1(Form("upperCut%s", centralityBin < 0 ? "Default" : Form("Bin%d", centralityBin)), "[0]", 0, 20); // Use constant upper cut
345 TF1 *lowerCut = new TF1(Form("lowerCut%s", centralityBin < 0 ? "Default" : Form("Bin%d", centralityBin)), lowerCutParam == NULL ? "[0] * TMath::Exp([1]*x) + [2]": lowerCutParam, 0, 20);
347 upperCut->SetParameter(0, upperTPCCut); // pp
350 for(Int_t ipar = 0; ipar < lowerCut->GetNpar(); ipar++)
352 lowerCut->SetParameter(ipar, params[ipar]);
353 // printf("printout %i %s %f \n", centralityBin, lowerCutParam, params[ipar]);
356 // Set default parameterization
357 if(HasMCData()) lowerCut->SetParameter(0, -2.5);
358 else lowerCut->SetParameter(0, -4.03); //pp
359 lowerCut->SetParameter(1, -0.22); // pp
361 if(HasMCData()) lowerCut->SetParameter(2, -2.2);
362 else lowerCut->SetParameter(2, 0.92); //pp
367 tpcpid->SetTPCnSigma(2);
368 if(centralityBin < 0){
369 tpcpid->SetUpperSigmaCutDefault(upperCut);
370 tpcpid->SetLowerSigmaCutDefault(lowerCut);
372 tpcpid->SetUpperSigmaCutCentrality(upperCut, centralityBin);
373 tpcpid->SetLowerSigmaCutCentrality(lowerCut, centralityBin);
376 AddCommonObject(upperCut);
377 AddCommonObject(lowerCut);
380 //____________________________________________________________
381 void AliHFEpid::PrintStatus() const {
383 // Print the PID configuration
385 printf("\n%s: Printing configuration\n", GetName());
386 printf("===============================================\n");
387 printf("PID Detectors: \n");
389 TString detectors[kNdetectorPID] = {"MC", "ESD", "ITS", "TPC", "TRD", "TOF", "EMCAL"};
390 for(Int_t idet = 0; idet < kNdetectorPID; idet++){
391 if(IsDetectorOn(idet)){
392 printf("\t%s\n", detectors[idet].Data());
396 if(!npid) printf("\tNone\n");