]>
Commit | Line | Data |
---|---|---|
4fb1e9d1 | 1 | //------------------------------------------------------------------------- |
2 | // Implementation of the ITS tracker class | |
3 | // The pattern recongintion based on the "cooked covariance" approach | |
4 | //------------------------------------------------------------------------- | |
5 | ||
4fb1e9d1 | 6 | #include <TTree.h> |
7 | #include <TClonesArray.h> | |
8 | ||
9 | #include "AliLog.h" | |
10 | #include "AliESDEvent.h" | |
11 | #include "AliITSUClusterPix.h" | |
12 | #include "AliITSUGeomTGeo.h" | |
13 | #include "AliITSUTrackerCooked.h" | |
14 | #include "AliITSUTrackCooked.h" | |
c64acbbe | 15 | #include "AliITSUReconstructor.h" |
4fb1e9d1 | 16 | |
17 | ClassImp(AliITSUTrackerCooked) | |
18 | ||
19 | //************************************************ | |
20 | // Constants hardcoded for the moment: | |
21 | //************************************************ | |
4fb1e9d1 | 22 | // seed "windows" in z and phi: MakeSeeds |
68ef59d1 | 23 | const Double_t kzWin=0.33; |
24 | const Double_t kminPt=0.05; | |
4fb1e9d1 | 25 | // Maximal accepted impact parameters for the seeds |
8cd74e8c | 26 | const Double_t kmaxDCAxy=3.; |
27 | const Double_t kmaxDCAz= 3.; | |
ae14edfe | 28 | // Layers for the seeding |
29 | const Int_t kSeedingLayer1=6, kSeedingLayer2=4, kSeedingLayer3=5; | |
66be9a4e | 30 | // Space point resolution |
31 | const Double_t kSigma2=0.0005*0.0005; | |
acfdd368 | 32 | // Max accepted chi2 |
33 | const Double_t kmaxChi2PerCluster=20.; | |
34 | const Double_t kmaxChi2PerTrack=30.; | |
a00907de | 35 | // Tracking "road" from layer to layer |
36 | const Double_t kRoadY=0.7; | |
37 | const Double_t kRoadZ=0.7; | |
ae63ad3b | 38 | // Minimal number of attached clusters |
c64acbbe | 39 | const Int_t kminNumberOfClusters=4; |
4fb1e9d1 | 40 | |
41 | //************************************************ | |
42 | // TODO: | |
43 | //************************************************ | |
44 | // Seeding: | |
45 | // Precalculate cylidnrical (r,phi) for the clusters; | |
46 | // use exact r's for the clusters | |
47 | ||
48 | ||
49 | AliITSUTrackerCooked::AliITSUlayer | |
50 | AliITSUTrackerCooked::fgLayers[AliITSUTrackerCooked::kNLayers]; | |
51 | ||
c64acbbe | 52 | AliITSUTrackerCooked::AliITSUTrackerCooked(AliITSUReconstructor *rec): |
53 | AliITSUTrackerGlo(rec), | |
66be9a4e | 54 | fSeeds(0), |
55 | fI(kNLayers-1), | |
56 | fBestTrack(0), | |
c64acbbe | 57 | fTrackToFollow(0), |
58 | fSAonly(kTRUE) | |
4fb1e9d1 | 59 | { |
60 | //-------------------------------------------------------------------- | |
61 | // This default constructor needs to be provided | |
62 | //-------------------------------------------------------------------- | |
45eedad7 | 63 | const Double_t |
64 | klRadius[7]={2.34, 3.15, 3.93, 19.61, 24.55, 34.39, 39.34}; //tdr6 | |
4fb1e9d1 | 65 | |
66 | AliITSUGeomTGeo *gm = new AliITSUGeomTGeo(kTRUE,kTRUE); | |
67 | AliITSUClusterPix::SetGeom(gm); | |
68 | ||
69 | for (Int_t i=0; i<kNLayers; i++) fgLayers[i].SetR(klRadius[i]); | |
70 | ||
71 | // Some default primary vertex | |
72 | Double_t xyz[]={0.,0.,0.}; | |
73 | Double_t ers[]={2.,2.,2.}; | |
74 | ||
75 | SetVertex(xyz,ers); | |
76 | ||
77 | } | |
78 | ||
66be9a4e | 79 | void AliITSUTrackerCooked::ResetTrackToFollow(const AliITSUTrackCooked &t) { |
80 | //-------------------------------------------------------------------- | |
81 | // Prepare to follow a new track seed | |
82 | //-------------------------------------------------------------------- | |
83 | delete fTrackToFollow; | |
84 | fTrackToFollow = new AliITSUTrackCooked(t); | |
85 | } | |
86 | ||
87 | void AliITSUTrackerCooked::ResetBestTrack() { | |
88 | //-------------------------------------------------------------------- | |
89 | // Replace the best track branch | |
90 | //-------------------------------------------------------------------- | |
91 | delete fBestTrack; | |
92 | fBestTrack = new AliITSUTrackCooked(*fTrackToFollow); | |
93 | } | |
94 | ||
4fb1e9d1 | 95 | AliITSUTrackerCooked::~AliITSUTrackerCooked() |
96 | { | |
97 | //-------------------------------------------------------------------- | |
98 | // Virtual destructor | |
99 | //-------------------------------------------------------------------- | |
100 | ||
66be9a4e | 101 | if (fSeeds) fSeeds->Delete(); delete fSeeds; |
102 | delete fBestTrack; | |
103 | delete fTrackToFollow; | |
4fb1e9d1 | 104 | |
4fb1e9d1 | 105 | } |
106 | ||
107 | static Double_t | |
108 | f1(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3) | |
109 | { | |
110 | //----------------------------------------------------------------- | |
111 | // Initial approximation of the track curvature | |
112 | //----------------------------------------------------------------- | |
113 | Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1); | |
114 | Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)- | |
115 | (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2)); | |
116 | Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)- | |
117 | (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1)); | |
118 | ||
68ef59d1 | 119 | Double_t xr=TMath::Abs(d/(d*x1-a)), yr=TMath::Abs(d/(d*y1-b)); |
4fb1e9d1 | 120 | |
68ef59d1 | 121 | Double_t crv=xr*yr/sqrt(xr*xr+yr*yr); |
122 | if (d>0) crv=-crv; | |
123 | ||
124 | return crv; | |
4fb1e9d1 | 125 | } |
126 | ||
127 | static Double_t | |
128 | f2(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3) | |
129 | { | |
130 | //----------------------------------------------------------------- | |
68ef59d1 | 131 | // Initial approximation of the x-coordinate of the center of curvature |
4fb1e9d1 | 132 | //----------------------------------------------------------------- |
68ef59d1 | 133 | |
134 | Double_t k1=(y2-y1)/(x2-x1), k2=(y3-y2)/(x3-x2); | |
135 | Double_t x0=0.5*(k1*k2*(y1-y3) + k2*(x1+x2) - k1*(x2+x3))/(k2-k1); | |
136 | ||
137 | return x0; | |
4fb1e9d1 | 138 | } |
139 | ||
140 | static Double_t | |
141 | f3(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t z1, Double_t z2) | |
142 | { | |
143 | //----------------------------------------------------------------- | |
144 | // Initial approximation of the tangent of the track dip angle | |
145 | //----------------------------------------------------------------- | |
146 | return (z1 - z2)/sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2)); | |
147 | } | |
148 | ||
149 | Bool_t AliITSUTrackerCooked:: | |
150 | AddCookedSeed(const Float_t r1[3], Int_t l1, Int_t i1, | |
151 | const Float_t r2[3], Int_t l2, Int_t i2, | |
152 | const AliCluster *c3,Int_t l3, Int_t i3) | |
153 | { | |
154 | //-------------------------------------------------------------------- | |
155 | // This is the main cooking function. | |
156 | // Creates seed parameters out of provided clusters. | |
157 | //-------------------------------------------------------------------- | |
158 | Float_t x,a; | |
159 | if (!c3->GetXAlphaRefPlane(x,a)) return kFALSE; | |
160 | ||
8cd74e8c | 161 | Double_t ca=TMath::Cos(a), sa=TMath::Sin(a); |
162 | Double_t x1 = r1[0]*ca + r1[1]*sa, | |
163 | y1 =-r1[0]*sa + r1[1]*ca, z1 = r1[2]; | |
164 | Double_t x2 = r2[0]*ca + r2[1]*sa, | |
165 | y2 =-r2[0]*sa + r2[1]*ca, z2 = r2[2]; | |
166 | Double_t x3 = x, y3 = c3->GetY(), z3 = c3->GetZ(); | |
167 | ||
4fb1e9d1 | 168 | Double_t par[5]; |
8cd74e8c | 169 | par[0]=y3; |
170 | par[1]=z3; | |
171 | Double_t crv=f1(x1, y1, x2, y2, x3, y3); //curvature | |
68ef59d1 | 172 | Double_t x0 =f2(x1, y1, x2, y2, x3, y3); //x-coordinate of the center |
8cd74e8c | 173 | Double_t tgl12=f3(x1, y1, x2, y2, z1, z2); |
174 | Double_t tgl23=f3(x2, y2, x3, y3, z2, z3); | |
175 | ||
68ef59d1 | 176 | Double_t sf=crv*(x-x0); |
8cd74e8c | 177 | if (TMath::Abs(sf) >= kAlmost1) return kFALSE; |
178 | par[2]=sf; | |
179 | ||
4fb1e9d1 | 180 | par[3]=0.5*(tgl12 + tgl23); |
181 | Double_t bz=GetBz(); | |
ae14edfe | 182 | par[4]=(TMath::Abs(bz) < kAlmost0Field) ? kAlmost0 : crv/(bz*kB2C); |
4fb1e9d1 | 183 | |
184 | Double_t cov[15]; | |
acfdd368 | 185 | /* |
4fb1e9d1 | 186 | for (Int_t i=0; i<15; i++) cov[i]=0.; |
66be9a4e | 187 | cov[0] =kSigma2*10; |
188 | cov[2] =kSigma2*10; | |
4fb1e9d1 | 189 | cov[5] =0.007*0.007*10; //FIXME all these lines |
190 | cov[9] =0.007*0.007*10; | |
191 | cov[14]=0.1*0.1*10; | |
acfdd368 | 192 | */ |
193 | const Double_t dlt=0.0005; | |
194 | Double_t | |
195 | fy=1./(fgLayers[kSeedingLayer3].GetR() - fgLayers[kSeedingLayer2].GetR()); | |
196 | Double_t tz=fy; | |
197 | Double_t cy=(f1(x1, y1, x2, y2+dlt, x3, y3) - crv)/dlt/bz/kB2C; | |
198 | cy*=20; //FIXME: MS contribution to the cov[14] | |
199 | Double_t s2=kSigma2; | |
200 | ||
201 | cov[0]=s2; | |
202 | cov[1]=0.; cov[2]=s2; | |
203 | cov[3]=s2*fy; cov[4]=0.; cov[5]=s2*fy*fy; | |
204 | cov[6]=0.; cov[7]=s2*tz; cov[8]=0.; cov[9]=s2*tz*tz; | |
205 | cov[10]=s2*cy; cov[11]=0.; cov[12]=s2*fy*cy; cov[13]=0.; cov[14]=s2*cy*cy; | |
4fb1e9d1 | 206 | |
207 | AliITSUTrackCooked *seed=new AliITSUTrackCooked(); | |
208 | seed->Set(Double_t(x), Double_t(a), par, cov); | |
209 | ||
210 | Float_t dz[2]; | |
211 | seed->GetDZ(GetX(),GetY(),GetZ(),GetBz(),dz); | |
212 | if (TMath::Abs(dz[0]) > kmaxDCAxy) {delete seed; return kFALSE;} | |
213 | if (TMath::Abs(dz[1]) > kmaxDCAz ) {delete seed; return kFALSE;} | |
214 | ||
acfdd368 | 215 | Double_t xx0 = 0.008; // Rough layer thickness |
216 | Double_t radl= 9.36; // Radiation length of Si [cm] | |
217 | Double_t rho = 2.33; // Density of Si [g/cm^3] | |
218 | Double_t mass= 0.139;// Pion | |
219 | if (!seed->CorrectForMeanMaterial(xx0, xx0*radl*rho, mass, kTRUE)) { | |
220 | delete seed; return kFALSE; | |
221 | } | |
222 | ||
4fb1e9d1 | 223 | seed->SetClusterIndex(l1,i1); |
224 | seed->SetClusterIndex(l2,i2); | |
225 | seed->SetClusterIndex(l3,i3); | |
226 | ||
227 | fSeeds->AddLast(seed); | |
228 | ||
229 | return kTRUE; | |
230 | } | |
231 | ||
ae14edfe | 232 | Int_t AliITSUTrackerCooked::MakeSeeds() { |
4fb1e9d1 | 233 | //-------------------------------------------------------------------- |
234 | // This is the main pattern recongition function. | |
235 | // Creates seeds out of two clusters and another point. | |
236 | //-------------------------------------------------------------------- | |
4fb1e9d1 | 237 | if (fSeeds) {fSeeds->Delete(); delete fSeeds;} |
238 | fSeeds=new TObjArray(77777); | |
239 | ||
68ef59d1 | 240 | const Double_t zv=GetZ(); |
4fb1e9d1 | 241 | |
ae14edfe | 242 | AliITSUlayer &layer1=fgLayers[kSeedingLayer1]; |
243 | AliITSUlayer &layer2=fgLayers[kSeedingLayer2]; | |
244 | AliITSUlayer &layer3=fgLayers[kSeedingLayer3]; | |
4fb1e9d1 | 245 | Double_t r1=layer1.GetR(); |
246 | Double_t r2=layer2.GetR(); | |
247 | Double_t r3=layer3.GetR(); | |
248 | ||
68ef59d1 | 249 | const Double_t maxC = TMath::Abs(GetBz()*kB2C/kminPt); |
250 | const Double_t kpWin = TMath::ASin(0.5*maxC*r1) - TMath::ASin(0.5*maxC*r2); | |
251 | ||
4fb1e9d1 | 252 | Int_t nClusters1=layer1.GetNumberOfClusters(); |
253 | Int_t nClusters2=layer2.GetNumberOfClusters(); | |
254 | Int_t nClusters3=layer3.GetNumberOfClusters(); | |
255 | for (Int_t n1=0; n1<nClusters1; n1++) { | |
256 | AliCluster *c1=layer1.GetCluster(n1); | |
257 | // | |
258 | //Int_t lab=c1->GetLabel(0); | |
259 | // | |
260 | Double_t z1=c1->GetZ(); | |
261 | Float_t xyz1[3]; c1->GetGlobalXYZ(xyz1); | |
4cce6dec | 262 | Double_t phi1=layer1.GetClusterPhi(n1); |
4fb1e9d1 | 263 | Double_t zr2=zv + r2/r1*(z1-zv); |
264 | Int_t start2=layer2.FindClusterIndex(zr2-kzWin); | |
265 | for (Int_t n2=start2; n2<nClusters2; n2++) { | |
266 | AliCluster *c2=layer2.GetCluster(n2); | |
267 | // | |
268 | //if (c2->GetLabel(0)!=lab) continue; | |
269 | // | |
270 | Double_t z2=c2->GetZ(); | |
4fb1e9d1 | 271 | if (z2 > (zr2+kzWin)) break; //check in Z |
acfdd368 | 272 | |
4fb1e9d1 | 273 | Float_t xyz2[3]; c2->GetGlobalXYZ(xyz2); |
4cce6dec | 274 | Double_t phi2=layer2.GetClusterPhi(n2); |
4fb1e9d1 | 275 | if (TMath::Abs(phi2-phi1) > kpWin) continue; //check in Phi |
acfdd368 | 276 | |
4fb1e9d1 | 277 | Double_t zr3=z1 + (r3-r1)/(r2-r1)*(z2-z1); |
68ef59d1 | 278 | Double_t crv=f1(xyz1[0], xyz1[1], xyz2[0], xyz2[1], GetX(), GetY()); |
279 | Double_t phir3 = phi1 + 0.5*crv*(r3 - r1); | |
280 | ||
4fb1e9d1 | 281 | Int_t start3=layer3.FindClusterIndex(zr3-kzWin/2); |
282 | for (Int_t n3=start3; n3<nClusters3; n3++) { | |
283 | AliCluster *c3=layer3.GetCluster(n3); | |
284 | // | |
285 | //if (c3->GetLabel(0)!=lab) continue; | |
286 | // | |
287 | Double_t z3=c3->GetZ(); | |
4fb1e9d1 | 288 | if (z3 > (zr3+kzWin/2)) break; //check in Z |
68ef59d1 | 289 | |
4fb1e9d1 | 290 | Float_t xyz3[3]; c3->GetGlobalXYZ(xyz3); |
4cce6dec | 291 | Double_t phi3=layer3.GetClusterPhi(n3); |
58c05301 | 292 | if (TMath::Abs(phir3-phi3) > kpWin/100) continue; //check in Phi |
4fb1e9d1 | 293 | |
294 | AliITSUClusterPix cc(*((AliITSUClusterPix*)c2)); | |
295 | cc.GoToFrameTrk(); | |
ae14edfe | 296 | AddCookedSeed(xyz1, kSeedingLayer1, n1, |
297 | xyz3, kSeedingLayer3, n3, | |
298 | &cc, kSeedingLayer2, n2); | |
4fb1e9d1 | 299 | |
300 | } | |
301 | } | |
302 | } | |
ae63ad3b | 303 | |
304 | for (Int_t n1=0; n1<nClusters1; n1++) { | |
305 | AliCluster *c1=layer1.GetCluster(n1); | |
306 | ((AliITSUClusterPix*)c1)->GoToFrameTrk(); | |
307 | } | |
308 | for (Int_t n2=0; n2<nClusters2; n2++) { | |
309 | AliCluster *c2=layer2.GetCluster(n2); | |
310 | ((AliITSUClusterPix*)c2)->GoToFrameTrk(); | |
311 | } | |
312 | for (Int_t n3=0; n3<nClusters3; n3++) { | |
313 | AliCluster *c3=layer3.GetCluster(n3); | |
314 | ((AliITSUClusterPix*)c3)->GoToFrameTrk(); | |
315 | } | |
316 | ||
d5a6c710 | 317 | fSeeds->Sort(); |
4fb1e9d1 | 318 | return fSeeds->GetEntriesFast(); |
319 | } | |
320 | ||
321 | Int_t AliITSUTrackerCooked::Clusters2Tracks(AliESDEvent *event) { | |
322 | //-------------------------------------------------------------------- | |
323 | // This is the main tracking function | |
324 | // The clusters must already be loaded | |
325 | //-------------------------------------------------------------------- | |
4fb1e9d1 | 326 | |
c64acbbe | 327 | if (!fSAonly) AliITSUTrackerGlo::Clusters2Tracks(event); |
328 | ||
ae14edfe | 329 | Int_t nSeeds=MakeSeeds(); |
4fb1e9d1 | 330 | |
331 | // Possibly, icrement the seeds with additional clusters (Kalman) | |
332 | ||
333 | // Possibly, (re)fit the found tracks | |
334 | ||
ae63ad3b | 335 | Int_t ngood=0; |
336 | for (Int_t s=0; s<nSeeds; s++) { | |
337 | const AliITSUTrackCooked *track=(AliITSUTrackCooked*)fSeeds->At(s); | |
4cce6dec | 338 | |
339 | Double_t x=track->GetX(); | |
340 | Double_t y=track->GetY(); | |
341 | Double_t phi=track->GetAlpha() + TMath::ATan2(y,x); | |
342 | const Float_t pi2 = 2.*TMath::Pi(); | |
343 | if (phi<0.) phi+=pi2; | |
344 | else if (phi >= pi2) phi-=pi2; | |
345 | for (Int_t n=0; n<kNLayers-3; n++) { | |
346 | Double_t z; | |
347 | track->GetZAt(fgLayers[n].GetR(),GetBz(),z); | |
348 | fgLayers[n].SelectClusters(phi,kRoadY,z,kRoadZ); | |
349 | } | |
350 | ||
ae63ad3b | 351 | ResetTrackToFollow(*track); |
352 | ResetBestTrack(); | |
353 | fI=kSeedingLayer2; | |
354 | fgLayers[fI].ResetTrack(*track); | |
355 | ||
356 | for (FollowProlongation(); fI<kSeedingLayer2; fI++) { | |
357 | while (TakeNextProlongation()) FollowProlongation(); | |
358 | } | |
359 | ||
360 | if (fBestTrack->GetNumberOfClusters() < kminNumberOfClusters) continue; | |
361 | ||
362 | CookLabel(fBestTrack,0.); //For comparison only | |
363 | Int_t label=fBestTrack->GetLabel(); | |
364 | if (label>0) ngood++; | |
365 | ||
366 | AliESDtrack iotrack; | |
367 | iotrack.UpdateTrackParams(fBestTrack,AliESDtrack::kITSin); | |
368 | iotrack.SetLabel(label); | |
369 | event->AddTrack(&iotrack); | |
370 | UseClusters(fBestTrack); | |
371 | } | |
66be9a4e | 372 | |
ae63ad3b | 373 | Info("Clusters2Tracks","Seeds: %d",nSeeds); |
374 | if (nSeeds) | |
375 | Info("Clusters2Tracks","Good tracks/seeds: %f",Float_t(ngood)/nSeeds); | |
4fb1e9d1 | 376 | |
377 | if (fSeeds) {fSeeds->Delete(); delete fSeeds;} | |
378 | fSeeds=0; | |
379 | ||
380 | return 0; | |
381 | } | |
382 | ||
66be9a4e | 383 | void AliITSUTrackerCooked::FollowProlongation() { |
384 | //-------------------------------------------------------------------- | |
385 | // Push this track tree branch towards the primary vertex | |
386 | //-------------------------------------------------------------------- | |
387 | while (fI) { | |
388 | fI--; | |
4cce6dec | 389 | fgLayers[fI].ResetSelectedClusters(); |
66be9a4e | 390 | if (!TakeNextProlongation()) return; |
66be9a4e | 391 | } |
392 | ||
393 | //deal with the best track | |
394 | Int_t ncl=fTrackToFollow->GetNumberOfClusters(); | |
395 | Int_t nclb=fBestTrack->GetNumberOfClusters(); | |
66be9a4e | 396 | if (ncl >= nclb) { |
397 | Double_t chi2=fTrackToFollow->GetChi2(); | |
acfdd368 | 398 | if (chi2 < kmaxChi2PerTrack) { |
66be9a4e | 399 | if (ncl > nclb || chi2 < fBestTrack->GetChi2()) { |
400 | ResetBestTrack(); | |
401 | } | |
402 | } | |
403 | } | |
404 | ||
405 | } | |
406 | ||
407 | Int_t AliITSUTrackerCooked::TakeNextProlongation() { | |
408 | //-------------------------------------------------------------------- | |
409 | // Switch to the next track tree branch | |
410 | //-------------------------------------------------------------------- | |
411 | AliITSUlayer &layer=fgLayers[fI]; | |
412 | ||
413 | const AliCluster *c=0; Int_t ci=-1; | |
414 | const AliCluster *cc=0; Int_t cci=-1; | |
415 | UShort_t volId=-1; | |
867f050b | 416 | Double_t z=0., dz=0., y=0., dy=0., chi2=0.; |
66be9a4e | 417 | while ((c=layer.GetNextCluster(ci))!=0) { |
66be9a4e | 418 | Int_t id=c->GetVolumeId(); |
419 | if (id != volId) { | |
420 | volId=id; | |
66be9a4e | 421 | const AliITSUTrackCooked *t = fgLayers[fI+1].GetTrack(); |
422 | ResetTrackToFollow(*t); | |
4cce6dec | 423 | Double_t x=layer.GetXRef(ci); |
424 | Double_t alpha=layer.GetAlphaRef(ci); | |
66be9a4e | 425 | if (!fTrackToFollow->Propagate(alpha, x, GetBz())) { |
426 | //Warning("TakeNextProlongation","propagation failed !\n"); | |
427 | continue; | |
428 | } | |
429 | dz=7*TMath::Sqrt(fTrackToFollow->GetSigmaZ2() + kSigma2); | |
430 | dy=7*TMath::Sqrt(fTrackToFollow->GetSigmaY2() + kSigma2); | |
431 | z=fTrackToFollow->GetZ(); | |
432 | y=fTrackToFollow->GetY(); | |
433 | } | |
434 | ||
435 | //if (TMath::Abs(fTrackToFollow.GetZ()-GetZ())>layer.GetR()+dz) continue; | |
436 | ||
437 | if (TMath::Abs(z - c->GetZ()) > dz) continue; | |
438 | if (TMath::Abs(y - c->GetY()) > dy) continue; | |
439 | ||
440 | Double_t ch2=fTrackToFollow->GetPredictedChi2(c); | |
441 | if (ch2 > kmaxChi2PerCluster) continue; | |
442 | chi2=ch2; | |
443 | cc=c; cci=ci; | |
444 | break; | |
445 | } | |
446 | ||
447 | if (!cc) return 0; | |
448 | ||
449 | if (!fTrackToFollow->Update(cc,chi2,(fI<<28)+cci)) { | |
450 | //Warning("TakeNextProlongation","filtering failed !\n"); | |
451 | return 0; | |
452 | } | |
453 | Double_t xx0 = (fI > 2) ? 0.008 : 0.003; // Rough layer thickness | |
454 | Double_t x0 = 9.36; // Radiation length of Si [cm] | |
455 | Double_t rho = 2.33; // Density of Si [g/cm^3] | |
456 | Double_t mass = fTrackToFollow->GetMass(); | |
457 | fTrackToFollow->CorrectForMeanMaterial(xx0, xx0*x0*rho, mass, kTRUE); | |
458 | layer.ResetTrack(*fTrackToFollow); | |
459 | ||
4cce6dec | 460 | |
66be9a4e | 461 | return 1; |
462 | } | |
463 | ||
4fb1e9d1 | 464 | Int_t AliITSUTrackerCooked::PropagateBack(AliESDEvent *event) { |
465 | //-------------------------------------------------------------------- | |
466 | // Here, we implement the Kalman smoother ? | |
467 | // The clusters must already be loaded | |
468 | //-------------------------------------------------------------------- | |
ae63ad3b | 469 | Int_t n=event->GetNumberOfTracks(); |
470 | Int_t ntrk=0; | |
471 | Int_t ngood=0; | |
472 | for (Int_t i=0; i<n; i++) { | |
473 | AliESDtrack *esdTrack=event->GetTrack(i); | |
474 | ||
c64acbbe | 475 | if (!esdTrack->IsOn(AliESDtrack::kITSin)) continue; |
476 | if ( esdTrack->IsOn(AliESDtrack::kTPCin)) continue;//skip a TPC+ITS track | |
ae63ad3b | 477 | |
478 | AliITSUTrackCooked track(*esdTrack); | |
479 | ||
480 | ResetTrackToFollow(track); | |
481 | ||
482 | fTrackToFollow->ResetCovariance(10.); fTrackToFollow->ResetClusters(); | |
483 | if (RefitAt(40., fTrackToFollow, &track)) { | |
484 | ||
485 | CookLabel(fTrackToFollow, 0.); //For comparison only | |
486 | Int_t label=fTrackToFollow->GetLabel(); | |
487 | if (label>0) ngood++; | |
488 | ||
489 | esdTrack->UpdateTrackParams(fTrackToFollow,AliESDtrack::kITSout); | |
490 | //UseClusters(fTrackToFollow); | |
491 | ntrk++; | |
492 | } | |
493 | } | |
494 | ||
495 | Info("PropagateBack","Back propagated tracks: %d",ntrk); | |
496 | if (ntrk) | |
497 | Info("PropagateBack","Good tracks/back propagated: %f",Float_t(ngood)/ntrk); | |
c64acbbe | 498 | |
499 | if (!fSAonly) AliITSUTrackerGlo::PropagateBack(event); | |
500 | ||
4fb1e9d1 | 501 | return 0; |
502 | } | |
503 | ||
ae63ad3b | 504 | Bool_t AliITSUTrackerCooked:: |
505 | RefitAt(Double_t xx, AliITSUTrackCooked *t, const AliITSUTrackCooked *c) { | |
506 | //-------------------------------------------------------------------- | |
507 | // This function refits the track "t" at the position "x" using | |
508 | // the clusters from "c" | |
509 | //-------------------------------------------------------------------- | |
510 | Int_t index[kNLayers]; | |
511 | Int_t k; | |
512 | for (k=0; k<kNLayers; k++) index[k]=-1; | |
513 | Int_t nc=c->GetNumberOfClusters(); | |
514 | for (k=0; k<nc; k++) { | |
515 | Int_t idx=c->GetClusterIndex(k), nl=(idx&0xf0000000)>>28; | |
516 | index[nl]=idx; | |
517 | } | |
518 | ||
519 | Int_t from, to, step; | |
520 | if (xx > t->GetX()) { | |
521 | from=0; to=kNLayers; | |
522 | step=+1; | |
523 | } else { | |
524 | from=kNLayers-1; to=-1; | |
525 | step=-1; | |
526 | } | |
527 | ||
528 | for (Int_t i=from; i != to; i += step) { | |
529 | Int_t idx=index[i]; | |
530 | if (idx>=0) { | |
531 | const AliCluster *cl=GetCluster(idx); | |
532 | Float_t xr,ar; cl->GetXAlphaRefPlane(xr, ar); | |
533 | if (!t->Propagate(Double_t(ar), Double_t(xr), GetBz())) { | |
534 | //Warning("RefitAt","propagation failed !\n"); | |
535 | return kFALSE; | |
536 | } | |
537 | Double_t chi2=t->GetPredictedChi2(cl); | |
538 | if (chi2 < kmaxChi2PerCluster) t->Update(cl, chi2, idx); | |
539 | } else { | |
540 | Double_t r=fgLayers[i].GetR(); | |
541 | Double_t phi,z; | |
542 | if (!t->GetPhiZat(r,phi,z)) { | |
543 | //Warning("RefitAt","failed to estimate track !\n"); | |
544 | return kFALSE; | |
545 | } | |
fc239e33 | 546 | if (!t->Propagate(phi, r, GetBz())) { |
6890ade8 | 547 | //Warning("RefitAt","propagation failed !\n"); |
548 | return kFALSE; | |
549 | } | |
ae63ad3b | 550 | } |
551 | Double_t xx0 = (i > 2) ? 0.008 : 0.003; // Rough layer thickness | |
552 | Double_t x0 = 9.36; // Radiation length of Si [cm] | |
553 | Double_t rho = 2.33; // Density of Si [g/cm^3] | |
554 | Double_t mass = t->GetMass(); | |
555 | t->CorrectForMeanMaterial(xx0, -step*xx0*x0*rho, mass, kTRUE); | |
556 | } | |
557 | ||
558 | if (!t->PropagateTo(xx,0.,0.)) return kFALSE; | |
559 | return kTRUE; | |
560 | } | |
561 | ||
4fb1e9d1 | 562 | Int_t AliITSUTrackerCooked::RefitInward(AliESDEvent *event) { |
563 | //-------------------------------------------------------------------- | |
564 | // Some final refit, after the outliers get removed by the smoother ? | |
565 | // The clusters must be loaded | |
566 | //-------------------------------------------------------------------- | |
ae63ad3b | 567 | Int_t n=event->GetNumberOfTracks(); |
568 | Int_t ntrk=0; | |
569 | Int_t ngood=0; | |
570 | for (Int_t i=0; i<n; i++) { | |
571 | AliESDtrack *esdTrack=event->GetTrack(i); | |
572 | ||
c64acbbe | 573 | if (!esdTrack->IsOn(AliESDtrack::kITSout)) continue; |
574 | if ( esdTrack->IsOn(AliESDtrack::kTPCin)) continue;//skip a TPC+ITS track | |
ae63ad3b | 575 | |
576 | AliITSUTrackCooked track(*esdTrack); | |
577 | ResetTrackToFollow(track); | |
578 | ||
579 | fTrackToFollow->ResetCovariance(10.); fTrackToFollow->ResetClusters(); | |
fc239e33 | 580 | if (!RefitAt(2.1, fTrackToFollow, &track)) continue; |
581 | //Cross the beam pipe | |
582 | if (!fTrackToFollow->PropagateTo(1.8, 2.27e-3, 35.28*1.848)) continue; | |
ae63ad3b | 583 | |
fc239e33 | 584 | CookLabel(fTrackToFollow, 0.); //For comparison only |
585 | Int_t label=fTrackToFollow->GetLabel(); | |
586 | if (label>0) ngood++; | |
ae63ad3b | 587 | |
fc239e33 | 588 | esdTrack->UpdateTrackParams(fTrackToFollow,AliESDtrack::kITSrefit); |
589 | //esdTrack->RelateToVertex(event->GetVertex(),GetBz(),33.); | |
590 | //UseClusters(fTrackToFollow); | |
591 | ntrk++; | |
ae63ad3b | 592 | } |
593 | ||
594 | Info("RefitInward","Refitted tracks: %d",ntrk); | |
595 | if (ntrk) | |
596 | Info("RefitInward","Good tracks/refitted: %f",Float_t(ngood)/ntrk); | |
4fb1e9d1 | 597 | |
c64acbbe | 598 | if (!fSAonly) AliITSUTrackerGlo::RefitInward(event); |
599 | ||
4fb1e9d1 | 600 | return 0; |
601 | } | |
602 | ||
603 | Int_t AliITSUTrackerCooked::LoadClusters(TTree *cTree) { | |
604 | //-------------------------------------------------------------------- | |
605 | // This function reads the ITSU clusters from the tree, | |
606 | // sort them, distribute over the internal tracker arrays, etc | |
607 | //-------------------------------------------------------------------- | |
608 | if (!cTree) { | |
609 | AliFatal("No cluster tree !"); | |
610 | return 1; | |
611 | } | |
612 | ||
c64acbbe | 613 | AliITSUTrackerGlo::LoadClusters(cTree); |
4fb1e9d1 | 614 | |
615 | for (Int_t i=0; i<kNLayers; i++) { | |
c64acbbe | 616 | TClonesArray *clusters=fReconstructor->GetClusters(i); |
ae14edfe | 617 | switch (i) { |
618 | case kSeedingLayer1: | |
619 | case kSeedingLayer2: | |
620 | case kSeedingLayer3: | |
c64acbbe | 621 | fgLayers[i].InsertClusters(clusters,kTRUE,fSAonly); |
ae14edfe | 622 | break; |
623 | default: | |
c64acbbe | 624 | fgLayers[i].InsertClusters(clusters,kFALSE,fSAonly); |
ae14edfe | 625 | break; |
626 | } | |
4fb1e9d1 | 627 | } |
628 | ||
629 | return 0; | |
630 | } | |
631 | ||
632 | void AliITSUTrackerCooked::UnloadClusters() { | |
633 | //-------------------------------------------------------------------- | |
634 | // This function unloads ITSU clusters from the RAM | |
635 | //-------------------------------------------------------------------- | |
c64acbbe | 636 | AliITSUTrackerGlo::UnloadClusters(); |
4fb1e9d1 | 637 | for (Int_t i=0; i<kNLayers; i++) fgLayers[i].DeleteClusters(); |
638 | } | |
639 | ||
640 | AliCluster *AliITSUTrackerCooked::GetCluster(Int_t index) const { | |
641 | //-------------------------------------------------------------------- | |
642 | // Return pointer to a given cluster | |
643 | //-------------------------------------------------------------------- | |
644 | Int_t l=(index & 0xf0000000) >> 28; | |
645 | Int_t c=(index & 0x0fffffff) >> 00; | |
646 | return fgLayers[l].GetCluster(c); | |
647 | } | |
66be9a4e | 648 | |
649 | AliITSUTrackerCooked::AliITSUlayer::AliITSUlayer(): | |
650 | fR(0), | |
651 | fN(0), | |
652 | fNsel(0), | |
4cce6dec | 653 | fI(0), |
66be9a4e | 654 | fTrack(0) |
655 | { | |
656 | //-------------------------------------------------------------------- | |
657 | // This default constructor needs to be provided | |
658 | //-------------------------------------------------------------------- | |
659 | for (Int_t i=0; i<kMaxClusterPerLayer; i++) fClusters[i]=0; | |
660 | for (Int_t i=0; i<kMaxSelected; i++) fIndex[i]=-1; | |
661 | } | |
662 | ||
663 | AliITSUTrackerCooked::AliITSUlayer::~AliITSUlayer() | |
664 | { | |
665 | //-------------------------------------------------------------------- | |
666 | // Simple destructor | |
667 | //-------------------------------------------------------------------- | |
66be9a4e | 668 | delete fTrack; |
669 | } | |
670 | ||
671 | void | |
672 | AliITSUTrackerCooked::AliITSUlayer::ResetTrack(const AliITSUTrackCooked &t) { | |
673 | //-------------------------------------------------------------------- | |
674 | // Replace the track estimate at this layer | |
675 | //-------------------------------------------------------------------- | |
676 | delete fTrack; | |
677 | fTrack=new AliITSUTrackCooked(t); | |
678 | } | |
679 | ||
c64acbbe | 680 | void AliITSUTrackerCooked::AliITSUlayer:: |
681 | InsertClusters(TClonesArray *clusters, Bool_t seedingLayer, Bool_t saOnly) | |
66be9a4e | 682 | { |
683 | //-------------------------------------------------------------------- | |
684 | // Load clusters to this layer | |
685 | //-------------------------------------------------------------------- | |
686 | Int_t ncl=clusters->GetEntriesFast(); | |
45eedad7 | 687 | Double_t r=0.; |
688 | for (Int_t i=0; i<ncl; i++) { | |
689 | AliITSUClusterPix *c=(AliITSUClusterPix*)clusters->UncheckedAt(i); | |
c64acbbe | 690 | if (!saOnly) if (c->IsClusterUsed()) continue; |
45eedad7 | 691 | c->GoToFrameGlo(); |
692 | Double_t x=c->GetX(), y=c->GetY(); | |
693 | r += TMath::Sqrt(x*x + y*y); | |
694 | if (!seedingLayer) c->GoToFrameTrk(); | |
66be9a4e | 695 | //if (!c->Misalign()) AliWarning("Can't misalign this cluster !"); |
4cce6dec | 696 | if (InsertCluster(c)) break; |
697 | } | |
698 | if (fN) fR = r/fN; | |
699 | const Float_t pi2 = 2.*TMath::Pi(); | |
700 | for (Int_t i=0; i<fN; i++) { | |
701 | AliCluster *c=fClusters[i]; | |
702 | c->GetXAlphaRefPlane(fXRef[i],fAlphaRef[i]); | |
703 | Float_t xyz[3]; c->GetGlobalXYZ(xyz); | |
704 | Float_t phi=TMath::ATan2(xyz[1],xyz[0]); | |
705 | if (phi<0.) phi+=pi2; | |
706 | else if (phi >= pi2) phi-=pi2; | |
707 | fPhi[i]=phi; | |
66be9a4e | 708 | } |
66be9a4e | 709 | } |
710 | ||
711 | void AliITSUTrackerCooked::AliITSUlayer::DeleteClusters() | |
712 | { | |
713 | //-------------------------------------------------------------------- | |
714 | // Load clusters to this layer | |
715 | //-------------------------------------------------------------------- | |
c64acbbe | 716 | //for (Int_t i=0; i<fN; i++) {delete fClusters[i]; fClusters[i]=0;} |
66be9a4e | 717 | fN=0; |
718 | } | |
719 | ||
720 | Int_t | |
721 | AliITSUTrackerCooked::AliITSUlayer::InsertCluster(AliCluster *c) { | |
722 | //-------------------------------------------------------------------- | |
723 | // This function inserts a cluster to this layer in increasing | |
724 | // order of the cluster's fZ | |
725 | //-------------------------------------------------------------------- | |
726 | if (fN>=kMaxClusterPerLayer) { | |
727 | ::Error("InsertCluster","Too many clusters !\n"); | |
728 | return 1; | |
729 | } | |
730 | if (fN==0) fClusters[0]=c; | |
731 | else { | |
732 | Int_t i=FindClusterIndex(c->GetZ()); | |
733 | Int_t k=fN-i; | |
734 | memmove(fClusters+i+1 ,fClusters+i,k*sizeof(AliCluster*)); | |
735 | fClusters[i]=c; | |
736 | } | |
737 | fN++; | |
738 | return 0; | |
739 | } | |
740 | ||
741 | Int_t | |
742 | AliITSUTrackerCooked::AliITSUlayer::FindClusterIndex(Double_t z) const { | |
743 | //-------------------------------------------------------------------- | |
744 | // This function returns the index of the first | |
745 | // with its fZ >= "z". | |
746 | //-------------------------------------------------------------------- | |
747 | if (fN==0) return 0; | |
748 | ||
749 | Int_t b=0; | |
750 | if (z <= fClusters[b]->GetZ()) return b; | |
751 | ||
752 | Int_t e=b+fN-1; | |
753 | if (z > fClusters[e]->GetZ()) return e+1; | |
754 | ||
755 | Int_t m=(b+e)/2; | |
756 | for (; b<e; m=(b+e)/2) { | |
757 | if (z > fClusters[m]->GetZ()) b=m+1; | |
758 | else e=m; | |
759 | } | |
760 | return m; | |
761 | } | |
762 | ||
4cce6dec | 763 | void AliITSUTrackerCooked::AliITSUlayer:: |
764 | SelectClusters(Float_t phi, Float_t dy, Float_t z, Float_t dz) { | |
66be9a4e | 765 | //-------------------------------------------------------------------- |
766 | // This function selects clusters within the "road" | |
767 | //-------------------------------------------------------------------- | |
4cce6dec | 768 | fNsel=0; |
769 | ||
770 | Float_t dphi=dy/fR; | |
771 | Float_t phiMin=phi-dphi; | |
772 | Float_t phiMax=phi+dphi; | |
773 | Float_t zMin=z-dz; | |
774 | Float_t zMax=z+dz; | |
775 | ||
776 | Int_t i=FindClusterIndex(zMin), imax=FindClusterIndex(zMax); | |
777 | for (; i<imax; i++) { | |
778 | Float_t cphi=fPhi[i]; | |
779 | if (cphi <= phiMin) continue; | |
780 | if (cphi > phiMax) continue; | |
66be9a4e | 781 | AliCluster *c=fClusters[i]; |
a00907de | 782 | if (c->IsClusterUsed()) continue; |
4cce6dec | 783 | |
66be9a4e | 784 | fIndex[fNsel++]=i; |
4cce6dec | 785 | if (fNsel>=kMaxSelected) break; |
786 | } | |
66be9a4e | 787 | } |
788 | ||
789 | const AliCluster *AliITSUTrackerCooked::AliITSUlayer::GetNextCluster(Int_t &ci){ | |
790 | //-------------------------------------------------------------------- | |
791 | // This function returns clusters within the "road" | |
792 | //-------------------------------------------------------------------- | |
4cce6dec | 793 | if (fI<fNsel) { |
794 | ci=fIndex[fI++]; | |
795 | return fClusters[ci]; | |
66be9a4e | 796 | } |
4cce6dec | 797 | ci=-1; |
798 | return 0; | |
66be9a4e | 799 | } |
800 |