AliITStrackerSA.cxx: minor fix for cosmic reco + adding debug prints with AliLog...
[u/mrichter/AliRoot.git] / ITS / AliITStrackerSA.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2003, 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 /* $Id$ */
17
18 ////////////////////////////////////////////////////
19 //  Stand alone tracker class                     //
20 //  Origin:  Elisabetta Crescio                   //
21 //  e-mail:  crescio@to.infn.it                   //
22 //  tracks are saved as AliITStrackV2 objects     //
23 ////////////////////////////////////////////////////
24
25 #include <stdlib.h>
26
27 #include <TArrayI.h>
28 #include <TBranch.h>
29 #include <TObjArray.h>
30 #include <TTree.h>
31
32 #include "AliESDEvent.h"
33 #include "AliESDVertex.h"
34 #include "AliESDtrack.h"
35 #include "AliITSVertexer.h"
36 #include "AliITSclusterTable.h"
37 #include "AliITSRecPoint.h"
38 #include "AliITSgeomTGeo.h"
39 #include "AliITStrackSA.h"
40 #include "AliITStrackerSA.h"
41 #include "AliITSReconstructor.h"
42 #include "AliLog.h"
43 #include "AliRun.h"
44
45 ClassImp(AliITStrackerSA)
46
47 //____________________________________________________________________________
48 AliITStrackerSA::AliITStrackerSA():AliITStrackerMI(),
49 fPhiEstimate(0),
50 fITSStandAlone(0),
51 fLambdac(0),
52 fPhic(0),
53 fCoef1(0),
54 fCoef2(0),
55 fCoef3(0),
56 fNloop(0),
57 fPhiWin(0),
58 fLambdaWin(0),
59 fVert(0),
60 fVertexer(0),
61 fListOfTracks(0),
62 fListOfSATracks(0),
63 fITSclusters(0),
64 fInwardFlag(0),
65 fOuterStartLayer(0),
66 fInnerStartLayer(5),
67 fMinNPoints(0),
68 fMinQ(0.),
69 fCluLayer(0),
70 fCluCoord(0){
71   // Default constructor
72   Init();
73  
74 }
75 //____________________________________________________________________________
76 AliITStrackerSA::AliITStrackerSA(const Char_t *geom):AliITStrackerMI(0),
77 fPhiEstimate(0),
78 fITSStandAlone(0),
79 fLambdac(0),
80 fPhic(0),
81 fCoef1(0),
82 fCoef2(0),
83 fCoef3(0),
84 fNloop(0),
85 fPhiWin(0),
86 fLambdaWin(0),
87 fVert(0),
88 fVertexer(0),
89 fListOfTracks(0),
90 fListOfSATracks(0),
91 fITSclusters(0),
92 fInwardFlag(0),
93 fOuterStartLayer(0),
94 fInnerStartLayer(5),
95 fMinNPoints(0),
96 fMinQ(0.),
97 fCluLayer(0),
98 fCluCoord(0) 
99 {
100   // Standard constructor (Vertex is known and passed to this obj.)
101   if (geom) {
102     AliWarning("\"geom\" is actually a dummy argument !");
103   }
104
105   Init();
106   fVert = 0;
107  
108 }
109
110 //____________________________________________________________________________
111 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliESDVertex *vert):AliITStrackerMI(0),
112 fPhiEstimate(0),
113 fITSStandAlone(0),
114 fLambdac(0),
115 fPhic(0),
116 fCoef1(0),
117 fCoef2(0),
118 fCoef3(0),
119 fNloop(0),
120 fPhiWin(0),
121 fLambdaWin(0),
122 fVert(vert),
123 fVertexer(0),
124 fListOfTracks(0),
125 fListOfSATracks(0),
126 fITSclusters(0),
127 fInwardFlag(0),
128 fOuterStartLayer(0),
129 fInnerStartLayer(5),
130 fMinNPoints(0),
131 fMinQ(0.),
132 fCluLayer(0),
133 fCluCoord(0)
134 {
135   // Standard constructor (Vertex is known and passed to this obj.)
136   if (geom) {
137     AliWarning("\"geom\" is actually a dummy argument !");
138   }
139   Init();
140  
141 }
142
143 //____________________________________________________________________________
144 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliITSVertexer *vertexer):AliITStrackerMI(0),
145 fPhiEstimate(0),
146 fITSStandAlone(0),
147 fLambdac(0),
148 fPhic(0),
149 fCoef1(0),
150 fCoef2(0),
151 fCoef3(0),
152 fNloop(0),
153 fPhiWin(0),
154 fLambdaWin(0),
155 fVert(),
156 fVertexer(vertexer),
157 fListOfTracks(0),
158 fListOfSATracks(0),
159 fITSclusters(0),
160 fInwardFlag(0),
161 fOuterStartLayer(0),
162 fInnerStartLayer(5),
163 fMinNPoints(0),
164 fMinQ(0.),
165 fCluLayer(0),
166 fCluCoord(0)
167 {
168   // Standard constructor (Vertex is unknown - vertexer is passed to this obj)
169   if (geom) {
170     AliWarning("\"geom\" is actually a dummy argument !");
171   }
172   Init();
173   fVertexer = vertexer;
174  
175 }
176
177 //____________________________________________________________________________
178 AliITStrackerSA::AliITStrackerSA(const AliITStrackerSA& tracker):AliITStrackerMI(),
179 fPhiEstimate(tracker.fPhiEstimate),
180 fITSStandAlone(tracker.fITSStandAlone),
181 fLambdac(tracker.fLambdac),
182 fPhic(tracker.fPhic),
183 fCoef1(tracker.fCoef1),
184 fCoef2(tracker.fCoef2),
185 fCoef3(tracker.fCoef3),
186 fNloop(tracker.fNloop),
187 fPhiWin(tracker.fPhiWin),
188 fLambdaWin(tracker.fLambdaWin),
189 fVert(tracker.fVert),
190 fVertexer(tracker.fVertexer),
191 fListOfTracks(tracker.fListOfTracks),
192 fListOfSATracks(tracker.fListOfSATracks),
193 fITSclusters(tracker.fITSclusters),
194 fInwardFlag(tracker.fInwardFlag),
195 fOuterStartLayer(tracker.fOuterStartLayer),
196 fInnerStartLayer(tracker.fInnerStartLayer),
197 fMinNPoints(tracker.fMinNPoints),
198 fMinQ(tracker.fMinQ),
199 fCluLayer(tracker.fCluLayer),
200 fCluCoord(tracker.fCluCoord) {
201   // Copy constructor
202   for(Int_t i=0;i<2;i++){
203     fPoint1[i]=tracker.fPoint1[i];
204     fPoint2[i]=tracker.fPoint2[i];
205     fPoint3[i]=tracker.fPoint3[i];
206     fPointc[i]=tracker.fPointc[i];
207   }
208   if(tracker.fVertexer && tracker.fVert){
209     fVert = new AliESDVertex(*tracker.fVert);
210   }
211   else {
212     fVert = tracker.fVert;
213   }
214   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
215     fCluLayer[i] = tracker.fCluLayer[i];
216     fCluCoord[i] = tracker.fCluCoord[i];
217   } 
218 }
219 //______________________________________________________________________
220 AliITStrackerSA& AliITStrackerSA::operator=(const AliITStrackerSA& source){
221     // Assignment operator. 
222   this->~AliITStrackerSA();
223   new(this) AliITStrackerSA(source);
224   return *this;
225  
226 }
227
228 //____________________________________________________________________________
229 AliITStrackerSA::~AliITStrackerSA(){
230   // destructor
231   // if fVertexer is not null, the AliESDVertex obj. is owned by this class
232   // and is deleted here
233   if(fVertexer){
234     if(fVert)delete fVert;
235   }
236   fVert = 0;
237   fVertexer = 0;
238  
239   if(fPhiWin)delete []fPhiWin;
240   if(fLambdaWin)delete []fLambdaWin;
241   fListOfTracks->Delete();
242   delete fListOfTracks;
243   fListOfSATracks->Delete();
244   delete fListOfSATracks;
245   if(fCluLayer){
246     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
247       if(fCluLayer[i]){
248         fCluLayer[i]->Delete();
249         delete fCluLayer[i];
250       }
251     }
252     delete [] fCluLayer;
253   }
254   if(fCluCoord){
255     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
256       if(fCluCoord[i]){
257         fCluCoord[i]->Delete();
258         delete fCluCoord[i];
259       }
260     }
261     delete [] fCluCoord;
262   }
263   
264 }
265
266 //____________________________________________________________________________
267 Int_t AliITStrackerSA::Clusters2Tracks(AliESDEvent *event){
268 // This method is used to find and fit the tracks. By default the corresponding
269 // method in the parent class is invoked. In this way a combined tracking
270 // TPC+ITS is performed. If the flag fITSStandAlone is true, the tracking
271 // is done in the ITS only. In the standard reconstruction chain this option
272 // can be set via AliReconstruction::SetOption("ITS","onlyITS")
273   Int_t rc=0;
274   if(!fITSStandAlone){
275     rc=AliITStrackerMI::Clusters2Tracks(event);
276   }
277   else {
278     AliDebug(1,"Stand Alone flag set: doing tracking in ITS alone\n");
279   }
280   if(!rc) rc=FindTracks(event);
281   return rc;
282 }
283
284 //____________________________________________________________________________
285 void AliITStrackerSA::Init(){
286   //  Reset all data members
287     fPhiEstimate=0;
288     for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
289     fLambdac=0;
290     fPhic=0;
291     fCoef1=0;
292     fCoef2=0;
293     fCoef3=0;
294     fPointc[0]=0;
295     fPointc[1]=0;
296     fVert = 0;
297     fVertexer = 0;
298     Int_t nLoops=AliITSReconstructor::GetRecoParam()->GetNLoopsSA();
299     if(nLoops==33){
300       SetFixedWindowSizes();
301     }else{
302       Double_t phimin=AliITSReconstructor::GetRecoParam()->GetMinPhiSA();
303       Double_t phimax=AliITSReconstructor::GetRecoParam()->GetMaxPhiSA();
304       Double_t lambmin=AliITSReconstructor::GetRecoParam()->GetMinLambdaSA();
305       Double_t lambmax=AliITSReconstructor::GetRecoParam()->GetMaxLambdaSA();
306       SetCalculatedWindowSizes(nLoops,phimin,phimax,lambmin,lambmax);
307     }
308     fMinQ=AliITSReconstructor::GetRecoParam()->GetSAMinClusterCharge();
309     fITSclusters = 0;
310     SetOuterStartLayer(1);
311     SetSAFlag(kFALSE);
312     fListOfTracks=new TObjArray(0,0);
313     fListOfSATracks=new TObjArray(0,0);
314     fCluLayer = 0;
315     fCluCoord = 0;
316     fMinNPoints = 3;
317  }
318 //_______________________________________________________________________
319 void AliITStrackerSA::ResetForFinding(){
320   //  Reset data members used in all loops during track finding
321     fPhiEstimate=0;
322     for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
323     fLambdac=0;
324     fPhic=0;
325     fCoef1=0;
326     fCoef2=0;
327     fCoef3=0;
328     fPointc[0]=0;
329     fPointc[1]=0;
330     fListOfTracks->Delete();
331     fListOfSATracks->Delete();
332 }
333
334  
335
336 //______________________________________________________________________
337 Int_t AliITStrackerSA::FindTracks(AliESDEvent* event){
338
339 // Track finder using the ESD object
340
341   AliDebug(2,Form(" field is %f",event->GetMagneticField()));
342   
343   if(!fITSclusters){
344     Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
345     return -1;
346   }
347   //Reads event and mark clusters of traks already found, with flag kITSin
348   Int_t nentr=event->GetNumberOfTracks();
349   if(AliITSReconstructor::GetRecoParam()->GetSAUseAllClusters()==kFALSE) {
350     while (nentr--) {
351       AliESDtrack *track=event->GetTrack(nentr);
352       if ((track->GetStatus()&AliESDtrack::kITSin) == AliESDtrack::kITSin){
353         Int_t idx[12];
354         Int_t ncl = track->GetITSclusters(idx);
355         for(Int_t k=0;k<ncl;k++){
356           AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
357           cll->SetBit(kSAflag);
358         }
359       }
360     }
361   }
362   //Get primary vertex
363   Double_t primaryVertex[3];
364   event->GetVertex()->GetXYZ(primaryVertex);
365   //Creates TClonesArray with clusters for each layer. The clusters already used
366   //by AliITStrackerMI are not considered
367   Int_t nclusters[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
368   Int_t dmar[AliITSgeomTGeo::kNLayers]={0,0,0,0,0,0};
369   if (fCluLayer == 0) {
370     fCluLayer = new TClonesArray*[AliITSgeomTGeo::kNLayers];
371     fCluCoord = new TClonesArray*[AliITSgeomTGeo::kNLayers];
372     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
373       fCluLayer[i]=0;
374       fCluCoord[i]=0;
375     }
376   }
377   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
378     AliITSlayer &layer=fgLayers[i];
379     if (!AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) {
380       for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
381         AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
382         if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
383         if(cls->GetQ()==0) continue; //fake clusters dead zones
384         if(i>1 && cls->GetQ()<=fMinQ) continue; // cut on SDD and SSD cluster charge
385         nclusters[i]++;
386       }
387     }
388     dmar[i]=0;
389     delete fCluLayer[i];
390     fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
391     delete fCluCoord[i];
392     fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
393   }
394
395   for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
396     TClonesArray &clulay = *fCluLayer[ilay];
397     TClonesArray &clucoo = *fCluCoord[ilay];
398     AliITSlayer &layer=fgLayers[ilay];
399     if (!AliITSReconstructor::GetRecoParam()->GetLayersToSkip(ilay)) {
400       for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
401         AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
402         if(cls->TestBit(kSAflag)==kTRUE) continue;
403         if(cls->GetQ()==0) continue;
404         if(ilay>1 && cls->GetQ()<=fMinQ) continue; 
405         Double_t phi=0;Double_t lambda=0;
406         Float_t x=0;Float_t y=0;Float_t z=0;
407         Float_t sx=0;Float_t sy=0;Float_t sz=0;
408         GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
409         GetCoorErrors(cls,sx,sy,sz);
410         new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
411         new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
412         dmar[ilay]++;
413       }
414     }
415   }
416    
417   // track counter
418   Int_t ntrack=0;
419
420   static Int_t nClusLay[AliITSgeomTGeo::kNLayers];//counter for clusters on each layer
421
422   if(AliITSReconstructor::GetRecoParam()->GetAllowProlongationWithEmptyRoad()) fMinNPoints=2;
423
424  
425   // loop on minimum number of points
426   for(Int_t iMinNPoints=AliITSgeomTGeo::GetNLayers(); iMinNPoints>=fMinNPoints; iMinNPoints--) {
427
428     if(!fInwardFlag){ // Tracking outwards from the inner layers
429       // loop on starting layer for track finding 
430       for(Int_t innLay=0; innLay<=fOuterStartLayer; innLay++) {
431         if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(innLay)) continue; 
432         Int_t minNPoints=iMinNPoints-innLay;
433         for(Int_t i=innLay+1;i<AliITSgeomTGeo::GetNLayers();i++)
434           if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) 
435             minNPoints--;
436         if(minNPoints<fMinNPoints) continue;
437
438         // loop on phi and lambda window size
439         for(Int_t nloop=0;nloop<fNloop;nloop++){
440           Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
441           while(nclInnLay--){ 
442             ResetForFinding();
443             Bool_t useRP=SetFirstPoint(innLay,nclInnLay,primaryVertex);
444             if(!useRP) continue;            
445             AliITStrackSA* trs = new AliITStrackSA(); 
446             
447             Int_t pflag=0;          
448             Int_t kk;
449             for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
450             
451             kk=0;
452             nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
453                                           trs,primaryVertex[2],pflag);
454             for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
455               kk++;
456               nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
457                                             trs,primaryVertex[2],pflag);
458               if(nClusLay[kk]!=0){
459                 pflag=1;
460                 if(kk==1) {
461                   fPoint3[0]=fPointc[0];
462                   fPoint3[1]=fPointc[1];
463                 } else {
464                   UpdatePoints();
465                 }
466               }
467             }
468             
469             Int_t layOK=0;
470             for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
471               if(nClusLay[nnp]!=0) layOK+=1;
472             }
473             if(layOK>=minNPoints){ 
474               AliDebug(2,Form("---NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]));
475               AliITStrackV2* tr2 = 0;
476               tr2 = FitTrack(trs,primaryVertex);
477               if(!tr2) continue;
478               AliDebug(2,Form("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
479               
480               StoreTrack(tr2,event);
481               ntrack++;
482               
483             }   
484             
485             delete trs;
486           }//end loop on clusters of innLay
487         } //end loop on window sizes
488       } //end loop on innLay
489     }else{// Tracking inwards from the outer layers
490     // loop on starting layer for track finding 
491       for(Int_t outLay=AliITSgeomTGeo::GetNLayers()-1; outLay>=fInnerStartLayer; outLay--) {
492
493         if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(outLay)) continue;
494         Int_t minNPoints=iMinNPoints-(AliITSgeomTGeo::GetNLayers()-1-outLay);
495         for(Int_t i=0;i<outLay-1;i++)
496           if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) 
497             minNPoints--;
498         if(minNPoints<fMinNPoints) continue;
499         
500         // loop on phi and lambda window size
501         for(Int_t nloop=0;nloop<fNloop;nloop++){
502           Int_t nclOutLay=fCluLayer[outLay]->GetEntries();
503           while(nclOutLay--){ 
504             ResetForFinding();
505             Bool_t useRP=SetFirstPoint(outLay,nclOutLay,primaryVertex);
506             if(!useRP) continue;            
507             AliITStrackSA* trs = new AliITStrackSA(); 
508             
509             Int_t pflag=0;          
510             Int_t kk;
511             for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
512             
513             kk=0;
514             nClusLay[kk] = SearchClusters(outLay,fPhiWin[nloop],fLambdaWin[nloop],
515                                   trs,primaryVertex[2],pflag);
516             for(Int_t nextLay=outLay-1; nextLay>=0; nextLay--) {
517               kk++;
518               nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
519                                             trs,primaryVertex[2],pflag);
520               if(nClusLay[kk]!=0){
521                 pflag=1;
522                 if(kk==1) {
523                   fPoint3[0]=fPointc[0];
524                   fPoint3[1]=fPointc[1];
525                 } else {
526                   UpdatePoints();
527                 }
528               }
529             }
530           
531             Int_t layOK=0;
532             for(Int_t nnp=outLay; nnp>=0; nnp--){
533               if(nClusLay[nnp]!=0) layOK+=1;
534             }
535             if(layOK>=minNPoints){ 
536               AliDebug(2,Form("---NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]));
537               AliITStrackV2* tr2 = 0;
538               tr2 = FitTrack(trs,primaryVertex);
539               if(!tr2) continue;
540               AliDebug(2,Form("---NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
541               
542               StoreTrack(tr2,event);
543               ntrack++;
544               
545             }   
546             
547             delete trs;
548         }//end loop on clusters of innLay
549         } //end loop on window sizes
550       } //end loop on innLay
551     } // end if (fInwardFlag)
552   }//end loop on min points
553
554   // search for 1-point tracks in SPD, only for cosmics
555   // (A.Dainese 21.03.08)
556   if(AliITSReconstructor::GetRecoParam()->GetSAOnePointTracks() && 
557      TMath::Abs(event->GetMagneticField())<0.01) {
558     Int_t outerLayer=1; // only SPD
559     for(Int_t innLay=0; innLay<=TMath::Min(1,fOuterStartLayer); innLay++) {
560       //   counter for clusters on each layer  
561
562       for(Int_t nloop=0;nloop<fNloop;nloop++){
563         Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
564         while(nclInnLay--){ //loop starting from layer innLay
565           ResetForFinding();
566           Bool_t useRP=SetFirstPoint(innLay,nclInnLay,primaryVertex);
567           if(!useRP) continue;
568           AliITStrackSA* trs = new AliITStrackSA(); 
569             
570           Int_t pflag=0;            
571           Int_t kk;
572           for(kk=0;kk<AliITSgeomTGeo::GetNLayers();kk++) nClusLay[kk] = 0;
573           
574           kk=0;
575           nClusLay[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
576                                   trs,primaryVertex[2],pflag);
577           for(Int_t nextLay=innLay+1; nextLay<=outerLayer; nextLay++) {
578             kk++;
579             nClusLay[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
580                                     trs,primaryVertex[2],pflag);
581             if(nClusLay[kk]!=0){
582               pflag=1;
583               if(kk==1) {
584                 fPoint3[0]=fPointc[0];
585                 fPoint3[1]=fPointc[1];
586               } else {
587                 UpdatePoints();
588               }
589             }
590           }
591           
592           Int_t layOK=0;
593           for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
594             if(nClusLay[nnp]!=0) layOK+=1;
595           }
596           if(layOK==1) {
597             AliDebug(2,Form("----NPOINTS: %d; MAP: %d %d %d %d %d %d\n",layOK,nClusLay[0],nClusLay[1],nClusLay[2],nClusLay[3],nClusLay[4],nClusLay[5]));
598             AliITStrackV2* tr2 = 0;
599             Bool_t onePoint = kTRUE;
600             tr2 = FitTrack(trs,primaryVertex,onePoint);
601             if(!tr2) continue;
602             AliDebug(2,Form("----NPOINTS fit: %d\n",tr2->GetNumberOfClusters()));
603             
604             StoreTrack(tr2,event);
605             ntrack++;
606             
607           }   
608           
609           delete trs;
610         }//end loop on clusters of innLay
611       } //end loop on window sizes
612       
613     } //end loop on innLay
614   } // end search 1-point tracks
615   
616   Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
617   return 0;
618
619 }
620  
621 //________________________________________________________________________
622
623 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex,Bool_t onePoint) {
624   //fit of the found track (most general case, also <6 points, layers missing)
625   // A.Dainese 16.11.07 
626
627   
628   const Int_t kMaxClu=AliITStrackSA::kMaxNumberOfClusters;
629
630   static Int_t firstmod[AliITSgeomTGeo::kNLayers];
631   
632   static Int_t clind[AliITSgeomTGeo::kNLayers][kMaxClu];
633   static Int_t clmark[AliITSgeomTGeo::kNLayers][kMaxClu];
634   static Int_t nnn[AliITSgeomTGeo::kNLayers];
635   static Int_t kkk[AliITSgeomTGeo::kNLayers];
636   static Int_t end[AliITSgeomTGeo::kNLayers];
637   static AliITSRecPoint *listlayer[AliITSgeomTGeo::kNLayers][kMaxClu];
638
639   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
640     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
641     nnn[i]=0;
642     kkk[i]=0;
643     end[i]=0;
644     for(Int_t j=0;j<kMaxClu; j++){
645       clind[i][j]=0;
646       clmark[i][j]=0;
647       listlayer[i][j]=0;
648    }
649   }
650   
651
652   Int_t nclusters = tr->GetNumberOfClustersSA();
653   for(Int_t ncl=0;ncl<nclusters;ncl++){
654     Int_t index = tr->GetClusterIndexSA(ncl); 
655     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
656     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
657     Int_t lay = (index & 0xf0000000) >> 28;
658     Int_t nInLay=end[lay];
659     listlayer[lay][nInLay]=cl;
660     Int_t ind=nnn[lay];
661     clind[lay][ind]=index;
662     end[lay]++;
663   }
664
665   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
666     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
667       Int_t mark = tr->GetClusterMark(nlay,ncl);
668       Int_t ind=kkk[nlay];
669       clmark[nlay][ind]=mark;
670     }
671   }
672
673
674   Int_t firstLay=-1,secondLay=-1;
675   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
676     if(end[i]==0) {
677       end[i]=1;
678     }else{
679       if(firstLay==-1) {
680         firstLay=i;
681       } else if(secondLay==-1) {
682         secondLay=i;
683       }
684     }
685   }
686
687   if(firstLay==-1 || (secondLay==-1 && !onePoint)) return 0;
688   
689   for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
690     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0][l0];
691     for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
692       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1][l1];
693       for(Int_t l2=0;l2<end[2];l2++){  //loop on layer 3
694         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2][l2];
695         for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4   
696           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3][l3];
697           for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
698             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4][l4];
699             for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6  
700               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5][l5];
701
702
703               Double_t x1,y1,z1,sx1,sy1,sz1;
704               Double_t x2,y2,z2,sx2,sy2,sz2;
705               AliITSRecPoint* p1=0;
706               AliITSRecPoint* p2=0;
707               Int_t index1=0,index2=0;
708               Int_t mrk1=0,mrk2=0;
709
710               switch(firstLay) {
711               case 0:
712                 p1=cl0;
713                 index1=clind[0][l0];mrk1=clmark[0][l0];
714                 break;
715               case 1:
716                 p1=cl1;
717                 index1=clind[1][l1];mrk1=clmark[1][l1];
718                 break;
719               case 2:
720                 p1=cl2;
721                 index1=clind[2][l2];mrk1=clmark[2][l2];
722                 break;
723               case 3:
724                 p1=cl3;
725                 index1=clind[3][l3];mrk1=clmark[3][l3];
726                 break;
727               case 4:
728                 p1=cl4;
729                 index1=clind[4][l4];mrk1=clmark[4][l4];
730                 break;
731               }
732
733               switch(secondLay) {
734               case 1:
735                 p2=cl1;
736                 index2=clind[1][l1];mrk2=clmark[1][l1];
737                 break;
738               case 2:
739                 p2=cl2;
740                 index2=clind[2][l2];mrk2=clmark[2][l2];
741                 break;
742               case 3:
743                 p2=cl3;
744                 index2=clind[3][l3];mrk2=clmark[3][l3];
745                 break;
746               case 4:
747                 p2=cl4;
748                 index2=clind[4][l4];mrk2=clmark[4][l4];
749                 break;
750               case 5:
751                 p2=cl5;
752                 index2=clind[5][l5];mrk2=clmark[5][l5];
753                 break;
754               default:
755                 p2=0;
756                 index2=-1;mrk2=-1;
757                 break;
758               }
759
760               Int_t module1 = p1->GetDetectorIndex()+firstmod[firstLay]; 
761               Int_t layer,ladder,detector;
762               AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
763               Float_t yclu1 = p1->GetY();
764               Float_t zclu1 = p1->GetZ();
765               Double_t cv=0,tgl2=0,phi2=0;
766               
767               Int_t cln1=mrk1;
768               AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(firstLay,cln1);
769               x1 = arr1->GetX();
770               y1 = arr1->GetY();
771               z1 = arr1->GetZ();
772               sx1 = arr1->GetSx();
773               sy1 = arr1->GetSy();
774               sz1 = arr1->GetSz();
775
776               if(secondLay>0) {
777                 Int_t cln2=mrk2;
778                 AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(secondLay,cln2);
779                 x2 = arr2->GetX();
780                 y2 = arr2->GetY();
781                 z2 = arr2->GetZ();
782                 sx2 = arr2->GetSx();
783                 sy2 = arr2->GetSy();
784                 sz2 = arr2->GetSz();
785                 cv = Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
786                 tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
787                 phi2 = TMath::ATan2((y2-y1),(x2-x1));
788               } else { // special case of 1-point tracks, only for cosmics (B=0)
789                 x2 = primaryVertex[0];
790                 y2 = primaryVertex[1];
791                 z2 = primaryVertex[2];
792                 cv = 0;
793                 tgl2 = (z1-z2)/TMath::Sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
794                 phi2 = TMath::ATan2((y1-y2),(x1-x2));
795               }
796
797
798               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
799
800
801               if(cl5!=0) {
802                 trac->AddClusterV2(5,(clind[5][l5] & 0x0fffffff)>>0);
803                 trac->AddClusterMark(5,clmark[5][l5]);
804               }
805               if(cl4!=0){
806                 trac->AddClusterV2(4,(clind[4][l4] & 0x0fffffff)>>0);
807                 trac->AddClusterMark(4,clmark[4][l4]);
808               }
809               if(cl3!=0){
810                 trac->AddClusterV2(3,(clind[3][l3] & 0x0fffffff)>>0);
811                 trac->AddClusterMark(3,clmark[3][l3]);
812               }
813               if(cl2!=0){
814                 trac->AddClusterV2(2,(clind[2][l2] & 0x0fffffff)>>0);
815                 trac->AddClusterMark(2,clmark[2][l2]);
816               }
817               if(cl1!=0){
818                 trac->AddClusterV2(1,(clind[1][l1] & 0x0fffffff)>>0);
819                 trac->AddClusterMark(1,clmark[1][l1]);
820               }
821               if(cl0!=0){
822                 trac->AddClusterV2(0,(clind[0][l0] & 0x0fffffff)>>0);
823                 trac->AddClusterMark(0,clmark[0][l0]);
824               }
825
826               //fit with Kalman filter using AliITStrackerMI::RefitAt()
827               AliITStrackMI* ot = new AliITStrackSA(*trac);
828               
829               ot->ResetCovariance(10.);
830               ot->ResetClusters();
831               
832               // Propagate inside the innermost layer with a cluster 
833               if(ot->Propagate(ot->GetX()-0.1*ot->GetX())) {
834
835                 if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
836                   AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
837                   otrack2->ResetCovariance(10.); 
838                   otrack2->ResetClusters();
839                   //fit from layer 6 to layer 1
840                   if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
841                     fListOfTracks->AddLast(otrack2);
842                     fListOfSATracks->AddLast(trac);
843                   } else {
844                     delete otrack2;
845                     delete trac;
846                   }
847                               
848                 }       
849               }
850               delete ot;
851             }//end loop layer 6
852           }//end loop layer 5
853         }//end loop layer 4        
854       }//end loop layer 3
855     }//end loop layer 2 
856   }//end loop layer 1
857
858
859
860
861   if(fListOfTracks->GetEntries()==0) return 0;
862
863   Int_t lowchi2 = FindTrackLowChiSquare();
864   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
865   AliITStrackSA* trsa = (AliITStrackSA*)fListOfSATracks->At(lowchi2);
866  
867   if(otrack==0) return 0;
868
869   Int_t indexc[AliITSgeomTGeo::kNLayers];
870   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
871   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
872     indexc[nind] = otrack->GetClusterIndex(nind);
873   }      
874   Int_t labl[6][3];
875   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
876     if(i<otrack->GetNumberOfClusters()) {
877       AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
878       labl[i][0]=cl->GetLabel(0);
879       labl[i][1]=cl->GetLabel(1);
880       labl[i][2]=cl->GetLabel(2);
881     } else {
882       labl[i][0]=-1;
883       labl[i][1]=-1;
884       labl[i][2]=-1;
885     }
886   }
887
888   CookLabel(otrack,0.); //MI change - to see fake ratio
889  
890   Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
891   Int_t lflag=0;
892   for(Int_t i=0;i<otrack->GetNumberOfClusters();i++)
893     if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
894   
895   if(lflag<otrack->GetNumberOfClusters()) label = -label;
896   otrack->SetLabel(label);  
897
898   //remove clusters of found track
899   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
900     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
901       Int_t index = trsa->GetClusterMark(nlay,cln);
902       fCluLayer[nlay]->RemoveAt(index);
903       RemoveClusterCoord(nlay,index);
904       fCluLayer[nlay]->Compress();
905     }    
906   }
907
908   return otrack;
909
910 }
911
912 //_______________________________________________________
913 void AliITStrackerSA::StoreTrack(AliITStrackV2 *t,AliESDEvent *event) const 
914 {
915   //
916   // Add new track to the ESD
917   //
918   AliESDtrack outtrack;
919   outtrack.UpdateTrackParams(t,AliESDtrack::kITSin);
920   for(Int_t i=0;i<12;i++) {
921     outtrack.SetITSModuleIndex(i,t->GetModuleIndex(i));
922   }
923   Double_t sdedx[4]={0.,0.,0.,0.};
924   for(Int_t i=0; i<4; i++) sdedx[i]=t->GetSampledEdx(i);
925   outtrack.SetITSdEdxSamples(sdedx);
926   event->AddTrack(&outtrack);
927
928   return;
929 }
930
931
932 //_______________________________________________________
933 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
934   //function used to to find the clusters associated to the track
935
936   if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(layer)) return 0;
937
938   Int_t nc=0;
939   AliITSlayer &lay = fgLayers[layer];
940   Double_t r=lay.GetR();
941   if(pflag==1){      
942     Float_t cx1,cx2,cy1,cy2;
943     FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
944     if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
945        return 0;
946     Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
947     Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
948     fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
949   }
950
951  
952   Int_t ncl = fCluLayer[layer]->GetEntries();
953   for (Int_t index=0; index<ncl; index++) {
954     AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
955     if (!c) continue;
956     if (c->GetQ()<=0) continue;
957     if(layer>1 && c->GetQ()<=fMinQ) continue;
958     
959      AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
960      Double_t phi = arr->GetPhi();
961      if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
962
963      Double_t lambda = arr->GetLambda();
964      if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
965
966      if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
967      if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
968      Int_t orind = arr->GetOrInd();
969      trs->AddClusterSA(layer,orind);
970      trs->AddClusterMark(layer,index);
971        
972      nc++;
973      fLambdac=lambda;
974      fPhiEstimate=phi;
975
976      fPointc[0]=arr->GetX();
977      fPointc[1]=arr->GetY();
978
979   }
980   return nc;
981 }
982
983 //________________________________________________________________
984 Bool_t AliITStrackerSA::SetFirstPoint(Int_t lay, Int_t clu, Double_t* primaryVertex){
985   // Sets the first point (seed) for tracking
986
987   AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[lay]->At(clu);
988   if(!cl) return kFALSE;
989   if (cl->GetQ()<=0) return kFALSE;
990   if(lay>1 && cl->GetQ()<=fMinQ) return kFALSE;
991
992   AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(lay,clu);
993   fPhic = arr->GetPhi();
994   fLambdac = arr->GetLambda();
995   fPhiEstimate = fPhic;
996   fPoint1[0]=primaryVertex[0];
997   fPoint1[1]=primaryVertex[1];
998   fPoint2[0]=arr->GetX();
999   fPoint2[1]=arr->GetY();
1000   return kTRUE; 
1001 }
1002
1003 //________________________________________________________________
1004 void AliITStrackerSA::UpdatePoints(){
1005   //update of points for the estimation of the curvature  
1006
1007   fPoint2[0]=fPoint3[0];
1008   fPoint2[1]=fPoint3[1];
1009   fPoint3[0]=fPointc[0];
1010   fPoint3[1]=fPointc[1];
1011
1012   
1013 }
1014
1015 //___________________________________________________________________
1016 Int_t AliITStrackerSA::FindEquation(Float_t x1, Float_t y1, Float_t x2, Float_t y2, Float_t x3, Float_t y3,Float_t& a, Float_t& b, Float_t& c){
1017
1018    //given (x,y) of three recpoints (in global coordinates) 
1019    //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1020
1021    Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1022    if(den==0) return 0;
1023    a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1024    b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1025    c = -x1*x1-y1*y1-a*x1-b*y1;
1026    return 1;
1027  }
1028 //__________________________________________________________________________
1029  Int_t AliITStrackerSA::FindIntersection(Float_t a1, Float_t b1, Float_t c1, Float_t c2,Float_t& x1,Float_t& y1, Float_t& x2, Float_t& y2){
1030  
1031  //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1032  //c2 is -rlayer*rlayer
1033
1034   if(a1==0) return 0;
1035  Double_t m = c2-c1; 
1036  Double_t aA = (b1*b1)/(a1*a1)+1;
1037  Double_t bB = (-2*m*b1/(a1*a1));
1038  Double_t cC = c2+(m*m)/(a1*a1);
1039  Double_t dD = bB*bB-4*aA*cC;
1040  if(dD<0) return 0;
1041  
1042  y1 = (-bB+TMath::Sqrt(dD))/(2*aA); 
1043  y2 = (-bB-TMath::Sqrt(dD))/(2*aA); 
1044  x1 = (c2-c1-b1*y1)/a1;
1045  x2 = (c2-c1-b1*y2)/a1;
1046
1047  return 1; 
1048 }
1049 //____________________________________________________________________
1050 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t 
1051 x2,Double_t y2,Double_t x3,Double_t y3){
1052
1053   //calculates the curvature of track  
1054   Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1055   if(den==0) return 0;
1056   Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1057   Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1058   Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1059   Double_t xc=-a/2.;
1060
1061   if((a*a+b*b-4*c)<0) return 0;
1062   Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1063   if(rad==0) return 0;
1064   
1065   if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1066   if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1067   //  if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1068   // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1069   
1070   return 1/rad;
1071  
1072 }
1073
1074
1075 //____________________________________________________________________
1076 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1077
1078   //Returns the point closest to pp
1079
1080   Double_t diff1 = p1-pp;
1081   Double_t diff2 = p2-pp;
1082   
1083   if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1084   else fPhiEstimate=p2;  
1085   return fPhiEstimate;
1086   
1087 }
1088
1089
1090 //_________________________________________________________________
1091 Int_t AliITStrackerSA::FindTrackLowChiSquare() const {
1092   // returns track with lowest chi square  
1093   Int_t dim=fListOfTracks->GetEntries();
1094   if(dim<=1) return 0;
1095   AliITStrackV2* trk = (AliITStrackV2*)fListOfTracks->At(0);
1096   Double_t minChi2=trk->GetChi2();
1097   Int_t index=0;
1098   for(Int_t i=1;i<dim;i++){
1099     trk = (AliITStrackV2*)fListOfTracks->At(i);
1100     Double_t chi2=trk->GetChi2();
1101     if(chi2<minChi2){
1102       minChi2=chi2;
1103       index=i;
1104     }
1105   }
1106   return index;
1107 }
1108
1109 //__________________________________________________________
1110 Int_t AliITStrackerSA::FindLabel(Int_t l0, Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5){
1111
1112   //function used to determine the track label
1113   
1114   Int_t lb[6] = {l0,l1,l2,l3,l4,l5};
1115   Int_t aa[6]={1,1,1,1,1,1};
1116   Int_t ff=0; 
1117   Int_t ll=0;
1118   Int_t k=0;Int_t w=0;Int_t num=6;
1119   for(Int_t i=5;i>=0;i--) if(lb[i]==-1) num=i;
1120   
1121   while(k<num){
1122   
1123     for(Int_t i=k+1;i<num;i++){
1124     
1125       if(lb[k]==lb[i] && aa[k]!=0){
1126       
1127         aa[k]+=1;
1128         aa[i]=0;
1129       }
1130     }
1131     k++;
1132   }
1133
1134   while(w<num){
1135  
1136     for(Int_t j=0;j<6;j++){
1137       if(aa[w]<aa[j]){
1138       ff=aa[w];
1139       aa[w]=aa[j];
1140       aa[j]=ff;
1141       ll=lb[w];
1142       lb[w]=lb[j];
1143       lb[j]=ll;
1144      }
1145     }
1146     w++;
1147   }
1148   
1149   if(num<1) return -1; 
1150   return lb[num-1];
1151 }
1152
1153 //_____________________________________________________________________________
1154 Int_t AliITStrackerSA::Label(Int_t gl1, Int_t gl2, Int_t gl3, Int_t gl4, Int_t gl5, Int_t gl6,Int_t gl7, Int_t gl8, Int_t gl9, Int_t gl10,Int_t gl11,
1155 Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t minNPoints){
1156
1157  
1158   //function used to assign label to the found track. If track is fake, the label is negative
1159
1160   Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1161   Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1162   Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1163   Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1164   Int_t lflag=0;Int_t num=6;
1165   if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1166
1167   for(Int_t i=0;i<num;i++){
1168     if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1169   }
1170
1171   if(lflag>=minNPoints) return ll;
1172   else return -ll;
1173
1174   
1175 }
1176 //_____________________________________________________________________________
1177 void AliITStrackerSA::SetCalculatedWindowSizes(Int_t n, Float_t phimin, Float_t phimax, Float_t lambdamin, Float_t lambdamax){
1178   // Set sizes of the phi and lambda windows used for track finding
1179   fNloop = n;
1180   if(fPhiWin) delete [] fPhiWin;
1181   if(fLambdaWin) delete [] fLambdaWin;
1182   fPhiWin = new Double_t[fNloop];
1183   fLambdaWin = new Double_t[fNloop];
1184   Float_t stepPhi=(phimax-phimin)/(Float_t)(fNloop-1);
1185   Float_t stepLambda=(lambdamax-lambdamin)/(Float_t)(fNloop-1);
1186   for(Int_t k=0;k<fNloop;k++){
1187     Float_t phi=phimin+k*stepPhi;
1188     Float_t lam=lambdamin+k*stepLambda;
1189     fPhiWin[k]=phi;
1190     fLambdaWin[k]=lam;
1191   }
1192 }
1193 //_____________________________________________________________________________
1194 void AliITStrackerSA::SetFixedWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1195   // Set sizes of the phi and lambda windows used for track finding
1196   fNloop = n;
1197   if(phi){ // user defined values
1198     fPhiWin = new Double_t[fNloop];
1199     fLambdaWin = new Double_t[fNloop];
1200     for(Int_t k=0;k<fNloop;k++){
1201       fPhiWin[k]=phi[k];
1202       fLambdaWin[k]=lam[k];
1203     }
1204   }
1205   else {  // default values
1206             
1207     Double_t phid[33]   = {0.002,0.003,0.004,0.0045,0.0047,
1208                            0.005,0.0053,0.0055,
1209                            0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1210                            0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1211                            0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1212     Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1213                             0.005,0.005,0.006,
1214                             0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1215                             0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1216                             0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1217     
1218     if(fNloop!=33){
1219       fNloop = 33;
1220     }
1221     
1222     
1223     fPhiWin = new Double_t[fNloop];
1224     fLambdaWin = new Double_t[fNloop];
1225
1226     Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
1227   
1228     for(Int_t k=0;k<fNloop;k++){
1229       fPhiWin[k]=phid[k]*factor;
1230       fLambdaWin[k]=lambdad[k]*factor;
1231     }
1232   
1233   }
1234
1235 }
1236 //_______________________________________________________________________
1237 void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
1238   //Returns values of phi (azimuthal) and lambda angles for a given cluster
1239 /*  
1240   Double_t rot[9];     fGeom->GetRotMatrix(module,rot);
1241   Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1242   Float_t tx,ty,tz;  fGeom->GetTrans(lay,lad,det,tx,ty,tz);     
1243
1244   Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1245   Double_t phi1=TMath::Pi()/2+alpha;
1246   if (lay==1) phi1+=TMath::Pi();
1247
1248   Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1249   Float_t r=tx*cp+ty*sp;
1250
1251   xyz= r*cp - cl->GetY()*sp;
1252   y= r*sp + cl->GetY()*cp;
1253   z=cl->GetZ();
1254 */
1255   Float_t xyz[3];
1256   cl->GetGlobalXYZ(xyz);
1257   x=xyz[0];
1258   y=xyz[1];
1259   z=xyz[2];
1260  
1261   phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
1262   lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1263 }
1264
1265 //________________________________________________________________________
1266 void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
1267
1268   //returns sigmax, y, z of cluster in global coordinates
1269 /*
1270   Double_t rot[9];     fGeom->GetRotMatrix(module,rot);
1271   Int_t lay,lad,det; 
1272   AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
1273  
1274   Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1275   Double_t phi=TMath::Pi()/2+alpha;
1276   if (lay==1) phi+=TMath::Pi();
1277
1278   Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
1279 */
1280   Float_t covm[6];
1281   cl->GetGlobalCov(covm);
1282   sx=TMath::Sqrt(covm[0]);
1283   sy=TMath::Sqrt(covm[3]);
1284   sz=TMath::Sqrt(covm[5]);
1285 /*
1286   sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1287   sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1288   sz = TMath::Sqrt(cl->GetSigmaZ2());
1289 */
1290 }
1291