Updated version of the PID code (from B. Batyunya)
[u/mrichter/AliRoot.git] / ITS / AliITStrackerV2.cxx
CommitLineData
006b5f7f 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-------------------------------------------------------------------------
17// Implementation of the ITS tracker class
18//
19// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
20//-------------------------------------------------------------------------
21#include <TFile.h>
22#include <TTree.h>
23#include <TRandom.h>
24#include <iostream.h>
25
26#include "AliITSgeom.h"
27#include "AliITSRecPoint.h"
517b130f 28#include "AliTPCtrack.h"
006b5f7f 29#include "AliITSclusterV2.h"
30#include "AliITStrackerV2.h"
31
32//#define DEBUG
33
34#ifdef DEBUG
7f6ddf58 35Int_t LAB=7;
006b5f7f 36#endif
37
006b5f7f 38AliITStrackerV2::AliITSlayer AliITStrackerV2::fLayers[kMaxLayer]; // ITS layers
39
40AliITStrackerV2::
7f6ddf58 41AliITStrackerV2(const AliITSgeom *geom, Int_t eventn) throw (const Char_t *) :
42AliTracker() {
006b5f7f 43 //--------------------------------------------------------------------
44 //This is the AliITStracker constructor
45 //It also reads clusters from a root file
46 //--------------------------------------------------------------------
7f6ddf58 47 fEventN=eventn; //MI change add event number - used to generate identifier
006b5f7f 48
49 AliITSgeom *g=(AliITSgeom*)geom;
50
51 Float_t x,y,z;
52 Int_t i;
53 for (i=1; i<kMaxLayer+1; i++) {
54 Int_t nlad=g->GetNladders(i);
55 Int_t ndet=g->GetNdetectors(i);
56
57 g->GetTrans(i,1,1,x,y,z);
58 Double_t r=TMath::Sqrt(x*x + y*y);
59 Double_t poff=TMath::ATan2(y,x);
60 Double_t zoff=z;
61
62 g->GetTrans(i,1,2,x,y,z);
63 r += TMath::Sqrt(x*x + y*y);
64 g->GetTrans(i,2,1,x,y,z);
65 r += TMath::Sqrt(x*x + y*y);
66 g->GetTrans(i,2,2,x,y,z);
67 r += TMath::Sqrt(x*x + y*y);
68 r*=0.25;
69
70 new (fLayers+i-1) AliITSlayer(r,poff,zoff,nlad,ndet);
71
72 for (Int_t j=1; j<nlad+1; j++) {
73 for (Int_t k=1; k<ndet+1; k++) { //Fill this layer with detectors
74 Float_t x,y,zshift; g->GetTrans(i,j,k,x,y,zshift);
75 Double_t rot[9]; g->GetRotMatrix(i,j,k,rot);
76
77 Double_t r =-x*rot[1] + y*rot[0]; if (i==1) r=-r;
78 Double_t phi=TMath::ATan2(rot[1],rot[0]); if (i==1) phi-=3.1415927;
79 phi+=0.5*TMath::Pi();
80 AliITSdetector &det=fLayers[i-1].GetDetector((j-1)*ndet + k-1);
81
82if (phi<0) phi += 2*TMath::Pi();
83
84 new(&det) AliITSdetector(r,phi);
85 }
86 }
87
88 }
89
90 try {
91 //Read clusters
c7ee21ca 92 //MI change
93 char cname[100];
94 sprintf(cname,"TreeC_ITS_%d",eventn);
743a19f2 95 TTree *cTree=(TTree*)gDirectory->Get(cname);
743a19f2 96
006b5f7f 97 if (!cTree) throw
98 ("AliITStrackerV2::AliITStrackerV2 can't get cTree !\n");
99
100 TBranch *branch=cTree->GetBranch("Clusters");
101 if (!branch) throw
102 ("AliITStrackerV2::AliITStrackerV2 can't get Clusters branch !\n");
103
104 TClonesArray dummy("AliITSclusterV2",10000), *clusters=&dummy;
105 branch->SetAddress(&clusters);
106
107 Int_t nentr=(Int_t)cTree->GetEntries();
108 for (i=0; i<nentr; i++) {
109 if (!cTree->GetEvent(i)) continue;
14825d5a 110 Int_t lay,lad,det; g->GetModuleId(i,lay,lad,det);
006b5f7f 111 Int_t ncl=clusters->GetEntriesFast();
112 while (ncl--) {
113 AliITSclusterV2 *c=(AliITSclusterV2*)clusters->UncheckedAt(ncl);
114
115#ifdef DEBUG
116if (c->GetLabel(2)!=LAB)
117if (c->GetLabel(1)!=LAB)
118if (c->GetLabel(0)!=LAB) continue;
119cout<<lay-1<<' '<<lad-1<<' '<<det-1<<' '<<c->GetY()<<' '<<c->GetZ()<<endl;
120#endif
121
7f6ddf58 122 fLayers[lay-1].InsertCluster(new AliITSclusterV2(*c));
006b5f7f 123 }
124 clusters->Delete();
125 }
c7ee21ca 126 delete cTree; //Thanks to Mariana Bondila
006b5f7f 127 }
128 catch (const Char_t *msg) {
129 cerr<<msg<<endl;
130 throw;
131 }
132
133 fI=kMaxLayer;
7f6ddf58 134
135 fPass=0;
136 fConstraint[0]=1; fConstraint[1]=0;
006b5f7f 137}
138
14825d5a 139#ifdef DEBUG
006b5f7f 140static Int_t lbl;
14825d5a 141#endif
006b5f7f 142
143Int_t AliITStrackerV2::Clusters2Tracks(const TFile *inp, TFile *out) {
144 //--------------------------------------------------------------------
145 //This functions reconstructs ITS tracks
146 //--------------------------------------------------------------------
147 TFile *in=(TFile*)inp;
148 TDirectory *savedir=gDirectory;
149
150 if (!in->IsOpen()) {
151 cerr<<"AliITStrackerV2::Clusters2Tracks(): ";
152 cerr<<"file with TPC tracks is not open !\n";
153 return 1;
154 }
155
156 if (!out->IsOpen()) {
157 cerr<<"AliITStrackerV2::Clusters2Tracks(): ";
158 cerr<<"file for ITS tracks is not open !\n";
159 return 2;
160 }
161
7f6ddf58 162 in->cd();
163
743a19f2 164 char tname[100];
7f6ddf58 165 Int_t nentr=0; TObjArray itsTracks(10000);
166
167 {/* Read TPC tracks */
168 sprintf(tname,"TreeT_TPC_%d",fEventN);
169 TTree *tpcTree=(TTree*)in->Get(tname);
170 if (!tpcTree) {
171 cerr<<"AliITStrackerV2::Clusters2Tracks(): "
172 "can't get a tree with TPC tracks !\n";
173 return 3;
174 }
175 AliTPCtrack *itrack=new AliTPCtrack;
176 tpcTree->SetBranchAddress("tracks",&itrack);
177 nentr=(Int_t)tpcTree->GetEntries();
178 for (Int_t i=0; i<nentr; i++) {
179 tpcTree->GetEvent(i);
180 itsTracks.AddLast(new AliITStrackV2(*itrack));
181 }
182 delete tpcTree; //Thanks to Mariana Bondila
183 delete itrack;
006b5f7f 184 }
7f6ddf58 185 itsTracks.Sort();
006b5f7f 186
187 out->cd();
7f6ddf58 188
189 sprintf(tname,"TreeT_ITS_%d",fEventN);
190 TTree itsTree(tname,"Tree with ITS tracks");
006b5f7f 191 AliITStrackV2 *otrack=&fBestTrack;
192 itsTree.Branch("tracks","AliITStrackV2",&otrack,32000,0);
193
7f6ddf58 194 for (fPass=0; fPass<2; fPass++) {
195 Int_t &constraint=fConstraint[fPass]; if (constraint<0) continue;
196 for (Int_t i=0; i<nentr; i++) {
197 if (i%10==0) cerr<<nentr-i<<" \r";
198 AliITStrackV2 *t=(AliITStrackV2*)itsTracks.UncheckedAt(i);
199 if (t==0) continue; //this track has been already tracked
200 Int_t tpcLabel=t->GetLabel(); //save the TPC track label
006b5f7f 201#ifdef DEBUG
14825d5a 202lbl=tpcLabel;
006b5f7f 203if (TMath::Abs(tpcLabel)!=LAB) continue;
204cout<<tpcLabel<<" *****************\n";
205#endif
7f6ddf58 206 try {
207 ResetTrackToFollow(*t);
208 } catch (const Char_t *msg) {
209 cerr<<msg<<endl;
210 continue;
211 }
212 ResetBestTrack();
006b5f7f 213
7f6ddf58 214 Double_t r2=fTrackToFollow.GetX()*fTrackToFollow.GetX();
215 Double_t x0=0.082/21.82*2.33*(45*45+40*40)/r2+2.000/41.28*0.03*75*75/r2;
216 if (constraint) fTrackToFollow.Improve(x0,GetY(),GetZ());
217
218 Double_t xk=80.;
219 fTrackToFollow.PropagateTo(xk,0.,0.); //Ne if it's still there
220
221 xk-=0.005;
222 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
223 xk-=0.02;
224 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
225 xk-=2.0;
226 fTrackToFollow.PropagateTo(xk, 2.0/41.28*0.029, 2.0*0.029);//Nomex
227 xk-=0.02;
228 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
229 xk-=0.005;
230 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
231
232 xk=61.;
233 fTrackToFollow.PropagateTo(xk,0.,0.); //C02
234
235 xk -=0.005;
236 fTrackToFollow.PropagateTo(xk, 0.005/24.01*2.7, 0.005*2.7); //Al
237 xk -=0.005;
238 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
239 xk -=0.02;
240 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
241 xk -=0.5;
242 fTrackToFollow.PropagateTo(xk, 0.5/41.28*.029, 0.5*0.029); //Nomex
243 xk -=0.02;
244 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
245 xk -=0.005;
246 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
247 xk -=0.005;
248 fTrackToFollow.PropagateTo(xk, 0.005/24.01*2.7, 0.005*2.7); //Al
249
250
251 for (FollowProlongation(); fI<kMaxLayer; fI++) {
252 while (TakeNextProlongation()) FollowProlongation();
253 }
006b5f7f 254
255#ifdef DEBUG
256cout<<fBestTrack.GetNumberOfClusters()<<" number of clusters\n\n";
257#endif
258
7f6ddf58 259 if (fBestTrack.GetNumberOfClusters() >= kMaxLayer-kLayersToSkip) {
260 fBestTrack.SetLabel(tpcLabel);
261 fBestTrack.CookdEdx();
262 CookLabel(&fBestTrack,0.); //For comparison only
263 itsTree.Fill();
264 UseClusters(&fBestTrack);
265 delete itsTracks.RemoveAt(i);
266 }
006b5f7f 267
7f6ddf58 268 }
006b5f7f 269 }
006b5f7f 270
7f6ddf58 271 itsTree.Write();
c7ee21ca 272
7f6ddf58 273 itsTracks.Delete();
274 savedir->cd();
275 cerr<<"Number of TPC tracks: "<<nentr<<endl;
276 cerr<<"Number of prolonged tracks: "<<itsTree.GetEntries()<<endl;
006b5f7f 277
278 return 0;
279}
280
14825d5a 281Int_t AliITStrackerV2::PropagateBack(const TFile *inp, TFile *out) {
282 //--------------------------------------------------------------------
283 //This functions propagates reconstructed ITS tracks back
284 //--------------------------------------------------------------------
285 TFile *in=(TFile*)inp;
286 TDirectory *savedir=gDirectory;
287
288 if (!in->IsOpen()) {
289 cerr<<"AliITStrackerV2::PropagateBack(): ";
290 cerr<<"file with ITS tracks is not open !\n";
291 return 1;
292 }
293
294 if (!out->IsOpen()) {
295 cerr<<"AliITStrackerV2::PropagateBack(): ";
296 cerr<<"file for back propagated ITS tracks is not open !\n";
297 return 2;
298 }
299
300 in->cd();
c7ee21ca 301 TTree *itsTree=(TTree*)in->Get("TreeT_ITS_0");
14825d5a 302 if (!itsTree) {
303 cerr<<"AliITStrackerV2::PropagateBack() ";
304 cerr<<"can't get a tree with ITS tracks !\n";
305 return 3;
306 }
307 AliITStrackV2 *itrack=new AliITStrackV2;
308 itsTree->SetBranchAddress("tracks",&itrack);
309
310 out->cd();
c7ee21ca 311 TTree backTree("TreeT_ITSb_0","Tree with back propagated ITS tracks");
14825d5a 312 AliTPCtrack *otrack=0;
313 backTree.Branch("tracks","AliTPCtrack",&otrack,32000,0);
314
315 Int_t ntrk=0;
316
317 Int_t nentr=(Int_t)itsTree->GetEntries();
318 for (Int_t i=0; i<nentr; i++) {
319 itsTree->GetEvent(i);
320 ResetTrackToFollow(*itrack);
321 fTrackToFollow.ResetCovariance(); fTrackToFollow.ResetClusters();
322 Int_t itsLabel=fTrackToFollow.GetLabel(); //save the ITS track label
323
324#ifdef DEBUG
325if (TMath::Abs(itsLabel)!=LAB) continue;
326cout<<itsLabel<<" *****************\n";
327#endif
328
329 try {
330 Int_t nc=itrack->GetNumberOfClusters();
331#ifdef DEBUG
332for (Int_t k=0; k<nc; k++) {
333 Int_t index=itrack->GetClusterIndex(k);
334 AliITSclusterV2 *c=(AliITSclusterV2*)GetCluster(index);
335 cout<<c->GetLabel(0)<<' '<<c->GetLabel(1)<<' '<<c->GetLabel(2)<<endl;
336}
337#endif
338 Int_t idx=0, l=0;
339 const AliITSclusterV2 *c=0;
340 if (nc--) {
341 idx=itrack->GetClusterIndex(nc); l=(idx&0xf0000000)>>28;
342 c=(AliITSclusterV2*)GetCluster(idx);
343 }
344 for (fI=0; fI<kMaxLayer; fI++) {
345 AliITSlayer &layer=fLayers[fI];
346 Double_t r=layer.GetR();
347 if (fI==2 || fI==4) {
348 Double_t rs=0.5*(fLayers[fI-1].GetR() + r);
349 Double_t ds=0.034; if (fI==4) ds=0.039;
350 Double_t dx0r=ds/21.82*2.33, dr=ds*2.33;
351 if (!fTrackToFollow.PropagateTo(rs,1*dx0r,dr)) throw "";
352 }
353
354 Double_t x,y,z;
355 if (!fTrackToFollow.GetGlobalXYZat(r,x,y,z))
356 throw "AliITStrackerV2::PropagateBack: failed to estimate track !";
357 Double_t phi=TMath::ATan2(y,x);
358 Double_t d=layer.GetThickness(phi,z);
359 Int_t idet=layer.FindDetectorIndex(phi,z);
360 if (idet<0)
361 throw "AliITStrackerV2::PropagateBack: failed to find a detector !\n";
362 const AliITSdetector &det=layer.GetDetector(idet);
363 r=det.GetR(); phi=det.GetPhi();
364 if (!fTrackToFollow.Propagate(phi,r,d/21.82*2.33,d*2.33)) throw "";
365
366 const AliITSclusterV2 *cl=0;
367 Int_t index=0;
368 Double_t maxchi2=kMaxChi2;
369
370 if (l==fI) {
371 if (idet != c->GetDetectorIndex()) {
372 idet=c->GetDetectorIndex();
373 const AliITSdetector &det=layer.GetDetector(idet);
374 r=det.GetR(); phi=det.GetPhi();
375 if (!fTrackToFollow.Propagate(phi,r,0.,0.)) throw "";
376 }
377 Double_t chi2=fTrackToFollow.GetPredictedChi2(c);
378 if (chi2<kMaxChi2) {
379 cl=c; maxchi2=chi2; index=idx;
380 }
381 if (nc--) {
382 idx=itrack->GetClusterIndex(nc); l=(idx&0xf0000000)>>28;
383 c=(AliITSclusterV2*)GetCluster(idx);
384 }
385 }
386
387 if (fTrackToFollow.GetNumberOfClusters()>2) {
388 Double_t dz=3*TMath::Sqrt(fTrackToFollow.GetSigmaZ2()+kSigmaZ2[fI]);
389 Double_t dy=3*TMath::Sqrt(fTrackToFollow.GetSigmaY2()+kSigmaY2[fI]);
390 Double_t zmin=fTrackToFollow.GetZ() - dz;
391 Double_t zmax=fTrackToFollow.GetZ() + dz;
392 Double_t ymin=fTrackToFollow.GetY() + phi*r - dy;
393 Double_t ymax=fTrackToFollow.GetY() + phi*r + dy;
394 layer.SelectClusters(zmin,zmax,ymin,ymax);
395
396 const AliITSclusterV2 *cc=0; Int_t ci;
397 while ((cc=layer.GetNextCluster(ci))!=0) {
398 if (idet != cc->GetDetectorIndex()) continue;
399 Double_t chi2=fTrackToFollow.GetPredictedChi2(cc);
400 if (chi2<maxchi2) {
401 cl=cc; index=(fI<<28)+ci; maxchi2=chi2;
402 }
403 }
404 }
405
406 if (cl) {
407 if (!fTrackToFollow.Update(cl,maxchi2,index))
408 cerr<<"AliITStrackerV2::PropagateBack: filtering failed !\n";
409 continue;
410 }
411 }
412
413 Double_t xk=61.;
414 fTrackToFollow.PropagateTo(xk,0.,0.); //Air
415
416 xk +=0.005;
417 fTrackToFollow.PropagateTo(xk, 0.005/24.01*2.7, 0.005*2.7); //Al
418 xk +=0.005;
419 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
420 xk +=0.02;
421 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
422 xk +=0.5;
423 fTrackToFollow.PropagateTo(xk, 0.5/41.28*.029, 0.5*0.029); //Nomex
424 xk +=0.02;
425 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
426 xk +=0.005;
427 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
428 xk +=0.005;
429 fTrackToFollow.PropagateTo(xk, 0.005/24.01*2.7, 0.005*2.7); //Al
430
431 xk=80.;
432 fTrackToFollow.PropagateTo(xk,0.,0.); //CO2
433
434 xk+=0.005;
435 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
436 xk+=0.02;
437 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
438 xk+=2.0;
439 fTrackToFollow.PropagateTo(xk, 2.0/41.28*0.029, 2.0*0.029);//Nomex
440 xk+=0.02;
441 fTrackToFollow.PropagateTo(xk, 0.02/44.86*1.45, 0.02*1.45); //Kevlar
442 xk+=0.005;
443 fTrackToFollow.PropagateTo(xk, 0.005/44.77*1.71, 0.005*1.71); //Tedlar
444
445 fTrackToFollow.SetLabel(itsLabel);
446 otrack=new AliTPCtrack(fTrackToFollow,fTrackToFollow.GetAlpha());
447 backTree.Fill(); delete otrack;
448 UseClusters(&fTrackToFollow);
449 cerr << ++ntrk << " \r";
450 }
451 catch (const Char_t *msg) {
452 cerr<<msg<<endl;
453 }
454 }
455
456 backTree.Write();
457 savedir->cd();
458 cerr<<"Number of ITS tracks: "<<nentr<<endl;
459 cerr<<"Number of back propagated ITS tracks: "<<ntrk<<endl;
460
461 delete itrack;
462
c7ee21ca 463 delete itsTree; //Thanks to Mariana Bondila
464
14825d5a 465 return 0;
466}
006b5f7f 467
468AliCluster *AliITStrackerV2::GetCluster(Int_t index) const {
469 //--------------------------------------------------------------------
470 // Return pointer to a given cluster
471 //--------------------------------------------------------------------
472 Int_t l=(index & 0xf0000000) >> 28;
473 Int_t c=(index & 0x0fffffff) >> 00;
474 return fLayers[l].GetCluster(c);
475}
476
477
478void AliITStrackerV2::FollowProlongation() {
479 //--------------------------------------------------------------------
480 //This function finds a track prolongation
481 //--------------------------------------------------------------------
482 Int_t tryAgain=kLayersToSkip;
483
484 while (fI) {
485 fI--;
006b5f7f 486#ifdef DEBUG
487cout<<fI<<' ';
488#endif
006b5f7f 489 AliITSlayer &layer=fLayers[fI];
490 AliITStrackV2 &track=fTracks[fI];
491
14825d5a 492 Double_t r=layer.GetR();
006b5f7f 493 if (fI==3 || fI==1) {
14825d5a 494 Double_t rs=0.5*(fLayers[fI+1].GetR() + r);
006b5f7f 495 Double_t ds=0.034; if (fI==3) ds=0.039;
496 Double_t dx0r=ds/21.82*2.33, dr=ds*2.33;
7f6ddf58 497 fTrackToFollow.Propagate(fTrackToFollow.GetAlpha(),rs,dx0r,dr);
006b5f7f 498 }
499
500 //find intersection
501 Double_t x,y,z;
14825d5a 502 if (!fTrackToFollow.GetGlobalXYZat(r,x,y,z)) {
7f6ddf58 503 //cerr<<"AliITStrackerV2::FollowProlongation: "
504 //"failed to estimate track !\n";
006b5f7f 505 break;
506 }
507 Double_t phi=TMath::ATan2(y,x);
508 Double_t d=layer.GetThickness(phi,z);
509 Int_t idet=layer.FindDetectorIndex(phi,z);
510 if (idet<0) {
7f6ddf58 511 //cerr<<"AliITStrackerV2::FollowProlongation: "
512 //"failed to find a detector !\n";
006b5f7f 513 break;
514 }
515
516 //propagate to the intersection
517 const AliITSdetector &det=layer.GetDetector(idet);
14825d5a 518 phi=det.GetPhi();
7f6ddf58 519 if (!fTrackToFollow.Propagate(phi,det.GetR(),d/21.82*2.33,d*2.33)) {
520 //cerr<<"AliITStrackerV2::FollowProlongation: propagation failed !\n";
006b5f7f 521 break;
522 }
7f6ddf58 523 if (TMath::Abs(fTrackToFollow.GetZ()-GetZ()) > r+1) break;
006b5f7f 524 fTrackToFollow.SetDetectorIndex(idet);
525
526 //Select possible prolongations and store the current track estimation
527 track.~AliITStrackV2(); new(&track) AliITStrackV2(fTrackToFollow);
528 Double_t dz=3*TMath::Sqrt(track.GetSigmaZ2() + kSigmaZ2[fI]);
14825d5a 529 if (dz < 0.5*TMath::Abs(track.GetTgl())) dz=0.5*TMath::Abs(track.GetTgl());
7f6ddf58 530 if (dz > kMaxRoad) {
006b5f7f 531 //cerr<<"AliITStrackerV2::FollowProlongation: too broad road in Z !\n";
532 break;
533 }
534 Double_t dy=4*TMath::Sqrt(track.GetSigmaY2() + kSigmaY2[fI]);
14825d5a 535 if (dy < 0.5*TMath::Abs(track.GetSnp())) dy=0.5*TMath::Abs(track.GetSnp());
7f6ddf58 536 if (dy > kMaxRoad) {
006b5f7f 537 //cerr<<"AliITStrackerV2::FollowProlongation: too broad road in Y !\n";
538 break;
539 }
7f6ddf58 540 Double_t zmin=track.GetZ() - dz;
006b5f7f 541 Double_t zmax=track.GetZ() + dz;
14825d5a 542 Double_t ymin=track.GetY() + r*phi - dy;
543 Double_t ymax=track.GetY() + r*phi + dy;
006b5f7f 544 layer.SelectClusters(zmin,zmax,ymin,ymax);
545
006b5f7f 546 //take another prolongation
547 if (!TakeNextProlongation()) if (!tryAgain--) break;
548 tryAgain=kLayersToSkip;
549
550 }
551
552 //deal with the best track
553 Int_t ncl=fTrackToFollow.GetNumberOfClusters();
554 Int_t nclb=fBestTrack.GetNumberOfClusters();
555 if (ncl)
556 if (ncl >= nclb) {
557 Double_t chi2=fTrackToFollow.GetChi2();
558 if (chi2/ncl < kChi2PerCluster) {
559 if (ncl > nclb || chi2 < fBestTrack.GetChi2()) {
560 ResetBestTrack();
561 }
562 }
563 }
564
565 if (fI) fI++;
566}
567
568
569Int_t AliITStrackerV2::TakeNextProlongation() {
570 //--------------------------------------------------------------------
571 //This function takes another track prolongation
572 //--------------------------------------------------------------------
573 //Double_t m[20];
574 Double_t d=GetEffectiveThickness(0,0); //Think of this !!!!
575
576 AliITSlayer &layer=fLayers[fI];
577 AliITStrackV2 &t=fTracks[fI];
7f6ddf58 578 Int_t &constraint=fConstraint[fPass];
006b5f7f 579
580 Double_t dz=4*TMath::Sqrt(t.GetSigmaZ2() + kSigmaZ2[fI]);
581 Double_t dy=4*TMath::Sqrt(t.GetSigmaY2() + kSigmaY2[fI]);
582
583 const AliITSclusterV2 *c=0; Int_t ci=-1;
584 Double_t chi2=12345.;
585 while ((c=layer.GetNextCluster(ci))!=0) {
14825d5a 586
587#ifdef DEBUG
588//if (fI==0)
589//if (c->GetLabel(0)!=TMath::Abs(lbl)) continue;
590#endif
591
006b5f7f 592 Int_t idet=c->GetDetectorIndex();
593
594 if (t.GetDetectorIndex()!=idet) {
595 const AliITSdetector &det=layer.GetDetector(idet);
596 if (!t.Propagate(det.GetPhi(),det.GetR(),0.,0.)) {
7f6ddf58 597 //cerr<<"AliITStrackerV2::TakeNextProlongation: "
598 //"propagation failed !\n";
006b5f7f 599 continue;
600 }
601 t.SetDetectorIndex(idet);
602
603#ifdef DEBUG
604cout<<fI<<" change detector !\n";
605#endif
606
607 }
608
609 if (TMath::Abs(t.GetZ() - c->GetZ()) > dz) continue;
610 if (TMath::Abs(t.GetY() - c->GetY()) > dy) continue;
611
612 //m[0]=fYV; m[1]=fZV;
613 //chi2=t.GetPredictedChi2(c,m,d/21.82*2.33);
614 chi2=t.GetPredictedChi2(c);
615
616 if (chi2<kMaxChi2) break;
617 }
618
619#ifdef DEBUG
620cout<<fI<<" chi2="<<chi2<<' '<<t.GetY()<<' '<<t.GetZ()<<' '<<dy<<' '<<dz<<endl;
621#endif
622
623 if (chi2>=kMaxChi2) return 0;
624 if (!c) return 0;
625
626 ResetTrackToFollow(t);
627
628 //if (!fTrackToFollow.Update(m,chi2,(fI<<28)+ci)) {
629 if (!fTrackToFollow.Update(c,chi2,(fI<<28)+ci)) {
7f6ddf58 630 //cerr<<"AliITStrackerV2::TakeNextProlongation: filtering failed !\n";
006b5f7f 631 return 0;
632 }
23efe5f1 633 // b.b. change for the dEdX
634 fTrackToFollow.SetSampledEdx(c->GetQ(),
635 fTrackToFollow.GetNumberOfClusters()-1);
636
7f6ddf58 637 if (constraint) fTrackToFollow.Improve(d/21.82*2.33,GetY(),GetZ());
006b5f7f 638
639#ifdef DEBUG
640cout<<"accepted lab="<<c->GetLabel(0)<<' '
641 <<fTrackToFollow.GetNumberOfClusters()<<' '
642 <<fTrackToFollow.GetY()<<' '<<fTrackToFollow.GetZ()<<endl<<endl;
643#endif
644
645 return 1;
646}
647
648
649
650
651AliITStrackerV2::AliITSlayer::AliITSlayer() {
652 //--------------------------------------------------------------------
653 //default AliITSlayer constructor
654 //--------------------------------------------------------------------
655 fN=0;
656 fDetectors=0;
657}
658
659AliITStrackerV2::AliITSlayer::
660AliITSlayer(Double_t r,Double_t p,Double_t z,Int_t nl,Int_t nd) {
661 //--------------------------------------------------------------------
662 //main AliITSlayer constructor
663 //--------------------------------------------------------------------
664 fR=r; fPhiOffset=p; fZOffset=z;
665 fNladders=nl; fNdetectors=nd;
666 fDetectors=new AliITSdetector[fNladders*fNdetectors];
667
668 fN=0;
669 fI=0;
670}
671
672AliITStrackerV2::AliITSlayer::~AliITSlayer() {
673 //--------------------------------------------------------------------
674 // AliITSlayer destructor
675 //--------------------------------------------------------------------
676 delete[] fDetectors;
677 for (Int_t i=0; i<fN; i++) delete fClusters[i];
678}
679
680Int_t AliITStrackerV2::AliITSlayer::InsertCluster(AliITSclusterV2 *c) {
681 //--------------------------------------------------------------------
682 //This function adds a cluster to this layer
683 //--------------------------------------------------------------------
684 if (fN==kMaxClusterPerLayer) {
685 cerr<<"AliITStrackerV2::AliITSlayer::InsertCluster(): Too many clusters !\n";
686 return 1;
687 }
688
689 if (fN==0) {fClusters[fN++]=c; return 0;}
690 Int_t i=FindClusterIndex(c->GetZ());
691 memmove(fClusters+i+1 ,fClusters+i,(fN-i)*sizeof(AliITSclusterV2*));
692 fClusters[i]=c; fN++;
693
694 return 0;
695}
696
697Int_t AliITStrackerV2::AliITSlayer::FindClusterIndex(Double_t z) const {
698 //--------------------------------------------------------------------
699 // This function returns the index of the nearest cluster
700 //--------------------------------------------------------------------
701 if (fN==0) return 0;
702 if (z <= fClusters[0]->GetZ()) return 0;
703 if (z > fClusters[fN-1]->GetZ()) return fN;
704 Int_t b=0, e=fN-1, m=(b+e)/2;
705 for (; b<e; m=(b+e)/2) {
706 if (z > fClusters[m]->GetZ()) b=m+1;
707 else e=m;
708 }
709 return m;
710}
711
712void AliITStrackerV2::AliITSlayer::
713SelectClusters(Double_t zmin,Double_t zmax,Double_t ymin, Double_t ymax) {
714 //--------------------------------------------------------------------
715 // This function sets the "window"
716 //--------------------------------------------------------------------
717 fI=FindClusterIndex(zmin); fZmax=zmax;
14825d5a 718 Double_t circle=2*TMath::Pi()*fR;
719 if (ymax>circle) { ymax-=circle; ymin-=circle; }
006b5f7f 720 fYmin=ymin; fYmax=ymax;
721}
722
723const AliITSclusterV2 *AliITStrackerV2::AliITSlayer::GetNextCluster(Int_t &ci){
724 //--------------------------------------------------------------------
725 // This function returns clusters within the "window"
726 //--------------------------------------------------------------------
727 const AliITSclusterV2 *cluster=0;
728 for (Int_t i=fI; i<fN; i++) {
729 const AliITSclusterV2 *c=fClusters[i];
730 if (c->GetZ() > fZmax) break;
731 if (c->IsUsed()) continue;
732 const AliITSdetector &det=GetDetector(c->GetDetectorIndex());
733 Double_t y=fR*det.GetPhi() + c->GetY();
734
735 if (y>2.*fR*TMath::Pi()) y -= 2*fR*TMath::Pi();
736 if (y>1.*fR*TMath::Pi() && fYmax<y) y -= 2*fR*TMath::Pi();
737
738 if (y<fYmin) continue;
739 if (y>fYmax) continue;
740 cluster=c; ci=i;
741 fI=i+1;
742 break;
743 }
7f6ddf58 744
006b5f7f 745 return cluster;
746}
747
748Int_t AliITStrackerV2::AliITSlayer::
749FindDetectorIndex(Double_t phi, Double_t z) const {
750 //--------------------------------------------------------------------
751 //This function finds the detector crossed by the track
752 //--------------------------------------------------------------------
753 Double_t dphi=phi-fPhiOffset;
754 if (dphi < 0) dphi += 2*TMath::Pi();
755 else if (dphi >= 2*TMath::Pi()) dphi -= 2*TMath::Pi();
756 Int_t np=Int_t(dphi*fNladders*0.5/TMath::Pi()+0.5);
757
758 Double_t dz=fZOffset-z;
759 Int_t nz=Int_t(dz*(fNdetectors-1)*0.5/fZOffset+0.5);
760
761 if (np>=fNladders) np-=fNladders;
762 if (np<0) np+=fNladders;
763
764#ifdef DEBUG
765cout<<np<<' '<<nz<<endl;
766#endif
767
768 return np*fNdetectors + nz;
769}
770
771Double_t
772AliITStrackerV2::AliITSlayer::GetThickness(Double_t phi, Double_t z) const {
773 //--------------------------------------------------------------------
774 //This function returns the thickness of this layer
775 //--------------------------------------------------------------------
776 //-pi<phi<+pi
7f6ddf58 777 if (3 <fR&&fR<8 ) return 2.5*0.096;
006b5f7f 778 if (13<fR&&fR<26) return 1.1*0.088;
779 if (37<fR&&fR<41) return 1.1*0.085;
780 return 1.1*0.081;
781}
782
783
784Double_t AliITStrackerV2::GetEffectiveThickness(Double_t phi,Double_t z) const
785{
786 //--------------------------------------------------------------------
787 //Returns the thickness between the current layer and the vertex
788 //--------------------------------------------------------------------
789 //-pi<phi<+pi
790 Double_t d=0.1;
791
792 Double_t xn=fLayers[fI].GetR();
793 for (Int_t i=0; i<fI; i++) {
794 Double_t xi=fLayers[i].GetR();
795 d+=fLayers[i].GetThickness(phi,z)*xi*xi;
796 }
797
798 if (fI>1) {
799 Double_t xi=0.5*(fLayers[1].GetR()+fLayers[2].GetR());
800 d+=0.034*xi*xi;
801 }
802
803 if (fI>3) {
804 Double_t xi=0.5*(fLayers[3].GetR()+fLayers[4].GetR());
805 d+=0.039*xi*xi;
806 }
807 return d/(xn*xn);
808}
809
810
811
812Int_t AliITStrackerV2::AliITSlayer::InRoad() const {
813 //--------------------------------------------------------------------
814 // This function returns number of clusters within the "window"
815 //--------------------------------------------------------------------
816 Int_t ncl=0;
817 for (Int_t i=fI; i<fN; i++) {
818 const AliITSclusterV2 *c=fClusters[i];
819 if (c->GetZ() > fZmax) break;
820 //if (c->IsUsed()) continue;
821 const AliITSdetector &det=GetDetector(c->GetDetectorIndex());
822 Double_t y=fR*det.GetPhi() + c->GetY();
823
824 if (y>2.*fR*TMath::Pi()) y -= 2*fR*TMath::Pi();
825 if (y>1.*fR*TMath::Pi() && fYmax<y) y -= 2*fR*TMath::Pi();
826
827 if (y<fYmin) continue;
828 if (y>fYmax) continue;
829 ncl++;
830 }
831 return ncl;
832}
833
834
835