]>
Commit | Line | Data |
---|---|---|
c630aafd | 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 | //----------------------------------------------------------------- | |
17 | // Implementation of the TOF PID class | |
18 | // Very naive one... Should be made better by the detector experts... | |
19 | // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch | |
20 | //----------------------------------------------------------------- | |
21 | #include "TFile.h" | |
22 | #include "TTree.h" | |
23 | #include "TClonesArray.h" | |
3f83f224 | 24 | #include "TError.h" |
c630aafd | 25 | |
26 | #include "AliTOFpidESD.h" | |
27 | #include "AliESD.h" | |
28 | #include "AliESDtrack.h" | |
29 | #include "AliTOFdigit.h" | |
8c29135e | 30 | #include "AliTOFGeometry.h" |
c630aafd | 31 | |
3f83f224 | 32 | #include <stdlib.h> |
33 | ||
c630aafd | 34 | |
35 | ClassImp(AliTOFpidESD) | |
36 | ||
0f4a7374 | 37 | //_________________________________________________________________________ |
8c29135e | 38 | AliTOFpidESD::AliTOFpidESD(Double_t *param) { |
0f4a7374 | 39 | // |
40 | // The main constructor | |
41 | // | |
0f4a7374 | 42 | fR=378.; |
8c29135e | 43 | fDy=AliTOFGeometry::XPad(); fDz=AliTOFGeometry::ZPad(); |
0f4a7374 | 44 | fN=0; fEventN=0; |
45 | ||
46 | fSigma=param[0]; | |
47 | fRange=param[1]; | |
48 | ||
49 | } | |
50 | ||
0f4a7374 | 51 | //_________________________________________________________________________ |
8c29135e | 52 | Int_t AliTOFpidESD::LoadClusters(TTree *dTree, AliTOFGeometry *geom) { |
c630aafd | 53 | //-------------------------------------------------------------------- |
54 | //This function loads the TOF clusters | |
55 | //-------------------------------------------------------------------- | |
56 | TBranch *branch=dTree->GetBranch("TOF"); | |
57 | if (!branch) { | |
58 | Error("LoadClusters"," can't get the branch with the TOF digits !\n"); | |
59 | return 1; | |
60 | } | |
61 | ||
62 | TClonesArray dummy("AliTOFdigit",10000), *digits=&dummy; | |
63 | branch->SetAddress(&digits); | |
64 | ||
65 | dTree->GetEvent(0); | |
66 | Int_t nd=digits->GetEntriesFast(); | |
67 | Info("LoadClusters","number of digits: %d",nd); | |
68 | ||
69 | for (Int_t i=0; i<nd; i++) { | |
70 | AliTOFdigit *d=(AliTOFdigit*)digits->UncheckedAt(i); | |
71 | Int_t dig[5]; Float_t g[3]; | |
72 | dig[0]=d->GetSector(); | |
73 | dig[1]=d->GetPlate(); | |
74 | dig[2]=d->GetStrip(); | |
75 | dig[3]=d->GetPadz(); | |
76 | dig[4]=d->GetPadx(); | |
77 | ||
8c29135e | 78 | geom->GetPos(dig,g); |
c630aafd | 79 | |
80 | Double_t h[5]; | |
81 | h[0]=TMath::Sqrt(g[0]*g[0]+g[1]*g[1]); | |
82 | h[1]=TMath::ATan2(g[1],g[0]); h[2]=g[2]; | |
83 | h[3]=d->GetTdc(); h[4]=d->GetAdc(); | |
84 | ||
3f83f224 | 85 | AliTOFcluster *cl=new AliTOFcluster(h,d->GetTracks(),i); |
c630aafd | 86 | InsertCluster(cl); |
87 | } | |
88 | ||
89 | return 0; | |
90 | } | |
91 | ||
0f4a7374 | 92 | //_________________________________________________________________________ |
c630aafd | 93 | void AliTOFpidESD::UnloadClusters() { |
94 | //-------------------------------------------------------------------- | |
95 | //This function unloads TOF clusters | |
96 | //-------------------------------------------------------------------- | |
97 | for (Int_t i=0; i<fN; i++) delete fClusters[i]; | |
98 | fN=0; | |
99 | } | |
100 | ||
0f4a7374 | 101 | //_________________________________________________________________________ |
c630aafd | 102 | Int_t AliTOFpidESD::InsertCluster(AliTOFcluster *c) { |
103 | //-------------------------------------------------------------------- | |
104 | //This function adds a cluster to the array of clusters sorted in Z | |
105 | //-------------------------------------------------------------------- | |
106 | if (fN==kMaxCluster) { | |
107 | Error("InsertCluster","Too many clusters !\n"); | |
108 | return 1; | |
109 | } | |
110 | ||
111 | if (fN==0) {fClusters[fN++]=c; return 0;} | |
112 | Int_t i=FindClusterIndex(c->GetZ()); | |
113 | memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliTOFcluster*)); | |
114 | fClusters[i]=c; fN++; | |
115 | ||
116 | return 0; | |
117 | } | |
118 | ||
0f4a7374 | 119 | //_________________________________________________________________________ |
c630aafd | 120 | Int_t AliTOFpidESD::FindClusterIndex(Double_t z) const { |
121 | //-------------------------------------------------------------------- | |
122 | // This function returns the index of the nearest cluster | |
123 | //-------------------------------------------------------------------- | |
124 | if (fN==0) return 0; | |
125 | if (z <= fClusters[0]->GetZ()) return 0; | |
126 | if (z > fClusters[fN-1]->GetZ()) return fN; | |
127 | Int_t b=0, e=fN-1, m=(b+e)/2; | |
128 | for (; b<e; m=(b+e)/2) { | |
129 | if (z > fClusters[m]->GetZ()) b=m+1; | |
130 | else e=m; | |
131 | } | |
132 | return m; | |
133 | } | |
134 | ||
3f83f224 | 135 | static int cmp(const void *p1, const void *p2) { |
136 | AliESDtrack *t1=*((AliESDtrack**)p1); | |
137 | AliESDtrack *t2=*((AliESDtrack**)p2); | |
138 | Double_t c1[15]; t1->GetExternalCovariance(c1); | |
139 | Double_t c2[15]; t2->GetExternalCovariance(c2); | |
140 | if (c1[0]*c1[2] <c2[0]*c2[2]) return -1; | |
141 | if (c1[0]*c1[2]==c2[0]*c2[2]) return 0; | |
142 | return 1; | |
143 | } | |
144 | ||
c630aafd | 145 | //_________________________________________________________________________ |
146 | Int_t AliTOFpidESD::MakePID(AliESD *event) | |
147 | { | |
148 | // | |
149 | // This function calculates the "detector response" PID probabilities | |
150 | // Just for a bare hint... | |
151 | ||
3f83f224 | 152 | static const Double_t kMasses[]={ |
c630aafd | 153 | 0.000511, 0.105658, 0.139570, 0.493677, 0.938272, 1.875613 |
154 | }; | |
155 | ||
156 | Int_t ntrk=event->GetNumberOfTracks(); | |
3f83f224 | 157 | AliESDtrack **tracks=new AliESDtrack*[ntrk]; |
158 | ||
159 | Int_t i; | |
160 | for (i=0; i<ntrk; i++) { | |
c630aafd | 161 | AliESDtrack *t=event->GetTrack(i); |
3f83f224 | 162 | tracks[i]=t; |
163 | } | |
164 | qsort(tracks,ntrk,sizeof(AliESDtrack*),cmp); | |
165 | ||
166 | Int_t nmatch=0; | |
167 | for (i=0; i<ntrk; i++) { | |
168 | AliESDtrack *t=tracks[i]; | |
c630aafd | 169 | |
170 | if ((t->GetStatus()&AliESDtrack::kTRDout)==0) continue; | |
b861b2b9 | 171 | if ((t->GetStatus()&AliESDtrack::kTRDStop)!=0) continue; |
c630aafd | 172 | |
c5957288 | 173 | Double_t time[10]; t->GetIntegratedTimes(time); |
174 | ||
c630aafd | 175 | Double_t x,par[5]; t->GetExternalParameters(x,par); |
176 | Double_t cov[15]; t->GetExternalCovariance(cov); | |
177 | ||
3f83f224 | 178 | Double_t dphi=(5*TMath::Sqrt(cov[0]) + 0.5*fDy + 2.5*TMath::Abs(par[2]))/fR; |
179 | Double_t dz=5*TMath::Sqrt(cov[2]) + 0.5*fDz + 2.5*TMath::Abs(par[3]); | |
c630aafd | 180 | |
181 | Double_t phi=TMath::ATan2(par[0],x) + t->GetAlpha(); | |
182 | if (phi<-TMath::Pi()) phi+=2*TMath::Pi(); | |
183 | if (phi>=TMath::Pi()) phi-=2*TMath::Pi(); | |
184 | Double_t z=par[1]; | |
185 | ||
186 | Double_t d2max=1000.; | |
187 | Int_t index=-1; | |
188 | for (Int_t k=FindClusterIndex(z-dz); k<fN; k++) { | |
189 | AliTOFcluster *c=fClusters[k]; | |
190 | ||
191 | if (c->GetZ() > z+dz) break; | |
192 | if (c->IsUsed()) continue; | |
193 | ||
c5957288 | 194 | Double_t tof=50*c->GetTDC()+32; |
195 | if (t->GetIntegratedLength()/tof > 0.031) continue; | |
196 | if (tof>35000) continue; | |
197 | ||
3f83f224 | 198 | Double_t dph=TMath::Abs(c->GetPhi()-phi); |
c5957288 | 199 | if (dph>TMath::Pi()) dph=2*TMath::Pi()-dph; //Thanks to B.Zagreev |
3f83f224 | 200 | if (dph>dphi) continue; |
c630aafd | 201 | |
3f83f224 | 202 | Double_t d2=dph*dph*fR*fR + (c->GetZ()-z)*(c->GetZ()-z); |
c630aafd | 203 | if (d2 > d2max) continue; |
204 | ||
205 | d2max=d2; | |
206 | index=k; | |
207 | } | |
208 | ||
209 | if (index<0) { | |
210 | //Info("MakePID","matching failed ! %d",TMath::Abs(t->GetLabel())); | |
211 | continue; | |
212 | } | |
213 | ||
214 | nmatch++; | |
215 | ||
216 | AliTOFcluster *c=fClusters[index]; | |
3f83f224 | 217 | c->Use(); |
218 | ||
219 | Double_t tof=50*c->GetTDC()+32; // in ps | |
220 | t->SetTOFsignal(tof); | |
221 | t->SetTOFcluster(c->GetIndex()); | |
222 | ||
223 | if ((t->GetStatus()&AliESDtrack::kTIME)==0) continue; | |
c630aafd | 224 | |
3f83f224 | 225 | //track length correction |
226 | Double_t rc=TMath::Sqrt(c->GetR()*c->GetR() + c->GetZ()*c->GetZ()); | |
227 | Double_t rt=TMath::Sqrt(x*x + par[0]*par[0] + par[1]*par[1]); | |
228 | Double_t dlt=rc-rt; | |
c630aafd | 229 | |
230 | Double_t p[10]; | |
231 | Double_t mom=t->GetP(); | |
232 | for (Int_t j=0; j<AliESDtrack::kSPECIES; j++) { | |
3f83f224 | 233 | Double_t mass=kMasses[j]; |
234 | Double_t dpp=0.01; //mean relative pt resolution; | |
235 | if (mom>0.5) dpp=0.01*mom; | |
236 | Double_t sigma=dpp*time[j]/(1.+ mom*mom/(mass*mass)); | |
c630aafd | 237 | sigma=TMath::Sqrt(sigma*sigma + fSigma*fSigma); |
3f83f224 | 238 | |
239 | time[j]+=dlt/3e-2*TMath::Sqrt(mom*mom+mass*mass)/mom; | |
240 | ||
c630aafd | 241 | if (TMath::Abs(tof-time[j]) > fRange*sigma) { |
431be10d | 242 | p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma; |
c630aafd | 243 | continue; |
244 | } | |
431be10d | 245 | p[j]=TMath::Exp(-0.5*(tof-time[j])*(tof-time[j])/(sigma*sigma))/sigma; |
c630aafd | 246 | } |
3f83f224 | 247 | /* |
248 | if (c->GetLabel(0)!=TMath::Abs(t->GetLabel())) { | |
249 | cerr<<"Wrong matching: "<<t->GetLabel()<<endl; | |
250 | continue; | |
251 | } | |
252 | */ | |
c630aafd | 253 | t->SetTOFpid(p); |
254 | ||
c630aafd | 255 | } |
256 | ||
257 | Info("MakePID","Number of matched ESD track: %d",nmatch); | |
258 | ||
3f83f224 | 259 | delete[] tracks; |
260 | ||
c630aafd | 261 | return 0; |
262 | } | |
0f4a7374 | 263 |