]>
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 | //-------------------------------------------------------------------- | |
68ef59d1 | 237 | const Double_t zv=GetZ(); |
4fb1e9d1 | 238 | |
ae14edfe | 239 | AliITSUlayer &layer1=fgLayers[kSeedingLayer1]; |
240 | AliITSUlayer &layer2=fgLayers[kSeedingLayer2]; | |
241 | AliITSUlayer &layer3=fgLayers[kSeedingLayer3]; | |
4fb1e9d1 | 242 | Double_t r1=layer1.GetR(); |
243 | Double_t r2=layer2.GetR(); | |
244 | Double_t r3=layer3.GetR(); | |
245 | ||
68ef59d1 | 246 | const Double_t maxC = TMath::Abs(GetBz()*kB2C/kminPt); |
247 | const Double_t kpWin = TMath::ASin(0.5*maxC*r1) - TMath::ASin(0.5*maxC*r2); | |
248 | ||
4fb1e9d1 | 249 | Int_t nClusters1=layer1.GetNumberOfClusters(); |
250 | Int_t nClusters2=layer2.GetNumberOfClusters(); | |
251 | Int_t nClusters3=layer3.GetNumberOfClusters(); | |
252 | for (Int_t n1=0; n1<nClusters1; n1++) { | |
253 | AliCluster *c1=layer1.GetCluster(n1); | |
254 | // | |
255 | //Int_t lab=c1->GetLabel(0); | |
256 | // | |
257 | Double_t z1=c1->GetZ(); | |
258 | Float_t xyz1[3]; c1->GetGlobalXYZ(xyz1); | |
4cce6dec | 259 | Double_t phi1=layer1.GetClusterPhi(n1); |
4fb1e9d1 | 260 | Double_t zr2=zv + r2/r1*(z1-zv); |
261 | Int_t start2=layer2.FindClusterIndex(zr2-kzWin); | |
262 | for (Int_t n2=start2; n2<nClusters2; n2++) { | |
263 | AliCluster *c2=layer2.GetCluster(n2); | |
264 | // | |
265 | //if (c2->GetLabel(0)!=lab) continue; | |
266 | // | |
267 | Double_t z2=c2->GetZ(); | |
4fb1e9d1 | 268 | if (z2 > (zr2+kzWin)) break; //check in Z |
acfdd368 | 269 | |
4fb1e9d1 | 270 | Float_t xyz2[3]; c2->GetGlobalXYZ(xyz2); |
4cce6dec | 271 | Double_t phi2=layer2.GetClusterPhi(n2); |
4fb1e9d1 | 272 | if (TMath::Abs(phi2-phi1) > kpWin) continue; //check in Phi |
acfdd368 | 273 | |
4fb1e9d1 | 274 | Double_t zr3=z1 + (r3-r1)/(r2-r1)*(z2-z1); |
68ef59d1 | 275 | Double_t crv=f1(xyz1[0], xyz1[1], xyz2[0], xyz2[1], GetX(), GetY()); |
276 | Double_t phir3 = phi1 + 0.5*crv*(r3 - r1); | |
277 | ||
4fb1e9d1 | 278 | Int_t start3=layer3.FindClusterIndex(zr3-kzWin/2); |
279 | for (Int_t n3=start3; n3<nClusters3; n3++) { | |
280 | AliCluster *c3=layer3.GetCluster(n3); | |
281 | // | |
282 | //if (c3->GetLabel(0)!=lab) continue; | |
283 | // | |
284 | Double_t z3=c3->GetZ(); | |
4fb1e9d1 | 285 | if (z3 > (zr3+kzWin/2)) break; //check in Z |
68ef59d1 | 286 | |
4fb1e9d1 | 287 | Float_t xyz3[3]; c3->GetGlobalXYZ(xyz3); |
4cce6dec | 288 | Double_t phi3=layer3.GetClusterPhi(n3); |
58c05301 | 289 | if (TMath::Abs(phir3-phi3) > kpWin/100) continue; //check in Phi |
4fb1e9d1 | 290 | |
291 | AliITSUClusterPix cc(*((AliITSUClusterPix*)c2)); | |
292 | cc.GoToFrameTrk(); | |
ae14edfe | 293 | AddCookedSeed(xyz1, kSeedingLayer1, n1, |
294 | xyz3, kSeedingLayer3, n3, | |
295 | &cc, kSeedingLayer2, n2); | |
4fb1e9d1 | 296 | |
297 | } | |
298 | } | |
299 | } | |
ae63ad3b | 300 | |
301 | for (Int_t n1=0; n1<nClusters1; n1++) { | |
302 | AliCluster *c1=layer1.GetCluster(n1); | |
303 | ((AliITSUClusterPix*)c1)->GoToFrameTrk(); | |
304 | } | |
305 | for (Int_t n2=0; n2<nClusters2; n2++) { | |
306 | AliCluster *c2=layer2.GetCluster(n2); | |
307 | ((AliITSUClusterPix*)c2)->GoToFrameTrk(); | |
308 | } | |
309 | for (Int_t n3=0; n3<nClusters3; n3++) { | |
310 | AliCluster *c3=layer3.GetCluster(n3); | |
311 | ((AliITSUClusterPix*)c3)->GoToFrameTrk(); | |
312 | } | |
313 | ||
4fb1e9d1 | 314 | return fSeeds->GetEntriesFast(); |
315 | } | |
316 | ||
317 | Int_t AliITSUTrackerCooked::Clusters2Tracks(AliESDEvent *event) { | |
318 | //-------------------------------------------------------------------- | |
319 | // This is the main tracking function | |
320 | // The clusters must already be loaded | |
321 | //-------------------------------------------------------------------- | |
4fb1e9d1 | 322 | |
c64acbbe | 323 | if (!fSAonly) AliITSUTrackerGlo::Clusters2Tracks(event); |
324 | ||
3bc79d8d | 325 | if (fSeeds) {fSeeds->Delete(); delete fSeeds;} |
326 | fSeeds=new TObjArray(77777); | |
327 | ||
328 | //Seeding with the triggered primary vertex | |
329 | Double_t xyz[3]; | |
330 | const AliESDVertex *vtx=0; | |
331 | vtx=event->GetPrimaryVertexSPD(); | |
332 | if (vtx->GetStatus()) { | |
333 | xyz[0]=vtx->GetXv(); xyz[1]=vtx->GetYv(); xyz[2]=vtx->GetZv(); | |
334 | SetVertex(xyz); | |
335 | MakeSeeds(); | |
336 | } | |
337 | //Seeding with the pileup primary vertices | |
338 | TClonesArray *verticesSPD=event->GetPileupVerticesSPD(); | |
339 | Int_t nfoundSPD=verticesSPD->GetEntries(); | |
340 | for (Int_t v=0; v<nfoundSPD; v++) { | |
341 | vtx=(AliESDVertex *)verticesSPD->UncheckedAt(v); | |
342 | if (!vtx->GetStatus()) continue; | |
343 | xyz[0]=vtx->GetXv(); xyz[1]=vtx->GetYv(); xyz[2]=vtx->GetZv(); | |
344 | SetVertex(xyz); | |
345 | MakeSeeds(); | |
346 | } | |
347 | fSeeds->Sort(); | |
348 | Int_t nSeeds=fSeeds->GetEntriesFast(); | |
4fb1e9d1 | 349 | |
350 | // Possibly, icrement the seeds with additional clusters (Kalman) | |
351 | ||
352 | // Possibly, (re)fit the found tracks | |
353 | ||
ae63ad3b | 354 | Int_t ngood=0; |
355 | for (Int_t s=0; s<nSeeds; s++) { | |
356 | const AliITSUTrackCooked *track=(AliITSUTrackCooked*)fSeeds->At(s); | |
4cce6dec | 357 | |
358 | Double_t x=track->GetX(); | |
359 | Double_t y=track->GetY(); | |
360 | Double_t phi=track->GetAlpha() + TMath::ATan2(y,x); | |
361 | const Float_t pi2 = 2.*TMath::Pi(); | |
362 | if (phi<0.) phi+=pi2; | |
363 | else if (phi >= pi2) phi-=pi2; | |
364 | for (Int_t n=0; n<kNLayers-3; n++) { | |
365 | Double_t z; | |
366 | track->GetZAt(fgLayers[n].GetR(),GetBz(),z); | |
367 | fgLayers[n].SelectClusters(phi,kRoadY,z,kRoadZ); | |
368 | } | |
369 | ||
ae63ad3b | 370 | ResetTrackToFollow(*track); |
371 | ResetBestTrack(); | |
372 | fI=kSeedingLayer2; | |
373 | fgLayers[fI].ResetTrack(*track); | |
374 | ||
375 | for (FollowProlongation(); fI<kSeedingLayer2; fI++) { | |
376 | while (TakeNextProlongation()) FollowProlongation(); | |
377 | } | |
378 | ||
379 | if (fBestTrack->GetNumberOfClusters() < kminNumberOfClusters) continue; | |
380 | ||
381 | CookLabel(fBestTrack,0.); //For comparison only | |
382 | Int_t label=fBestTrack->GetLabel(); | |
383 | if (label>0) ngood++; | |
384 | ||
385 | AliESDtrack iotrack; | |
386 | iotrack.UpdateTrackParams(fBestTrack,AliESDtrack::kITSin); | |
387 | iotrack.SetLabel(label); | |
388 | event->AddTrack(&iotrack); | |
389 | UseClusters(fBestTrack); | |
390 | } | |
66be9a4e | 391 | |
ae63ad3b | 392 | Info("Clusters2Tracks","Seeds: %d",nSeeds); |
393 | if (nSeeds) | |
394 | Info("Clusters2Tracks","Good tracks/seeds: %f",Float_t(ngood)/nSeeds); | |
4fb1e9d1 | 395 | |
396 | if (fSeeds) {fSeeds->Delete(); delete fSeeds;} | |
397 | fSeeds=0; | |
398 | ||
399 | return 0; | |
400 | } | |
401 | ||
66be9a4e | 402 | void AliITSUTrackerCooked::FollowProlongation() { |
403 | //-------------------------------------------------------------------- | |
404 | // Push this track tree branch towards the primary vertex | |
405 | //-------------------------------------------------------------------- | |
406 | while (fI) { | |
407 | fI--; | |
4cce6dec | 408 | fgLayers[fI].ResetSelectedClusters(); |
66be9a4e | 409 | if (!TakeNextProlongation()) return; |
66be9a4e | 410 | } |
411 | ||
412 | //deal with the best track | |
413 | Int_t ncl=fTrackToFollow->GetNumberOfClusters(); | |
414 | Int_t nclb=fBestTrack->GetNumberOfClusters(); | |
66be9a4e | 415 | if (ncl >= nclb) { |
416 | Double_t chi2=fTrackToFollow->GetChi2(); | |
acfdd368 | 417 | if (chi2 < kmaxChi2PerTrack) { |
66be9a4e | 418 | if (ncl > nclb || chi2 < fBestTrack->GetChi2()) { |
419 | ResetBestTrack(); | |
420 | } | |
421 | } | |
422 | } | |
423 | ||
424 | } | |
425 | ||
426 | Int_t AliITSUTrackerCooked::TakeNextProlongation() { | |
427 | //-------------------------------------------------------------------- | |
428 | // Switch to the next track tree branch | |
429 | //-------------------------------------------------------------------- | |
430 | AliITSUlayer &layer=fgLayers[fI]; | |
431 | ||
432 | const AliCluster *c=0; Int_t ci=-1; | |
433 | const AliCluster *cc=0; Int_t cci=-1; | |
434 | UShort_t volId=-1; | |
867f050b | 435 | Double_t z=0., dz=0., y=0., dy=0., chi2=0.; |
66be9a4e | 436 | while ((c=layer.GetNextCluster(ci))!=0) { |
66be9a4e | 437 | Int_t id=c->GetVolumeId(); |
438 | if (id != volId) { | |
439 | volId=id; | |
66be9a4e | 440 | const AliITSUTrackCooked *t = fgLayers[fI+1].GetTrack(); |
441 | ResetTrackToFollow(*t); | |
4cce6dec | 442 | Double_t x=layer.GetXRef(ci); |
443 | Double_t alpha=layer.GetAlphaRef(ci); | |
66be9a4e | 444 | if (!fTrackToFollow->Propagate(alpha, x, GetBz())) { |
445 | //Warning("TakeNextProlongation","propagation failed !\n"); | |
446 | continue; | |
447 | } | |
448 | dz=7*TMath::Sqrt(fTrackToFollow->GetSigmaZ2() + kSigma2); | |
449 | dy=7*TMath::Sqrt(fTrackToFollow->GetSigmaY2() + kSigma2); | |
450 | z=fTrackToFollow->GetZ(); | |
451 | y=fTrackToFollow->GetY(); | |
452 | } | |
453 | ||
454 | //if (TMath::Abs(fTrackToFollow.GetZ()-GetZ())>layer.GetR()+dz) continue; | |
455 | ||
456 | if (TMath::Abs(z - c->GetZ()) > dz) continue; | |
457 | if (TMath::Abs(y - c->GetY()) > dy) continue; | |
458 | ||
459 | Double_t ch2=fTrackToFollow->GetPredictedChi2(c); | |
460 | if (ch2 > kmaxChi2PerCluster) continue; | |
461 | chi2=ch2; | |
462 | cc=c; cci=ci; | |
463 | break; | |
464 | } | |
465 | ||
466 | if (!cc) return 0; | |
467 | ||
468 | if (!fTrackToFollow->Update(cc,chi2,(fI<<28)+cci)) { | |
469 | //Warning("TakeNextProlongation","filtering failed !\n"); | |
470 | return 0; | |
471 | } | |
472 | Double_t xx0 = (fI > 2) ? 0.008 : 0.003; // Rough layer thickness | |
473 | Double_t x0 = 9.36; // Radiation length of Si [cm] | |
474 | Double_t rho = 2.33; // Density of Si [g/cm^3] | |
475 | Double_t mass = fTrackToFollow->GetMass(); | |
476 | fTrackToFollow->CorrectForMeanMaterial(xx0, xx0*x0*rho, mass, kTRUE); | |
477 | layer.ResetTrack(*fTrackToFollow); | |
478 | ||
4cce6dec | 479 | |
66be9a4e | 480 | return 1; |
481 | } | |
482 | ||
4fb1e9d1 | 483 | Int_t AliITSUTrackerCooked::PropagateBack(AliESDEvent *event) { |
484 | //-------------------------------------------------------------------- | |
485 | // Here, we implement the Kalman smoother ? | |
486 | // The clusters must already be loaded | |
487 | //-------------------------------------------------------------------- | |
ae63ad3b | 488 | Int_t n=event->GetNumberOfTracks(); |
489 | Int_t ntrk=0; | |
490 | Int_t ngood=0; | |
491 | for (Int_t i=0; i<n; i++) { | |
492 | AliESDtrack *esdTrack=event->GetTrack(i); | |
493 | ||
c64acbbe | 494 | if (!esdTrack->IsOn(AliESDtrack::kITSin)) continue; |
495 | if ( esdTrack->IsOn(AliESDtrack::kTPCin)) continue;//skip a TPC+ITS track | |
ae63ad3b | 496 | |
497 | AliITSUTrackCooked track(*esdTrack); | |
498 | ||
499 | ResetTrackToFollow(track); | |
500 | ||
501 | fTrackToFollow->ResetCovariance(10.); fTrackToFollow->ResetClusters(); | |
502 | if (RefitAt(40., fTrackToFollow, &track)) { | |
503 | ||
504 | CookLabel(fTrackToFollow, 0.); //For comparison only | |
505 | Int_t label=fTrackToFollow->GetLabel(); | |
506 | if (label>0) ngood++; | |
507 | ||
508 | esdTrack->UpdateTrackParams(fTrackToFollow,AliESDtrack::kITSout); | |
509 | //UseClusters(fTrackToFollow); | |
510 | ntrk++; | |
511 | } | |
512 | } | |
513 | ||
514 | Info("PropagateBack","Back propagated tracks: %d",ntrk); | |
515 | if (ntrk) | |
516 | Info("PropagateBack","Good tracks/back propagated: %f",Float_t(ngood)/ntrk); | |
c64acbbe | 517 | |
518 | if (!fSAonly) AliITSUTrackerGlo::PropagateBack(event); | |
519 | ||
4fb1e9d1 | 520 | return 0; |
521 | } | |
522 | ||
ae63ad3b | 523 | Bool_t AliITSUTrackerCooked:: |
524 | RefitAt(Double_t xx, AliITSUTrackCooked *t, const AliITSUTrackCooked *c) { | |
525 | //-------------------------------------------------------------------- | |
526 | // This function refits the track "t" at the position "x" using | |
527 | // the clusters from "c" | |
528 | //-------------------------------------------------------------------- | |
529 | Int_t index[kNLayers]; | |
530 | Int_t k; | |
531 | for (k=0; k<kNLayers; k++) index[k]=-1; | |
532 | Int_t nc=c->GetNumberOfClusters(); | |
533 | for (k=0; k<nc; k++) { | |
534 | Int_t idx=c->GetClusterIndex(k), nl=(idx&0xf0000000)>>28; | |
535 | index[nl]=idx; | |
536 | } | |
537 | ||
538 | Int_t from, to, step; | |
539 | if (xx > t->GetX()) { | |
540 | from=0; to=kNLayers; | |
541 | step=+1; | |
542 | } else { | |
543 | from=kNLayers-1; to=-1; | |
544 | step=-1; | |
545 | } | |
546 | ||
547 | for (Int_t i=from; i != to; i += step) { | |
548 | Int_t idx=index[i]; | |
549 | if (idx>=0) { | |
550 | const AliCluster *cl=GetCluster(idx); | |
551 | Float_t xr,ar; cl->GetXAlphaRefPlane(xr, ar); | |
552 | if (!t->Propagate(Double_t(ar), Double_t(xr), GetBz())) { | |
553 | //Warning("RefitAt","propagation failed !\n"); | |
554 | return kFALSE; | |
555 | } | |
556 | Double_t chi2=t->GetPredictedChi2(cl); | |
557 | if (chi2 < kmaxChi2PerCluster) t->Update(cl, chi2, idx); | |
558 | } else { | |
559 | Double_t r=fgLayers[i].GetR(); | |
560 | Double_t phi,z; | |
561 | if (!t->GetPhiZat(r,phi,z)) { | |
562 | //Warning("RefitAt","failed to estimate track !\n"); | |
563 | return kFALSE; | |
564 | } | |
fc239e33 | 565 | if (!t->Propagate(phi, r, GetBz())) { |
6890ade8 | 566 | //Warning("RefitAt","propagation failed !\n"); |
567 | return kFALSE; | |
568 | } | |
ae63ad3b | 569 | } |
570 | Double_t xx0 = (i > 2) ? 0.008 : 0.003; // Rough layer thickness | |
571 | Double_t x0 = 9.36; // Radiation length of Si [cm] | |
572 | Double_t rho = 2.33; // Density of Si [g/cm^3] | |
573 | Double_t mass = t->GetMass(); | |
574 | t->CorrectForMeanMaterial(xx0, -step*xx0*x0*rho, mass, kTRUE); | |
575 | } | |
576 | ||
577 | if (!t->PropagateTo(xx,0.,0.)) return kFALSE; | |
578 | return kTRUE; | |
579 | } | |
580 | ||
4fb1e9d1 | 581 | Int_t AliITSUTrackerCooked::RefitInward(AliESDEvent *event) { |
582 | //-------------------------------------------------------------------- | |
583 | // Some final refit, after the outliers get removed by the smoother ? | |
584 | // The clusters must be loaded | |
585 | //-------------------------------------------------------------------- | |
ae63ad3b | 586 | Int_t n=event->GetNumberOfTracks(); |
587 | Int_t ntrk=0; | |
588 | Int_t ngood=0; | |
589 | for (Int_t i=0; i<n; i++) { | |
590 | AliESDtrack *esdTrack=event->GetTrack(i); | |
591 | ||
c64acbbe | 592 | if (!esdTrack->IsOn(AliESDtrack::kITSout)) continue; |
593 | if ( esdTrack->IsOn(AliESDtrack::kTPCin)) continue;//skip a TPC+ITS track | |
ae63ad3b | 594 | |
595 | AliITSUTrackCooked track(*esdTrack); | |
596 | ResetTrackToFollow(track); | |
597 | ||
598 | fTrackToFollow->ResetCovariance(10.); fTrackToFollow->ResetClusters(); | |
fc239e33 | 599 | if (!RefitAt(2.1, fTrackToFollow, &track)) continue; |
600 | //Cross the beam pipe | |
601 | if (!fTrackToFollow->PropagateTo(1.8, 2.27e-3, 35.28*1.848)) continue; | |
ae63ad3b | 602 | |
fc239e33 | 603 | CookLabel(fTrackToFollow, 0.); //For comparison only |
604 | Int_t label=fTrackToFollow->GetLabel(); | |
605 | if (label>0) ngood++; | |
ae63ad3b | 606 | |
fc239e33 | 607 | esdTrack->UpdateTrackParams(fTrackToFollow,AliESDtrack::kITSrefit); |
608 | //esdTrack->RelateToVertex(event->GetVertex(),GetBz(),33.); | |
609 | //UseClusters(fTrackToFollow); | |
610 | ntrk++; | |
ae63ad3b | 611 | } |
612 | ||
613 | Info("RefitInward","Refitted tracks: %d",ntrk); | |
614 | if (ntrk) | |
615 | Info("RefitInward","Good tracks/refitted: %f",Float_t(ngood)/ntrk); | |
4fb1e9d1 | 616 | |
c64acbbe | 617 | if (!fSAonly) AliITSUTrackerGlo::RefitInward(event); |
618 | ||
4fb1e9d1 | 619 | return 0; |
620 | } | |
621 | ||
622 | Int_t AliITSUTrackerCooked::LoadClusters(TTree *cTree) { | |
623 | //-------------------------------------------------------------------- | |
624 | // This function reads the ITSU clusters from the tree, | |
625 | // sort them, distribute over the internal tracker arrays, etc | |
626 | //-------------------------------------------------------------------- | |
627 | if (!cTree) { | |
628 | AliFatal("No cluster tree !"); | |
629 | return 1; | |
630 | } | |
631 | ||
c64acbbe | 632 | AliITSUTrackerGlo::LoadClusters(cTree); |
4fb1e9d1 | 633 | |
634 | for (Int_t i=0; i<kNLayers; i++) { | |
c64acbbe | 635 | TClonesArray *clusters=fReconstructor->GetClusters(i); |
ae14edfe | 636 | switch (i) { |
637 | case kSeedingLayer1: | |
638 | case kSeedingLayer2: | |
639 | case kSeedingLayer3: | |
c64acbbe | 640 | fgLayers[i].InsertClusters(clusters,kTRUE,fSAonly); |
ae14edfe | 641 | break; |
642 | default: | |
c64acbbe | 643 | fgLayers[i].InsertClusters(clusters,kFALSE,fSAonly); |
ae14edfe | 644 | break; |
645 | } | |
4fb1e9d1 | 646 | } |
647 | ||
648 | return 0; | |
649 | } | |
650 | ||
651 | void AliITSUTrackerCooked::UnloadClusters() { | |
652 | //-------------------------------------------------------------------- | |
653 | // This function unloads ITSU clusters from the RAM | |
654 | //-------------------------------------------------------------------- | |
c64acbbe | 655 | AliITSUTrackerGlo::UnloadClusters(); |
4fb1e9d1 | 656 | for (Int_t i=0; i<kNLayers; i++) fgLayers[i].DeleteClusters(); |
657 | } | |
658 | ||
659 | AliCluster *AliITSUTrackerCooked::GetCluster(Int_t index) const { | |
660 | //-------------------------------------------------------------------- | |
661 | // Return pointer to a given cluster | |
662 | //-------------------------------------------------------------------- | |
663 | Int_t l=(index & 0xf0000000) >> 28; | |
664 | Int_t c=(index & 0x0fffffff) >> 00; | |
665 | return fgLayers[l].GetCluster(c); | |
666 | } | |
66be9a4e | 667 | |
668 | AliITSUTrackerCooked::AliITSUlayer::AliITSUlayer(): | |
669 | fR(0), | |
670 | fN(0), | |
671 | fNsel(0), | |
4cce6dec | 672 | fI(0), |
66be9a4e | 673 | fTrack(0) |
674 | { | |
675 | //-------------------------------------------------------------------- | |
676 | // This default constructor needs to be provided | |
677 | //-------------------------------------------------------------------- | |
678 | for (Int_t i=0; i<kMaxClusterPerLayer; i++) fClusters[i]=0; | |
679 | for (Int_t i=0; i<kMaxSelected; i++) fIndex[i]=-1; | |
680 | } | |
681 | ||
682 | AliITSUTrackerCooked::AliITSUlayer::~AliITSUlayer() | |
683 | { | |
684 | //-------------------------------------------------------------------- | |
685 | // Simple destructor | |
686 | //-------------------------------------------------------------------- | |
66be9a4e | 687 | delete fTrack; |
688 | } | |
689 | ||
690 | void | |
691 | AliITSUTrackerCooked::AliITSUlayer::ResetTrack(const AliITSUTrackCooked &t) { | |
692 | //-------------------------------------------------------------------- | |
693 | // Replace the track estimate at this layer | |
694 | //-------------------------------------------------------------------- | |
695 | delete fTrack; | |
696 | fTrack=new AliITSUTrackCooked(t); | |
697 | } | |
698 | ||
c64acbbe | 699 | void AliITSUTrackerCooked::AliITSUlayer:: |
700 | InsertClusters(TClonesArray *clusters, Bool_t seedingLayer, Bool_t saOnly) | |
66be9a4e | 701 | { |
702 | //-------------------------------------------------------------------- | |
703 | // Load clusters to this layer | |
704 | //-------------------------------------------------------------------- | |
705 | Int_t ncl=clusters->GetEntriesFast(); | |
45eedad7 | 706 | Double_t r=0.; |
707 | for (Int_t i=0; i<ncl; i++) { | |
708 | AliITSUClusterPix *c=(AliITSUClusterPix*)clusters->UncheckedAt(i); | |
c64acbbe | 709 | if (!saOnly) if (c->IsClusterUsed()) continue; |
45eedad7 | 710 | c->GoToFrameGlo(); |
711 | Double_t x=c->GetX(), y=c->GetY(); | |
712 | r += TMath::Sqrt(x*x + y*y); | |
713 | if (!seedingLayer) c->GoToFrameTrk(); | |
66be9a4e | 714 | //if (!c->Misalign()) AliWarning("Can't misalign this cluster !"); |
4cce6dec | 715 | if (InsertCluster(c)) break; |
716 | } | |
717 | if (fN) fR = r/fN; | |
718 | const Float_t pi2 = 2.*TMath::Pi(); | |
719 | for (Int_t i=0; i<fN; i++) { | |
720 | AliCluster *c=fClusters[i]; | |
721 | c->GetXAlphaRefPlane(fXRef[i],fAlphaRef[i]); | |
722 | Float_t xyz[3]; c->GetGlobalXYZ(xyz); | |
723 | Float_t phi=TMath::ATan2(xyz[1],xyz[0]); | |
724 | if (phi<0.) phi+=pi2; | |
725 | else if (phi >= pi2) phi-=pi2; | |
726 | fPhi[i]=phi; | |
66be9a4e | 727 | } |
66be9a4e | 728 | } |
729 | ||
730 | void AliITSUTrackerCooked::AliITSUlayer::DeleteClusters() | |
731 | { | |
732 | //-------------------------------------------------------------------- | |
733 | // Load clusters to this layer | |
734 | //-------------------------------------------------------------------- | |
c64acbbe | 735 | //for (Int_t i=0; i<fN; i++) {delete fClusters[i]; fClusters[i]=0;} |
66be9a4e | 736 | fN=0; |
737 | } | |
738 | ||
739 | Int_t | |
740 | AliITSUTrackerCooked::AliITSUlayer::InsertCluster(AliCluster *c) { | |
741 | //-------------------------------------------------------------------- | |
742 | // This function inserts a cluster to this layer in increasing | |
743 | // order of the cluster's fZ | |
744 | //-------------------------------------------------------------------- | |
745 | if (fN>=kMaxClusterPerLayer) { | |
746 | ::Error("InsertCluster","Too many clusters !\n"); | |
747 | return 1; | |
748 | } | |
749 | if (fN==0) fClusters[0]=c; | |
750 | else { | |
751 | Int_t i=FindClusterIndex(c->GetZ()); | |
752 | Int_t k=fN-i; | |
753 | memmove(fClusters+i+1 ,fClusters+i,k*sizeof(AliCluster*)); | |
754 | fClusters[i]=c; | |
755 | } | |
756 | fN++; | |
757 | return 0; | |
758 | } | |
759 | ||
760 | Int_t | |
761 | AliITSUTrackerCooked::AliITSUlayer::FindClusterIndex(Double_t z) const { | |
762 | //-------------------------------------------------------------------- | |
763 | // This function returns the index of the first | |
764 | // with its fZ >= "z". | |
765 | //-------------------------------------------------------------------- | |
766 | if (fN==0) return 0; | |
767 | ||
768 | Int_t b=0; | |
769 | if (z <= fClusters[b]->GetZ()) return b; | |
770 | ||
771 | Int_t e=b+fN-1; | |
772 | if (z > fClusters[e]->GetZ()) return e+1; | |
773 | ||
774 | Int_t m=(b+e)/2; | |
775 | for (; b<e; m=(b+e)/2) { | |
776 | if (z > fClusters[m]->GetZ()) b=m+1; | |
777 | else e=m; | |
778 | } | |
779 | return m; | |
780 | } | |
781 | ||
4cce6dec | 782 | void AliITSUTrackerCooked::AliITSUlayer:: |
783 | SelectClusters(Float_t phi, Float_t dy, Float_t z, Float_t dz) { | |
66be9a4e | 784 | //-------------------------------------------------------------------- |
785 | // This function selects clusters within the "road" | |
786 | //-------------------------------------------------------------------- | |
4cce6dec | 787 | fNsel=0; |
788 | ||
789 | Float_t dphi=dy/fR; | |
790 | Float_t phiMin=phi-dphi; | |
791 | Float_t phiMax=phi+dphi; | |
792 | Float_t zMin=z-dz; | |
793 | Float_t zMax=z+dz; | |
794 | ||
795 | Int_t i=FindClusterIndex(zMin), imax=FindClusterIndex(zMax); | |
796 | for (; i<imax; i++) { | |
797 | Float_t cphi=fPhi[i]; | |
798 | if (cphi <= phiMin) continue; | |
799 | if (cphi > phiMax) continue; | |
66be9a4e | 800 | AliCluster *c=fClusters[i]; |
a00907de | 801 | if (c->IsClusterUsed()) continue; |
4cce6dec | 802 | |
66be9a4e | 803 | fIndex[fNsel++]=i; |
4cce6dec | 804 | if (fNsel>=kMaxSelected) break; |
805 | } | |
66be9a4e | 806 | } |
807 | ||
808 | const AliCluster *AliITSUTrackerCooked::AliITSUlayer::GetNextCluster(Int_t &ci){ | |
809 | //-------------------------------------------------------------------- | |
810 | // This function returns clusters within the "road" | |
811 | //-------------------------------------------------------------------- | |
4cce6dec | 812 | if (fI<fNsel) { |
813 | ci=fIndex[fI++]; | |
814 | return fClusters[ci]; | |
66be9a4e | 815 | } |
4cce6dec | 816 | ci=-1; |
817 | return 0; | |
66be9a4e | 818 | } |
819 |