Tracker update
[u/mrichter/AliRoot.git] / ITS / UPGRADE / AliITSUTrackerGlo.cxx
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   TO CHECK
18   GetBz usage
19
20  */
21
22 //-------------------------------------------------------------------------
23 //  Implementation of the ITS Upgrade tracker mother class.              
24 //-------------------------------------------------------------------------
25 #include <TTree.h>
26 #include <Riostream.h> 
27 #include <TMath.h>
28
29 #include "AliITSUTrackerGlo.h"
30 #include "AliESDEvent.h"
31 #include "AliESDtrack.h"
32 #include "AliITSURecoDet.h"
33 #include "AliITSURecoSens.h"
34 #include "AliITSUReconstructor.h"
35 #include "AliITSReconstructor.h"
36 #include "AliITSUSeed.h"
37 #include "AliITSUAux.h"
38 #include "AliITSUClusterPix.h"
39 using namespace AliITSUAux;
40 using namespace TMath;
41
42
43 ClassImp(AliITSUTrackerGlo)
44 //_________________________________________________________________________
45 AliITSUTrackerGlo::AliITSUTrackerGlo(AliITSUReconstructor* rec)
46 :  fReconstructor(rec)
47   ,fITS(0)
48   ,fCurrMass(kPionMass)
49   ,fSeedsLr(0)
50   ,fSeedsPool("AliITSUSeed",0)
51 {
52   // Default constructor
53   if (rec) Init(rec);
54 }
55
56 //_________________________________________________________________________
57 AliITSUTrackerGlo::~AliITSUTrackerGlo()
58 {
59  // Default destructor
60  //  
61   delete fITS;
62   delete[] fSeedsLr;
63   //
64 }
65
66 //_________________________________________________________________________
67 void AliITSUTrackerGlo::Init(AliITSUReconstructor* rec)
68 {
69   // init with external reconstructor
70   //
71   fITS = new AliITSURecoDet(rec->GetGeom(),"ITSURecoInterface");
72   for (int ilr=fITS->GetNLayersActive();ilr--;) {
73     fITS->GetLayerActive(ilr)->SetClusters(rec->GetClusters(ilr));
74   }
75   //
76   fSeedsPool.ExpandCreateFast(1000); // RS TOCHECK
77   int n = fITS->GetNLayersActive()+1;
78   fSeedsLr = new TObjArray[n];
79   //
80
81 }
82
83 //_________________________________________________________________________
84 Int_t AliITSUTrackerGlo::Clusters2Tracks(AliESDEvent *esdEv)
85 {
86   //
87   //
88   fITS->ProcessClusters();
89   // select ESD tracks to propagate
90   int nTrESD = esdEv->GetNumberOfTracks();
91   for (int itr=0;itr<nTrESD;itr++) {
92     AliESDtrack *esdTr = esdEv->GetTrack(itr);
93     FindTrack(esdTr);
94   }
95
96   return 0;
97 }
98
99 //_________________________________________________________________________
100 Int_t AliITSUTrackerGlo::PropagateBack(AliESDEvent * /*event*/)
101 {
102   //
103   // To be implemented 
104   //
105   
106  Info("PropagateBack","To be implemented");
107   return 0;
108 }
109
110 //_________________________________________________________________________
111 Int_t AliITSUTrackerGlo::RefitInward(AliESDEvent * /*event*/)
112 {
113   //
114   // To be implemented 
115   //
116   
117   Info("RefitInward","To be implemented");
118   return 0;
119 }
120
121 //_________________________________________________________________________
122 Int_t AliITSUTrackerGlo::LoadClusters(TTree * treeRP)
123 {
124   // read from tree (if pointer provided) or directly from the ITS reco interface
125   //
126   return fReconstructor->LoadClusters(treeRP);
127
128
129 //_________________________________________________________________________
130 void AliITSUTrackerGlo::UnloadClusters()
131 {
132   //
133   // To be implemented 
134   //
135   
136   Info("UnloadClusters","To be implemented");
137
138 //_________________________________________________________________________
139 AliCluster * AliITSUTrackerGlo::GetCluster(Int_t /*index*/) const
140 {
141   //
142   // To be implemented 
143   //
144   
145   Info("GetCluster","To be implemented");
146   return 0x0;
147
148
149 //_________________________________________________________________________
150 Bool_t AliITSUTrackerGlo::NeedToProlong(AliESDtrack* esdTr)
151 {
152   // do we need to match this track to ITS?
153   //
154   static double bz = GetBz();
155   if (!esdTr->IsOn(AliESDtrack::kTPCin) ||
156       esdTr->IsOn(AliESDtrack::kTPCout) ||
157       esdTr->IsOn(AliESDtrack::kITSin)  ||
158       esdTr->GetKinkIndex(0)>0) return kFALSE;
159   //
160   if (esdTr->Pt()<AliITSUReconstructor::GetRecoParam()->GetMinPtForProlongation()) return kFALSE;
161   //
162   float dtz[2];
163   esdTr->GetDZ(GetX(),GetY(),GetZ(),bz,dtz); 
164   // if track is not V0 candidata but has large offset wrt IP, reject it. RS TOCHECK
165   if ( !(esdTr->GetV0Index(0)>0 && dtz[0]>AliITSUReconstructor::GetRecoParam()->GetMaxDforV0dghtrForProlongation()) 
166        && (Abs(dtz[0])>AliITSUReconstructor::GetRecoParam()->GetMaxDForProlongation() ||
167            Abs(dtz[1])>AliITSUReconstructor::GetRecoParam()->GetMaxDZForProlongation())) return kFALSE;
168   //
169   return kTRUE;
170 }
171
172 //_________________________________________________________________________
173 void AliITSUTrackerGlo::FindTrack(AliESDtrack* esdTr)
174 {
175   // find prolongaion candidates finding for single seed
176   //
177   const double kTolerX = 5e-4;
178   //
179   if (!NeedToProlong(esdTr)) return;  // are we interested in this track?
180   if (!InitSeed(esdTr))      return;  // initialize prolongations hypotheses tree
181   //
182   AliITSURecoSens *hitSens[AliITSURecoSens::kNNeighbors];
183   AliITSUSeed seed0;
184   TObjArray clArr; // container for transfer of clusters matching to seed
185   //
186   for (int ila=fITS->GetNLayersActive();ila--;) {
187     int ilaUp = ila+1;                         // prolong seeds from layer above
188     int nSeedsUp = GetNSeeds(ilaUp);
189     for (int isd=0;isd<nSeedsUp;isd++) {
190       AliITSUSeed* seedU = GetSeed(ilaUp,isd);  // seed on prev.active layer to prolong
191       seed0 = *seedU;
192       // go till next active layer
193       if (!TransportToLayer(&seed0, fITS->GetLrIDActive(ilaUp), fITS->GetLrIDActive(ila)) ) {
194         //
195         // Check if the seed satisfies to track definition
196         if (NeedToKill(&seed0,kTransportFailed)) KillSeed(ilaUp,isd); 
197         continue; // RS TODO: decide what to do with tracks stopped on higher layers w/o killing
198       }
199       AliITSURecoLayer* lrA = fITS->GetLayerActive(ila);
200       if (!GetRoadWidth(&seed0, ila)) { // failed to find road width on the layer
201         if (NeedToKill(&seed0,kRWCheckFailed)) KillSeed(ilaUp,isd); 
202         continue;
203       }
204       int nsens = lrA->FindSensors(&fTrImpData[kTrPhi0], hitSens);  // find detectors which may be hit by the track (max 4)
205       //
206       for (int isn=nsens;isn--;) {
207         AliITSURecoSens* sens = hitSens[isn];
208         AliITSUSeed* seedT = NewSeedFromPool(&seed0);
209         seedT->SetParent(seedU);
210         if (!seedT->Propagate(sens->GetPhiTF(),sens->GetXTF(),GetBz())) {DeleteLastSeedFromPool(); continue;}
211         int clID    = sens->GetFirstClusterId();
212         for (int icl=sens->GetNClusters();icl--;) {
213           AliITSUClusterPix* cl = (AliITSUClusterPix*)lrA->GetCluster(clID); // cluster data is in the tracking frame
214           if (TMath::Abs(cl->GetX())>kTolerX) { // if due to the misalingment X is large, propagate track only
215             if (!seedT->PropagateParamOnlyTo(seedT->GetX()+cl->GetX(),GetBz())) {DeleteLastSeedFromPool(); continue;}
216           }
217           
218
219         }
220       }
221       
222     }
223   }
224   //
225   ResetSeedTree();
226 }
227
228 //_________________________________________________________________________
229 Bool_t AliITSUTrackerGlo::InitSeed(AliESDtrack *esdTr)
230 {
231   // init prolongaion candidates finding for single seed
232   fCurrMass = esdTr->GetMass();
233   if (fCurrMass<kPionMass*0.9) fCurrMass = kPionMass; // don't trust to mu, e identification from TPCin
234   //
235   AliITSUSeed* seed = NewSeedFromPool();
236   seed->AliExternalTrackParam::operator=(*esdTr);
237   seed->SetParent(esdTr);
238   fSeedsLr[fITS->GetNLayersActive()].AddLast(seed);
239   return kTRUE;
240   // TO DO
241 }
242
243 //_________________________________________________________________________
244 void AliITSUTrackerGlo::ResetSeedTree()
245 {
246   // reset current hypotheses tree
247   for (int i=fITS->GetNLayersActive()+1;i--;) fSeedsLr[fITS->GetNLayersActive()].Clear();
248 }
249
250 //_________________________________________________________________________
251 Bool_t AliITSUTrackerGlo::TransportToLayer(AliITSUSeed* seed, Int_t lFrom, Int_t lTo)
252 {
253   // transport seed from layerFrom to the entrance of layerTo
254   //  
255   const double kToler = 1e-6; // tolerance for layer on-surface check
256   //
257   int dir = lTo > lFrom ? 1:-1;
258   AliITSURecoLayer* lrFr = fITS->GetLayer(lFrom); // this can be 0 when extrapolation from TPC to ITS is requested
259   Bool_t checkFirst = kTRUE;
260   while(lFrom!=lTo) {
261     double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
262     if (lrFr) {
263       Bool_t doLayer = kTRUE;
264       double xToGo = dir>0 ? lrFr->GetRMax() : lrFr->GetRMin();
265       if (checkFirst) { // do we need to track till the surface of the current layer ?
266         checkFirst = kFALSE;
267         if      (dir>0) { if (curR2-xToGo*xToGo>kToler) doLayer = kFALSE; } // on the surface or outside of the layer
268         else if (dir<0) { if (xToGo*xToGo-curR2>kToler) doLayer = kFALSE; } // on the surface or outside of the layer
269       }
270       if (doLayer) {
271         if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
272         // go via layer to its boundary, applying material correction.
273         if (!PropagateTrackTo(seed,xToGo,fCurrMass, lrFr->GetMaxStep(), kFALSE, -1, 0, kTRUE)) return kFALSE;
274       }
275     }
276     AliITSURecoLayer* lrTo =  fITS->GetLayer( (lFrom+=dir) );
277     if (!lrTo) AliFatal(Form("Layer %d does not exist",lFrom));
278     //
279     // go the entrance of the layer, assuming no materials in between
280     double xToGo = dir>0 ? lrTo->GetRMin() : lrTo->GetRMax();
281     if (seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
282     if (!seed->PropagateTo(xToGo, GetBz())) return kFALSE; // RS: do we need BxByBz?
283     lrFr = lrTo;
284   }
285   return kTRUE;
286   //
287 }
288
289 //_________________________________________________________________________
290 Bool_t AliITSUTrackerGlo::GetRoadWidth(AliITSUSeed* seed, int ilrA)
291 {
292   // calculate road width in terms of phi and z for the track which MUST be on the external radius of the layer
293   // as well as some aux info
294   double bz = GetBz();
295   AliITSURecoLayer* lrA = fITS->GetLayerActive(ilrA);
296   seed->GetXYZ(&fTrImpData[kTrXIn]);    // lab position at the entrance from above
297   //
298   fTrImpData[kTrPhiIn] = ATan2(fTrImpData[kTrYIn],fTrImpData[kTrXIn]);
299   if (!seed->Rotate(fTrImpData[kTrPhiIn])) return kFALSE; // go to the frame of the entry point into the layer
300   double dr  = lrA->GetDR();                              // approximate X dist at the inner radius
301   if (!seed->GetXYZAt(seed->GetX()-dr, bz, fTrImpData + kTrXOut)) {
302     // special case: track does not reach inner radius, might be tangential
303     double r = seed->GetD(0,0,bz);
304     double x;
305     if (!seed->GetXatLabR(r,x,bz,-1)) {
306       AliError(Form("This should not happen: r=%f",r));
307       seed->Print();
308       return kFALSE;
309     }
310     dr = Abs(seed->GetX() - x);
311     if (!seed->GetXYZAt(x, bz, fTrImpData + kTrXOut)) {
312       AliError(Form("This should not happen: x=%f",x));
313       seed->Print();
314       return kFALSE;      
315     }
316   }
317   //
318   fTrImpData[kTrPhiOut] = ATan2(fTrImpData[kTrYOut],fTrImpData[kTrXOut]);
319   double sgy = seed->GetSigmaY2() + dr*dr*seed->GetSigmaSnp2() + AliITSUReconstructor::GetRecoParam()->GetSigmaY2(ilrA);
320   double sgz = seed->GetSigmaZ2() + dr*dr*seed->GetSigmaTgl2() + AliITSUReconstructor::GetRecoParam()->GetSigmaZ2(ilrA);
321   sgy = Sqrt(sgy)*AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadY();
322   sgz = Sqrt(sgz)*AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadZ();
323   fTrImpData[kTrPhi0] = 0.5*(fTrImpData[kTrPhiOut]+fTrImpData[kTrPhiIn]);
324   fTrImpData[kTrZ0]   = 0.5*(fTrImpData[kTrZOut]+fTrImpData[kTrPhiIn]);
325   fTrImpData[kTrDPhi] = 0.5*Abs(fTrImpData[kTrPhiOut]-fTrImpData[kTrPhiIn]) + sgy/lrA->GetR();
326   fTrImpData[kTrDZ]   = 0.5*Abs(fTrImpData[kTrZOut]-fTrImpData[kTrPhiIn])   + sgz;
327   //  
328   return kTRUE;
329 }
330
331 //_________________________________________________________________________
332 AliITSUSeed* AliITSUTrackerGlo::NewSeedFromPool(const AliITSUSeed* src)
333 {
334   // create new seed, optionally copying from the source
335   return src ? 
336     new(fSeedsPool[fSeedsPool.GetEntriesFast()]) AliITSUSeed(*src) :
337     new(fSeedsPool[fSeedsPool.GetEntriesFast()]) AliITSUSeed();
338 }