]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWG2/FLOW/AliFlowTasks/AliFlowBayesianPID.h
coverty fix + dE/dx parameterization dependent on eta
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowTasks / AliFlowBayesianPID.h
1 #ifndef AliFlowBayesianPID_h
2 #define AliFlowBayesianPID_h
3
4 #include "TObject.h"
5 #include "AliESDpid.h"
6 #include "AliPIDResponse.h"
7
8 class TDatabasePDG;
9 class AliESDEvent;
10 class AliESDtrack;
11 class TH2D;
12 class TSpline3;
13 class TF1;
14 class AliTOFGeometry;
15 class AliTOFT0maker;
16
17 /*
18 HOW TO
19
20 1) in your main program:
21
22   AliFlowBayesianPID *mypid = new AliFlowBayesianPID();
23
24 or (if you want a global AliESDpid object)
25
26   AliESDpid *PIDesd = new AliESDpid();
27   ....
28   gROOT->LoadMacro("AliFlowBayesianPID.cxx++g");
29   AliFlowBayesianPID *mypid = new AliFlowBayesianPID(PIDesd);
30
31 for data starting from the PbPb pass2 reconstruction
32
33   mypid->SetNewTrackParam();
34
35
36 2) pass the pointer to your task
37
38 3) In your Task Exec:
39      mypid->SetDetResponse(esdEvent, centrality,AliESDpid::kTOF_T0,kFALSE); // centrality = PbPb centrality class (0-100%) or -1 for pp collisions
40      for(...){ // track loop
41      mypid->ComputeProb(track,centrality);
42      Float_t *prob = mypid->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mismatch parameterization
43
44      // for the single detector weights (no priors)
45      Float_t *tpcWeight = mypid->GetWeights(0); // TPC weights (equal weights in case of no kTPCpid)
46      Float_t *tofWeight = mypid->GetWeights(1); // TOF weights (equal weights in case of no kTOFpid)
47
48      // more tools
49      TSpline3 *mismatchShape = mypid->GetMismatch(); // time distribution of mismatched tracks - L/c (L is the distance of the pad from the PV)
50      TF1 *tpcprob = mypid->GetTPCprob(); // normalized TPC response function (in Nsigmas)  
51      TF1 *tofprob = mypid->GetTOFprob(); // normalized TOF response function (in Nsigmas)  
52
53      TH2D *hPr = mypid->GetHistoPriors(isp); // 2D (centrality - pT) histo for the priors of specie-isp (centrality < 0 means pp collisions)
54                                              // all the priors are normalized to the pion ones
55
56   }
57
58 */
59
60 class AliFlowBayesianPID : public AliPIDResponse {
61  public:
62   AliFlowBayesianPID(AliESDpid *esdpid=NULL); 
63   virtual ~AliFlowBayesianPID();
64   
65   // virtual method of AliPIDResponse
66   virtual Float_t NumberOfSigmasTOF(const AliVParticle *vtrack, AliPID::EParticleType type) const {if(vtrack) printf("Don't call AliFlowBayesianPID::NumberOfSigmasTOF method (%i)\n",type); return 0.0;} // do not use it
67
68   // setter
69   void SetMC(Bool_t flag = kTRUE){fIsMC=flag;} // actually do nothing
70   void SetDetResponse(AliESDEvent *esd,Float_t centrality=-1.0,EStartTimeType_t flagStart=AliESDpid::kTOF_T0,Bool_t recomputeT0TOF=kFALSE);
71   void SetNewTrackParam(Bool_t flag=kTRUE){fNewTrackParam=flag;};
72   void SetDetAND(Int_t idet){if(idet < fNdetectors && idet >= 0) fMaskAND[idet] = kTRUE;};
73   void SetDetOR(Int_t idet){if(idet < fNdetectors && idet >= 0) fMaskOR[idet] = kTRUE;};
74   void ResetDetAND(Int_t idet){if(idet < fNdetectors && idet >= 0) fMaskAND[idet] = kFALSE;};
75   void ResetDetOR(Int_t idet){if(idet < fNdetectors && idet >= 0) fMaskOR[idet] = kFALSE;};
76
77   // getter
78   AliESDpid* GetESDpid(){return fPIDesd;};
79   TH2D *GetHistoPriors(Int_t specie){if(specie >=0 && specie < fNspecies) return hPriors[specie]; else return NULL;};
80   TSpline3 *GetMismatch();  
81   TF1 *GetTOFprob(){return fTOFResponse;};
82   TF1 *GetTPCprob(){return fTPCResponse;};
83   Float_t *GetWeights(Int_t det){if(det < fNdetectors && det >= 0){return fWeights[det];} else{return NULL;}};
84   Float_t *GetProb(){return fProb;};
85   Float_t GetTOFMismWeight(){return fWTofMism;};
86   Float_t GetTOFMismProb(){return fProbTofMism;};
87   Float_t GetMassOverZ(){return fMassTOF;};
88   Float_t GetZ(){return fZ;};
89   Bool_t GetDetANDstatus(Int_t idet){if(idet < fNdetectors && idet >= 0){return fMaskAND[idet];} else{return kFALSE;} };
90   Bool_t GetDetORstatus(Int_t idet){if(idet < fNdetectors && idet >= 0){return fMaskOR[idet];} else{return kFALSE;} };
91   Bool_t GetCurrentMask(Int_t idet){if(idet < fNdetectors && idet >= 0){return fMaskCurrent[idet];} else{return kFALSE;} };
92   Float_t GetExpDeDx(const AliESDtrack *t,Int_t iS);
93
94   // methods for Bayesina Combined PID
95   void ComputeWeights(const AliESDtrack *t);
96   void ComputeProb(const AliESDtrack *t,Float_t); // obsolete method
97   void ComputeProb(const AliESDtrack *t){ComputeProb(t,0.0);}; 
98
99   void SetTOFres(Float_t res){fTOFres=res;};
100
101  private: 
102   void SetPriors();
103
104   static const Int_t fNdetectors = 2;
105   static const Int_t fNspecies = 8;// 0=el, 1=mu, 2=pi, 3=ka, 4=pr, 5=deuteron, 6=triton, 7=He3 
106   static TH2D* hPriors[fNspecies]; // histo with priors (hardcoded)
107   static TSpline3 *fMism; // function for mismatch
108   static AliTOFGeometry *fTofGeo; // TOF geometry needed to reproduce mismatch shape
109
110   AliESDpid *fPIDesd;//ESDpid object
111   TDatabasePDG *fDB; // Database pdg
112   Double_t fMass[fNspecies]; // mass for el(0),mu(1),pi(2),K(3),p(4)
113
114   Bool_t fNewTrackParam; // switch for new tracking resolution TOF parameterization
115   Bool_t fIsMC; // switch if MC data
116   Float_t fTOFres; // TOF res needed only if T0-TOF should be recomputed
117
118   AliFlowBayesianPID(const AliFlowBayesianPID&); // not implemented
119   AliFlowBayesianPID& operator=(const AliFlowBayesianPID&); // not implemented 
120
121   TF1 *fTOFResponse; // TOF Gaussian+tail response function (tail at 1.1 sigma)
122   TF1 *fTPCResponse; // TPC Gaussian+tail response function (tail at 1.8 sigma)
123
124   AliTOFT0maker *fTOFmaker; //TOF-T0 maker object
125
126   Float_t fWeights[fNdetectors][fNspecies]; // weights: 0=tpc,1=tof
127   Float_t fProb[fNspecies],fWTofMism,fProbTofMism; // Bayesian Combined PID + mismatch weights and probability 
128
129   Float_t fZ,fMassTOF; //measured charge(Z) and mass/Z
130   TF1 *fBBdata; // Bethe Bloch function (needed to compute the charge of the particle)
131
132   Bool_t fMaskAND[fNdetectors],fMaskOR[fNdetectors],fMaskCurrent[fNdetectors]; // mask detector should be used
133
134   Float_t fCurrCentrality; // Centrality in current event
135
136   ClassDef(AliFlowBayesianPID, 4); // example of analysis
137 };
138
139 #endif
140
141