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