1 #include "AliSpectraBothPID.h"
2 #include "AliAODEvent.h"
6 #include "AliAODTrack.h"
7 #include "AliAODMCParticle.h"
8 #include "AliPIDResponse.h"
9 #include "AliAnalysisManager.h"
10 #include "AliInputEventHandler.h"
11 #include "AliSpectraBothTrackCuts.h"
13 ClassImp(AliSpectraBothPID)
15 AliSpectraBothPID::AliSpectraBothPID() : TNamed("PID", "PID object"), fPIDType(kNSigmacircleTPCTOF), fNSigmaPID(3), fPIDResponse(0),fshiftTPC(0),fshiftTOF(0)
20 AliSpectraBothPID::AliSpectraBothPID(BothPIDType_t pidType) : TNamed("PID", "PID object"), fPIDType(pidType), fNSigmaPID(3), fPIDResponse(0), fshiftTPC(0),fshiftTOF(0)
29 void AliSpectraBothPID::FillQAHistos(AliSpectraBothHistoManager * hman, AliVTrack * track, AliSpectraBothTrackCuts * trackCuts)
32 // fill a bunch of QA histos
35 // Get PID response object, if needed
38 AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
39 AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
40 fPIDResponse = inputHandler->GetPIDResponse();
42 Bool_t ycut[3]={false,false,false};
43 ycut[0]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpPion);
44 ycut[1]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpKaon);
45 ycut[2]=trackCuts->CheckYCut ((BothParticleSpecies_t) kSpProton);
49 AliVParticle *inEvHMain = dynamic_cast<AliVParticle *>(track);
51 hman->GetPIDHistogram(kHistPIDTPC)->Fill(track->GetTPCmomentum(), track->GetTPCsignal()*track->Charge()); // PID histo
52 hman->GetPIDHistogram(kHistPIDTPCPion)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kPion)*track->Charge()); // Expected PIDPion histo
53 hman->GetPIDHistogram(kHistPIDTPCKaon)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kKaon)*track->Charge()); // Expected PIDKaon histo
54 hman->GetPIDHistogram(kHistPIDTPCProton)->Fill(track->GetTPCmomentum(),fPIDResponse->GetTPCResponse().GetExpectedSignal(track->GetTPCmomentum(),AliPID::kProton)*track->Charge()); // Expected PIDProton histo
56 Double_t nsigmaTPCkProton = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton)+fshiftTPC;
57 Double_t nsigmaTPCkKaon = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)+fshiftTPC;
58 Double_t nsigmaTPCkPion = fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)+fshiftTPC;
59 Double_t nsigmaTOFkProton=100.0,nsigmaTOFkKaon=100.0,nsigmaTOFkPion=100.0;
61 Double_t nsigmaTPCTOFkProton = TMath::Abs(nsigmaTPCkProton);
62 Double_t nsigmaTPCTOFkKaon = TMath::Abs(nsigmaTPCkKaon);
63 Double_t nsigmaTPCTOFkPion = TMath::Abs(nsigmaTPCkPion);
64 if(fPIDType==kNSigmaTOF)
66 nsigmaTPCTOFkProton = 100.0;
67 nsigmaTPCTOFkKaon = 100.0;
68 nsigmaTPCTOFkPion =100.0;
72 if(track->Pt()>trackCuts->GetPtTOFMatching())
75 hman->GetPIDHistogram(kHistPIDTOF)->Fill(track->P(),(track->GetTOFsignal()/100)*track->Charge()); // PID histo
77 nsigmaTOFkProton = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF;
78 nsigmaTOFkKaon = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF;
79 nsigmaTOFkPion = fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF;
84 hman->GetPtHistogram(kHistNSigPionTOF)->Fill(track->P(),nsigmaTOFkPion );
85 hman->GetPtHistogram(kHistNSigPionPtTOF)->Fill(track->Pt(),nsigmaTOFkPion );
89 hman->GetPtHistogram(kHistNSigKaonTOF)->Fill(track->P(),nsigmaTOFkKaon);
90 hman->GetPtHistogram(kHistNSigKaonPtTOF)->Fill(track->Pt(),nsigmaTOFkKaon );
94 hman->GetPtHistogram(kHistNSigProtonTOF)->Fill(track->P(),nsigmaTOFkProton );
95 hman->GetPtHistogram(kHistNSigProtonPtTOF)->Fill(track->Pt(),nsigmaTOFkProton );
99 nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProton*nsigmaTPCkProton+nsigmaTOFkProton*nsigmaTOFkProton)/2.0);
100 nsigmaTPCTOFkKaon = TMath::Sqrt((nsigmaTPCkKaon*nsigmaTPCkKaon+nsigmaTOFkKaon*nsigmaTOFkKaon)/2.0);
101 nsigmaTPCTOFkPion = TMath::Sqrt((nsigmaTPCkPion*nsigmaTPCkPion+nsigmaTOFkPion*nsigmaTOFkPion)/2.0);
107 hman->GetPtHistogram(kHistNSigPionTPC)->Fill(track->P(),nsigmaTPCkPion );
108 hman->GetPtHistogram(kHistNSigPionPtTPC)->Fill(track->Pt(),nsigmaTPCkPion );
109 hman->GetPtHistogram(kHistNSigPionTPCTOF)->Fill(track->P(),nsigmaTPCTOFkPion);
110 hman->GetPtHistogram(kHistNSigPionPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkPion);
114 hman->GetPtHistogram(kHistNSigKaonTPC)->Fill(track->P(),nsigmaTPCkKaon );
115 hman->GetPtHistogram(kHistNSigKaonPtTPC)->Fill(track->Pt(),nsigmaTPCkKaon );
116 hman->GetPtHistogram(kHistNSigKaonTPCTOF)->Fill(track->P(),nsigmaTPCTOFkKaon);
117 hman->GetPtHistogram(kHistNSigKaonPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkKaon);
122 hman->GetPtHistogram(kHistNSigProtonTPC)->Fill(track->P(),nsigmaTPCkProton );
123 hman->GetPtHistogram(kHistNSigProtonPtTPC)->Fill(track->Pt(),nsigmaTPCkProton );
124 hman->GetPtHistogram(kHistNSigProtonTPCTOF)->Fill(track->P(),nsigmaTPCTOFkProton);
125 hman->GetPtHistogram(kHistNSigProtonPtTPCTOF)->Fill(track->Pt(),nsigmaTPCTOFkProton);
130 Int_t AliSpectraBothPID::GetParticleSpecie(AliAODMCParticle * part)
132 // return PID according to MC truth
133 switch(TMath::Abs(part->PdgCode()))
149 Int_t AliSpectraBothPID::GetParticleSpecie(TParticle* part)
151 // return PID according to MC truth
152 switch(TMath::Abs(part->GetPdgCode()))
168 Int_t AliSpectraBothPID::GetParticleSpecie(AliSpectraBothHistoManager * hman,AliVTrack * trk, AliSpectraBothTrackCuts * trackCuts, Bool_t* rec)
170 // return PID according to detectors
171 // Get PID response object, if needed
173 // guess the particle based on the smaller nsigma
176 rec[kSpProton]=false;
180 AliAnalysisManager *man = AliAnalysisManager::GetAnalysisManager();
181 AliInputEventHandler* inputHandler = (AliInputEventHandler*)(man->GetInputEventHandler());
182 fPIDResponse = inputHandler->GetPIDResponse();
187 AliFatal("Cannot get pid response");
192 // Compute nsigma for each hypthesis
193 AliVParticle *inEvHMain = dynamic_cast<AliVParticle *>(trk);
196 Double_t nsigmaTPCkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kProton)+fshiftTPC);
197 Double_t nsigmaTPCkKaon = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kKaon)+fshiftTPC);
198 Double_t nsigmaTPCkPion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(inEvHMain, AliPID::kPion)+fshiftTPC);
202 Double_t nsigmaTOFkProton=0.0,nsigmaTOFkKaon=0.0,nsigmaTOFkPion=0.0;
206 Double_t nsigmaTPCTOFkProton = TMath::Abs(nsigmaTPCkProton);
207 Double_t nsigmaTPCTOFkKaon = TMath::Abs(nsigmaTPCkKaon);
208 Double_t nsigmaTPCTOFkPion = TMath::Abs(nsigmaTPCkPion);
209 if(fPIDType==kNSigmaTOF)
211 nsigmaTPCTOFkProton = 100.0;
212 nsigmaTPCTOFkKaon = 100.0;
213 nsigmaTPCTOFkPion =100.0;
214 nsigmaTOFkProton = 100.0;
215 nsigmaTOFkKaon = 100.0;
216 nsigmaTOFkPion =100.0;
222 if(trk->Pt()>trackCuts->GetPtTOFMatching())
224 nsigmaTOFkProton = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kProton)+fshiftTOF);
225 nsigmaTOFkKaon = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kKaon)+fshiftTOF);
226 nsigmaTOFkPion = TMath::Abs(fPIDResponse->NumberOfSigmasTOF(inEvHMain, AliPID::kPion)+fshiftTOF);
227 // the TOF info is used in combined
228 nsigmaTPCTOFkProton = TMath::Sqrt((nsigmaTPCkProton*nsigmaTPCkProton+nsigmaTOFkProton*nsigmaTOFkProton)/2.);
229 nsigmaTPCTOFkKaon = TMath::Sqrt((nsigmaTPCkKaon*nsigmaTPCkKaon+nsigmaTOFkKaon*nsigmaTOFkKaon)/2.);
230 nsigmaTPCTOFkPion = TMath::Sqrt((nsigmaTPCkPion*nsigmaTPCkPion+nsigmaTOFkPion*nsigmaTOFkPion)/2.);
234 if (fPIDType==kNSigmaTOF)
239 // select the nsigma to be used for the actual PID
240 Double_t nsigmaPion=100;
241 Double_t nsigmaKaon=100;
242 Double_t nsigmaProton=100;
247 nsigmaProton = nsigmaTPCkProton;
248 nsigmaKaon = nsigmaTPCkKaon ;
249 nsigmaPion = nsigmaTPCkPion ;
252 nsigmaProton = nsigmaTOFkProton;
253 nsigmaKaon = nsigmaTOFkKaon ;
254 nsigmaPion = nsigmaTOFkPion ;
256 case kNSigmacircleTPCTOF:
257 nsigmaProton = nsigmaTPCTOFkProton;
258 nsigmaKaon = nsigmaTPCTOFkKaon ;
259 nsigmaPion = nsigmaTPCTOFkPion ;
261 case kNSigmasquareTPCTOF:
262 nsigmaProton = TMath::Max(nsigmaTPCkProton,nsigmaTOFkProton);
263 nsigmaKaon = TMath::Max(nsigmaTPCkKaon,nsigmaTOFkKaon);
264 nsigmaPion = TMath::Max(nsigmaTPCkPion,nsigmaTOFkPion) ;
269 if(nsigmaPion < fNSigmaPID)
271 if(nsigmaKaon < fNSigmaPID)
273 if(nsigmaProton < fNSigmaPID)
276 if( ( nsigmaKaon==nsigmaPion ) && ( nsigmaKaon==nsigmaProton ))
277 return kSpUndefined;//if is the default value for the three
278 if( ( nsigmaKaon < nsigmaPion ) && ( nsigmaKaon < nsigmaProton ) && (nsigmaKaon < fNSigmaPID))
280 if(hman->GetPIDHistogram(kHistPIDTPCKaonRec))
281 hman->GetPIDHistogram(kHistPIDTPCKaonRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDKaon histo
284 if( ( nsigmaPion < nsigmaKaon ) && ( nsigmaPion < nsigmaProton ) && (nsigmaPion < fNSigmaPID))
286 if(hman->GetPIDHistogram(kHistPIDTPCPionRec))
287 hman->GetPIDHistogram(kHistPIDTPCPionRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDPion histo
290 if( ( nsigmaProton < nsigmaKaon ) && ( nsigmaProton < nsigmaPion ) && (nsigmaProton < fNSigmaPID))
292 if(hman->GetPIDHistogram(kHistPIDTPCProtonRec))
293 hman->GetPIDHistogram(kHistPIDTPCProtonRec)->Fill(trk->GetTPCmomentum(), trk->GetTPCsignal()*trk->Charge()); // Reconstructed PIDProton histo
296 // else, return undefined
301 Long64_t AliSpectraBothPID::Merge(TCollection* list)
303 // Merging interface.
304 // Returns the number of merged objects (including this).
314 TIterator* iter = list->MakeIterator();
317 // Actually, we don't do anything here...
318 // collections of all histograms
319 // TList collections;
323 while ((obj = iter->Next())) {
324 AliSpectraBothPID* entry = dynamic_cast<AliSpectraBothPID*> (obj);
328 // TH1I * histo = entry->GetHistoCuts();
329 // collections.Add(histo);
333 // fHistoCuts->Merge(&collections);