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 **************************************************************************/
16 // QA class for Bayes PID
17 // Plot likelihoods vs various observables
18 // combined likelihood as well as indiv. detector response
19 // also plot mass spectrum using TOF for prior determination
22 // Yvonne Pachmayer <Y.Pachmayer@gsi.de>
27 #include <THnSparse.h>
30 #include "AliAODTrack.h"
31 #include "AliESDtrack.h"
34 #include "AliPIDResponse.h"
36 #include "AliHFEcollection.h"
37 #include "AliHFEpidBase.h"
38 #include "AliHFEpidQAmanager.h"
39 #include "AliHFEpidTPC.h"
40 #include "AliHFEpidTOF.h"
41 #include "AliHFEpidTRD.h"
42 #include "AliHFEpidBayes.h"
43 #include "AliHFEtools.h"
44 #include "AliHFEbayesPIDqa.h"
45 #include "AliHFEdetPIDqa.h"
47 ClassImp(AliHFEbayesPIDqa)
49 //_________________________________________________________
50 AliHFEbayesPIDqa::AliHFEbayesPIDqa():
59 //_________________________________________________________
60 AliHFEbayesPIDqa::AliHFEbayesPIDqa(const char* name):
61 AliHFEdetPIDqa(name, "QA for Bayes")
65 // Default constructor
69 //____________________________________________________________
70 AliHFEbayesPIDqa::AliHFEbayesPIDqa(const AliHFEbayesPIDqa &o):
79 //____________________________________________________________
80 AliHFEbayesPIDqa &AliHFEbayesPIDqa::operator=(const AliHFEbayesPIDqa &o){
85 AliHFEdetPIDqa::operator=(o);
92 //_________________________________________________________
93 AliHFEbayesPIDqa::~AliHFEbayesPIDqa(){
97 if(fHistos) delete fHistos;
99 //_________________________________________________________
100 void AliHFEbayesPIDqa::Copy(TObject &o) const {
104 AliHFEbayesPIDqa &target = dynamic_cast<AliHFEbayesPIDqa &>(o);
106 delete target.fHistos;
107 target.fHistos = NULL;
109 if(fHistos) target.fHistos = new AliHFEcollection(*fHistos);
113 //_________________________________________________________
114 Long64_t AliHFEbayesPIDqa::Merge(TCollection *coll){
116 // Merge with other objects
119 if(coll->IsEmpty()) return 1;
122 AliHFEbayesPIDqa *refQA = NULL;
127 refQA = dynamic_cast<AliHFEbayesPIDqa *>(o);
130 listHistos.Add(refQA->fHistos);
133 fHistos->Merge(&listHistos);
136 //_________________________________________________________
137 void AliHFEbayesPIDqa::Initialize(){
142 fHistos = new AliHFEcollection("BAYESqahistos", "Collection of Bayes QA histograms");
144 CreateDetectorSignalHistograms();
146 // fBAYESpid = new AliHFEpidBayes("QAbayesPID");
151 void AliHFEbayesPIDqa::CreateDetectorSignalHistograms(){
153 // Create Histogram for Probability Studies
156 Int_t kPbins = 1000; //fQAmanager->HasHighResolutionHistos() ? 1000 : 100;
157 Int_t kSigmaBins = 300; //fQAmanager->HasHighResolutionHistos() ? 600 : 150;
158 Int_t trdLikelihoodBins = 100; // fQAmanager->HasHighResolutionHistos() ? 200 : 100;
159 const Int_t kPIDbins = AliPID::kSPECIES + 1;
160 const Int_t kProbbins = 100;
161 const Int_t kSteps = 2;
162 const Int_t kCentralityBins = 11;
163 const Double_t kMinPID = -1;
164 const Double_t kMinP = 0.;
165 const Double_t kMaxPID = (Double_t)AliPID::kSPECIES;
166 const Double_t kMaxP = 10.;
168 Int_t nBinsContr[8] = {kPIDbins, kPbins, kProbbins, kSigmaBins, kSigmaBins, trdLikelihoodBins, kSteps, kCentralityBins};
169 Double_t minContr[8] = {kMinPID, kMinP, 0, -10., -7., 0., 0., 0.};
170 Double_t maxContr[8] = {kMaxPID, kMaxP, 1, 5., 8., 1., 2., 11.};
171 fHistos->CreateTHnSparse("control", "Control; species; p [GeV/c]; Comb prob; TPC sigma; TOF sigma; TRD electron Likelihood; selection step; centrality", 8, nBinsContr, minContr,maxContr);
175 Int_t nBinsTOFmass[5] = {kPIDbins, kPbins, 150, kSteps, kCentralityBins};
176 Double_t minTOFmass[5] = {kMinPID, kMinP, 0., 0., 0.};
177 Double_t maxTOFmass[5] = {kMaxPID, kMaxP, 1.5, 2., 11.};
178 fHistos->CreateTHnSparse("tofmass", "TOF mass; species; p [GeV/c]; TOF mass; selection step", 5, nBinsTOFmass, minTOFmass, maxTOFmass);
182 //_________________________________________________________
183 void AliHFEbayesPIDqa::ProcessTrack(const AliHFEpidObject *track, AliHFEdetPIDqa::EStep_t step){
185 // Fill TPC histograms
187 AliDebug(1, Form("QA started for BAYES PID for step %d", (Int_t)step));
188 AliHFEpidObject::AnalysisType_t anatype = track->IsESDanalysis() ? AliHFEpidObject::kESDanalysis : AliHFEpidObject::kAODanalysis;
189 Float_t centrality = track->GetCentrality();
190 Int_t species = track->GetAbInitioPID();
191 if(species >= AliPID::kSPECIES) species = -1;
193 Double_t probComb[AliPID::kSPECIES]={0.};
194 AliHFEpidBayes *bayespid = dynamic_cast<AliHFEpidBayes *>(fQAmanager->GetDetectorPID(AliHFEpid::kBAYESpid));
195 const AliPIDResponse *pidResponseBayes = bayespid ? bayespid->GetPIDResponse() : NULL;
196 if(!pidResponseBayes){
197 AliError("No PID Response available");
200 bayespid->CalcCombProb(track,pidResponseBayes, probComb);
202 Double_t contentSignal[8];
203 contentSignal[0] = species;
204 contentSignal[1] = track->GetRecTrack()->P();
205 contentSignal[2] = probComb[AliPID::kElectron];
206 // contentSignal[2] = contentSignal[2] = pidResponse ? pidResponse->NumberOfSigmasTPC(track->GetRecTrack(), AliPID::kElectron) : 0.;
208 AliHFEpidTOF *tofpid = dynamic_cast<AliHFEpidTOF *>(fQAmanager->GetDetectorPID(AliHFEpid::kTOFpid));
209 const AliPIDResponse *pidResponseTOF = tofpid ? tofpid->GetPIDResponse() : NULL;
211 AliError("No PID response available");
214 contentSignal[4] = pidResponseTOF->NumberOfSigmasTOF(track->GetRecTrack(), AliPID::kElectron);
216 AliHFEpidTPC *tpcpid = dynamic_cast<AliHFEpidTPC *>(fQAmanager->GetDetectorPID(AliHFEpid::kTPCpid));
217 const AliPIDResponse *pidResponse = tpcpid ? tpcpid->GetPIDResponse() : NULL;
219 AliError("No PID response available");
222 contentSignal[3] = pidResponse->NumberOfSigmasTPC(track->GetRecTrack(), AliPID::kElectron);
224 AliHFEpidTRD *trdpid = dynamic_cast<AliHFEpidTRD *>(fQAmanager->GetDetectorPID(AliHFEpid::kTRDpid));
225 contentSignal[5] = trdpid ? trdpid->GetElectronLikelihood(static_cast<const AliVTrack*>(track->GetRecTrack()), anatype) : -10;
227 contentSignal[6] = step;
228 contentSignal[7] = centrality;
229 fHistos->Fill("control", contentSignal);
231 Double_t contentFill[5];
232 contentFill[0]=contentSignal[0];
233 contentFill[1]=contentSignal[1];
234 contentFill[3]=contentSignal[6];
235 contentFill[4]=contentSignal[7];
237 Double_t masscalcfromtof=CalcTOFMass(track);
238 contentFill[2]=masscalcfromtof;
239 fHistos->Fill("tofmass", contentFill);
242 Double_t AliHFEbayesPIDqa::CalcTOFMass(const AliHFEpidObject *track){
247 Double_t mass=-99; //GeV
248 Double_t length=((AliESDtrack*)track->GetRecTrack())->GetIntegratedLength();
249 if (length<=0) return mass;
250 Double_t tofTime=((AliESDtrack*)track->GetRecTrack())->GetTOFsignal();//in ps
251 Double_t tof= tofTime*1E-3; // ns, average T0 fill subtracted, no info from T0detector
252 if (tof<=0) return mass;
253 Double_t c=TMath::C()*1.E-9;// m/ns
254 length =length*0.01; // in meters
256 Double_t fact= (tof/length)*(tof/length) -1.;
257 Double_t mom=track->GetRecTrack()->P();
258 if(mom==0) return mass;
261 mass = -mom*TMath::Sqrt(-fact);
263 mass = mom*TMath::Sqrt(fact);