Possibility of skipping layers in trackers MI and SA (A. Dainese)
[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 "AliRun.h"
43
44 ClassImp(AliITStrackerSA)
45
46 //____________________________________________________________________________
47 AliITStrackerSA::AliITStrackerSA():AliITStrackerMI(),
48 fPhiEstimate(0),
49 fITSStandAlone(0),
50 fLambdac(0),
51 fPhic(0),
52 fCoef1(0),
53 fCoef2(0),
54 fCoef3(0),
55 fNloop(0),
56 fPhiWin(0),
57 fLambdaWin(0),
58 fVert(0),
59 fVertexer(0),
60 fListOfTracks(0),
61 fITSclusters(0),
62 fSixPoints(0),
63 fOuterStartLayer(0),
64 fCluLayer(0),
65 fCluCoord(0){
66   // Default constructor
67   Init();
68  
69 }
70 //____________________________________________________________________________
71 AliITStrackerSA::AliITStrackerSA(const Char_t *geom):AliITStrackerMI(0),
72 fPhiEstimate(0),
73 fITSStandAlone(0),
74 fLambdac(0),
75 fPhic(0),
76 fCoef1(0),
77 fCoef2(0),
78 fCoef3(0),
79 fNloop(0),
80 fPhiWin(0),
81 fLambdaWin(0),
82 fVert(0),
83 fVertexer(0),
84 fListOfTracks(0),
85 fITSclusters(0),
86 fSixPoints(0),
87 fOuterStartLayer(0),
88 fCluLayer(0),
89 fCluCoord(0) 
90 {
91   // Standard constructor (Vertex is known and passed to this obj.)
92   if (geom) {
93     AliWarning("\"geom\" is actually a dummy argument !");
94   }
95
96   Init();
97   fVert = 0;
98  
99 }
100
101 //____________________________________________________________________________
102 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliESDVertex *vert):AliITStrackerMI(0),
103 fPhiEstimate(0),
104 fITSStandAlone(0),
105 fLambdac(0),
106 fPhic(0),
107 fCoef1(0),
108 fCoef2(0),
109 fCoef3(0),
110 fNloop(0),
111 fPhiWin(0),
112 fLambdaWin(0),
113 fVert(vert),
114 fVertexer(0),
115 fListOfTracks(0),
116 fITSclusters(0),
117 fSixPoints(0),
118 fOuterStartLayer(0),
119 fCluLayer(0),
120 fCluCoord(0)
121 {
122   // Standard constructor (Vertex is known and passed to this obj.)
123   if (geom) {
124     AliWarning("\"geom\" is actually a dummy argument !");
125   }
126   Init();
127  
128 }
129
130 //____________________________________________________________________________
131 AliITStrackerSA::AliITStrackerSA(const Char_t *geom, AliITSVertexer *vertexer):AliITStrackerMI(0),
132 fPhiEstimate(0),
133 fITSStandAlone(0),
134 fLambdac(0),
135 fPhic(0),
136 fCoef1(0),
137 fCoef2(0),
138 fCoef3(0),
139 fNloop(0),
140 fPhiWin(0),
141 fLambdaWin(0),
142 fVert(),
143 fVertexer(vertexer),
144 fListOfTracks(0),
145 fITSclusters(0),
146 fSixPoints(0),
147 fOuterStartLayer(0),
148 fCluLayer(0),
149 fCluCoord(0)
150 {
151   // Standard constructor (Vertex is unknown - vertexer is passed to this obj)
152   if (geom) {
153     AliWarning("\"geom\" is actually a dummy argument !");
154   }
155   Init();
156   fVertexer = vertexer;
157  
158 }
159
160 //____________________________________________________________________________
161 AliITStrackerSA::AliITStrackerSA(const AliITStrackerSA& tracker):AliITStrackerMI(),
162 fPhiEstimate(tracker.fPhiEstimate),
163 fITSStandAlone(tracker.fITSStandAlone),
164 fLambdac(tracker.fLambdac),
165 fPhic(tracker.fPhic),
166 fCoef1(tracker.fCoef1),
167 fCoef2(tracker.fCoef2),
168 fCoef3(tracker.fCoef3),
169 fNloop(tracker.fNloop),
170 fPhiWin(tracker.fPhiWin),
171 fLambdaWin(tracker.fLambdaWin),
172 fVert(tracker.fVert),
173 fVertexer(tracker.fVertexer),
174 fListOfTracks(tracker.fListOfTracks),
175 fITSclusters(tracker.fITSclusters),
176 fSixPoints(tracker.fSixPoints),
177 fOuterStartLayer(tracker.fOuterStartLayer),
178 fCluLayer(tracker.fCluLayer),
179 fCluCoord(tracker.fCluCoord) {
180   // Copy constructor
181   for(Int_t i=0;i<2;i++){
182     fPoint1[i]=tracker.fPoint1[i];
183     fPoint2[i]=tracker.fPoint2[i];
184     fPoint3[i]=tracker.fPoint3[i];
185     fPointc[i]=tracker.fPointc[i];
186   }
187   if(tracker.fVertexer && tracker.fVert){
188     fVert = new AliESDVertex(*tracker.fVert);
189   }
190   else {
191     fVert = tracker.fVert;
192   }
193   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
194     fCluLayer[i] = tracker.fCluLayer[i];
195     fCluCoord[i] = tracker.fCluCoord[i];
196   } 
197 }
198 //______________________________________________________________________
199 AliITStrackerSA& AliITStrackerSA::operator=(const AliITStrackerSA& source){
200     // Assignment operator. 
201   this->~AliITStrackerSA();
202   new(this) AliITStrackerSA(source);
203   return *this;
204  
205 }
206
207 //____________________________________________________________________________
208 AliITStrackerSA::~AliITStrackerSA(){
209   // destructor
210   // if fVertexer is not null, the AliESDVertex obj. is owned by this class
211   // and is deleted here
212   if(fVertexer){
213     if(fVert)delete fVert;
214   }
215   fVert = 0;
216   fVertexer = 0;
217  
218   if(fPhiWin)delete []fPhiWin;
219   if(fLambdaWin)delete []fLambdaWin;
220   fListOfTracks->Delete();
221   delete fListOfTracks;
222   if(fCluLayer){
223     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
224       if(fCluLayer[i]){
225         fCluLayer[i]->Delete();
226         delete fCluLayer[i];
227       }
228     }
229     delete [] fCluLayer;
230   }
231   if(fCluCoord){
232     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
233       if(fCluCoord[i]){
234         fCluCoord[i]->Delete();
235         delete fCluCoord[i];
236       }
237     }
238     delete [] fCluCoord;
239   }
240   
241 }
242
243 //____________________________________________________________________________
244 Int_t AliITStrackerSA::Clusters2Tracks(AliESDEvent *event){
245 // This method is used to find and fit the tracks. By default the corresponding
246 // method in the parent class is invoked. In this way a combined tracking
247 // TPC+ITS is performed. If the flag fITSStandAlone is true, the tracking
248 // is done in the ITS only. In the standard reconstruction chain this option
249 // can be set via AliReconstruction::SetOption("ITS","onlyITS")
250   Int_t rc=0;
251   if(!fITSStandAlone){
252     rc=AliITStrackerMI::Clusters2Tracks(event);
253   }
254   else {
255     AliDebug(1,"Stand Alone flag set: doing tracking in ITS alone\n");
256   }
257   if(!rc) rc=FindTracks(event);
258   return rc;
259 }
260
261 //____________________________________________________________________________
262 void AliITStrackerSA::Init(){
263   //  Reset all data members
264     fPhiEstimate=0;
265     for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
266     fLambdac=0;
267     fPhic=0;
268     fCoef1=0;
269     fCoef2=0;
270     fCoef3=0;
271     fPointc[0]=0;
272     fPointc[1]=0;
273     fVert = 0;
274     fVertexer = 0;
275     SetWindowSizes();
276     fITSclusters = 0;
277     SetSixPoints();
278     SetOuterStartLayer(0);
279     SetSAFlag(kFALSE);
280     fListOfTracks=new TObjArray(0,0);
281     fCluLayer = 0;
282     fCluCoord = 0;
283  }
284 //_______________________________________________________________________
285 void AliITStrackerSA::ResetForFinding(){
286   //  Reset data members used in all loops during track finding
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     fListOfTracks->Delete();
297 }
298
299  
300
301 //______________________________________________________________________
302 Int_t AliITStrackerSA::FindTracks(AliESDEvent* event){
303
304 // Track finder using the ESD object
305
306
307   //controllare numero cluster sui layer1 e 2 (morti?)
308   //non trova tracce...controllare..
309
310   if(!fITSclusters){
311     Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
312     return -1;
313   }
314   
315    
316   //Reads event and mark clusters of traks already found, with flag kITSin
317    Int_t nentr=event->GetNumberOfTracks();
318    while (nentr--) {
319      AliESDtrack *track=event->GetTrack(nentr);
320      if (track->GetStatus()&AliESDtrack::kITSin==AliESDtrack::kITSin){
321        Int_t idx[12];
322        Int_t ncl = track->GetITSclusters(idx);
323        for(Int_t k=0;k<ncl;k++){
324          AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
325          cll->SetBit(kSAflag);
326        }
327      }
328    }
329
330    //Get primary vertex
331    Double_t primaryVertex[3];
332    event->GetVertex()->GetXYZ(primaryVertex);
333    //Creates TClonesArray with clusters for each layer. The clusters already used
334    //by AliITStrackerMI are not considered
335    
336    Int_t nclusters[6]={0,0,0,0,0,0};
337    Int_t dmar[6]={0,0,0,0,0,0};
338    if (fCluLayer == 0) {
339    fCluLayer = new TClonesArray*[AliITSgeomTGeo::GetNLayers()];
340    fCluCoord = new TClonesArray*[AliITSgeomTGeo::GetNLayers()];
341      for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
342         fCluLayer[i]=0;
343         fCluCoord[i]=0;
344      }
345    }
346    for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
347      AliITSlayer &layer=fgLayers[i];
348      for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
349        AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
350        if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
351        if(cls->GetQ()==0) continue; //fake clusters dead zones
352        nclusters[i]++;
353      }
354      dmar[i]=0;
355      delete fCluLayer[i];
356      fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
357      delete fCluCoord[i];
358      fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
359    }
360
361    
362    Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
363    for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
364      firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
365    }
366    
367    for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
368      TClonesArray &clulay = *fCluLayer[ilay];
369      TClonesArray &clucoo = *fCluCoord[ilay];
370      AliITSlayer &layer=fgLayers[ilay];
371      for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
372        AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
373        if(cls->TestBit(kSAflag)==kTRUE) continue;
374        if(cls->GetQ()==0) continue;
375        Double_t phi=0;Double_t lambda=0;
376        Float_t x=0;Float_t y=0;Float_t z=0;
377        Float_t sx=0;Float_t sy=0;Float_t sz=0;
378        GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
379        GetCoorErrors(cls,sx,sy,sz);
380        new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
381        new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
382        dmar[ilay]++;
383      }
384    }
385    
386  
387    Int_t minNPoints = (fSixPoints ? AliITSgeomTGeo::GetNLayers() : AliITSgeomTGeo::GetNLayers()-1);  
388    Bool_t jumpy = kFALSE;
389    for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
390      if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) {
391        minNPoints--;
392        jumpy = kTRUE;
393      }
394    }
395
396    Int_t ntrack=0;
397    //loop on the different windows
398    for(Int_t nloop=0;nloop<fNloop;nloop++){
399      for(Int_t ncl=0;ncl<fCluLayer[0]->GetEntries();ncl++){ //loop starting from layer 0
400        
401        ResetForFinding();
402        Int_t pflag=0;
403        
404        AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[0]->At(ncl);
405
406        if(!cl) continue;
407
408
409        if (cl->GetQ()<=0) continue;
410        
411        AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(0,ncl); 
412        fPhic = arr->GetPhi();
413        fLambdac = arr->GetLambda();
414        if (TMath::Abs(fLambdac)>0.26*TMath::Pi()) continue;
415        fPhiEstimate = fPhic;
416        AliITStrackSA* trs = new AliITStrackSA(); 
417        fPoint1[0]=primaryVertex[0];
418        fPoint1[1]=primaryVertex[1];
419
420
421        fPoint2[0]=arr->GetX();
422        fPoint2[1]=arr->GetY();
423        Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()];//counter for clusters on each layer
424        for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){ nn[i]=0;}
425        nn[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
426        
427        nn[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
428       if(nn[1]>0){
429         pflag=1;
430         fPoint3[0] = fPointc[0];
431         fPoint3[1] = fPointc[1];
432       }
433       nn[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
434       if(nn[1]==0 && nn[2]==0) pflag=0;
435       if(nn[2]!=0 && nn[1]!=0){ pflag=1; UpdatePoints();}
436       if(nn[2]!=0 && nn[1]==0){
437         pflag=1;
438         fPoint3[0]=fPointc[0];
439         fPoint3[1]=fPointc[1];
440       }
441
442       nn[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
443       pflag=1;
444       if(nn[3]!=0) UpdatePoints();
445       nn[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag); 
446       pflag=1;
447       if(nn[4]!=0) UpdatePoints();
448       nn[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag); 
449           
450
451       Int_t layOK=0;
452       //check of the candidate track
453       for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers();nnp++) {
454         if(nn[nnp]!=0) layOK+=1;
455       }
456
457       if(layOK>=minNPoints){
458         AliITStrackV2* tr2 = 0;
459         if(!jumpy) {
460           tr2 = FitTrack(trs,primaryVertex);
461         } else {
462           tr2 = FitJumpyTrack(trs,primaryVertex);
463         }
464         if(!tr2) continue;
465
466         AliESDtrack outtrack;
467         outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
468         event->AddTrack(&outtrack);
469         ntrack++;
470       }
471         
472
473       delete trs;
474       delete[] nn;
475       
476      }//end loop on clusters of layer1
477     
478      //end loop2
479    }
480
481    minNPoints--;
482
483    //if 5/6 points are required, second loop starting 
484    //from second layer (SPD2), to find tracks with point of 
485    //layer 1 missing   
486    if(!fSixPoints) {
487      //printf("looking from SPD2\n");
488      //   counter for clusters on each layer  
489      Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()-1];      
490      for(Int_t nloop=0;nloop<fNloop;nloop++){
491        Int_t ncl2=fCluLayer[1]->GetEntries();
492        while(ncl2--){ //loop starting from layer 2
493          ResetForFinding();
494          Int_t pflag=0;
495          AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[1]->At(ncl2);
496          
497          if(!cl) continue;
498          AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(1,ncl2);
499          fPhic = arr->GetPhi();
500          fLambdac = arr->GetLambda();
501          fPhiEstimate = fPhic;
502          
503          AliITStrackSA* trs = new AliITStrackSA(); 
504          fPoint1[0]=primaryVertex[0];
505          fPoint1[1]=primaryVertex[1];
506          
507          fPoint2[0]=arr->GetX();
508          fPoint2[1]=arr->GetY();
509          for(Int_t kk=0;kk<AliITSgeomTGeo::GetNLayers()-1;kk++)nn[kk] = 0;
510          nn[0] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],
511                                    trs,primaryVertex[2],pflag);
512          nn[1] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],
513                                    trs,primaryVertex[2],pflag);
514          if(nn[1]!=0){
515            pflag=1;
516            fPoint3[0]=fPointc[0];
517            fPoint3[1]=fPointc[1];
518          }
519          nn[2]= SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],
520                                    trs,primaryVertex[2],pflag);
521          if(nn[2]!=0){
522            pflag=1;
523            UpdatePoints();
524          }
525          nn[3]= SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],
526                                    trs,primaryVertex[2],pflag);
527          if(nn[3]!=0){
528            pflag=1;
529            UpdatePoints();
530          }
531          nn[4]=SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],
532                                    trs,primaryVertex[2],pflag);
533
534          Int_t fl=0;
535          for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-1;nnp++){
536            if(nn[nnp]!=0) fl+=1;
537          }
538          if(fl>=minNPoints){  // 5/6       
539            AliITStrackV2* tr2 = 0;
540            if(!jumpy) {
541              tr2 = FitTrack(trs,primaryVertex);
542            } else {
543              tr2 = FitJumpyTrack(trs,primaryVertex);
544            }
545            if(!tr2) continue;
546            
547            AliESDtrack outtrack;
548            outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
549            event->AddTrack(&outtrack);
550            ntrack++;
551
552          }   
553               
554         delete trs;
555       }//end loop on clusters of layer2
556      }
557    
558      delete [] nn;
559    } //end opt="5/6"  
560
561
562
563    // search for tracks starting from SPD2, SDD1, SDD2, SSD2
564    // for cosmics (A. Dainese 31.07.07)
565    if(fOuterStartLayer>0) {
566      for(Int_t innLay=1; innLay<=fOuterStartLayer; innLay++) {
567        //printf("Searching from layer %d outward\n",innLay);
568        minNPoints=AliITSgeomTGeo::GetNLayers()-innLay;
569          for(Int_t i=innLay;i<AliITSgeomTGeo::GetNLayers();i++) 
570            if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(i)) 
571              minNPoints--;
572        //   counter for clusters on each layer  
573        Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()-innLay];      
574        for(Int_t nloop=0;nloop<fNloop;nloop++){
575          Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
576          while(nclInnLay--){ //loop starting from layer innLay
577            ResetForFinding();
578            Int_t pflag=0;
579            AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
580          
581            if(!cl) continue;
582            AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
583            fPhic = arr->GetPhi();
584            fLambdac = arr->GetLambda();
585            fPhiEstimate = fPhic;
586          
587            AliITStrackSA* trs = new AliITStrackSA(); 
588            fPoint1[0]=primaryVertex[0];
589            fPoint1[1]=primaryVertex[1];
590            fPoint2[0]=arr->GetX();
591            fPoint2[1]=arr->GetY();
592
593            Int_t kk;
594            for(kk=0;kk<AliITSgeomTGeo::GetNLayers()-innLay;kk++) nn[kk] = 0;
595
596            kk=0;
597            nn[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
598                                   trs,primaryVertex[2],pflag);
599            for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
600              kk++;
601              nn[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
602                                      trs,primaryVertex[2],pflag);
603              if(nn[kk]!=0){
604                pflag=1;
605                if(kk==1) {
606                  fPoint3[0]=fPointc[0];
607                  fPoint3[1]=fPointc[1];
608                } else {
609                  UpdatePoints();
610                }
611              }
612            }
613
614            Int_t fl=0;
615            for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
616              if(nn[nnp]!=0) fl+=1;
617            }
618            if(fl>=minNPoints){ 
619              AliITStrackV2* tr2 = 0;
620              if(!jumpy) {
621                tr2 = FitShortTrack(trs,primaryVertex,innLay);
622              } else {
623                tr2 = FitJumpyTrack(trs,primaryVertex);
624              }
625              if(!tr2) continue;
626              
627              
628              AliESDtrack outtrack;
629              outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
630              event->AddTrack(&outtrack);
631              ntrack++;
632              
633            }   
634               
635            delete trs;
636          }//end loop on clusters of innLay
637        } //end loop on window sizes
638    
639        delete [] nn;
640      } //end loop on innLay
641    } //end if(fOuterStartLayer>0)
642
643    
644   delete [] firstmod;
645   Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
646   return 0;
647
648 }
649  
650
651 //________________________________________________________________________
652
653 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex){
654   //fit of the found track
655
656   
657   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
658   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
659     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
660   }  
661
662   Int_t nclusters = tr->GetNumberOfClustersSA();
663   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
664   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
665     listlayer[i] = new TObjArray(0,0);
666   }
667
668   TArrayI clind0(20);
669   TArrayI clind1(20);
670   TArrayI clind2(20);
671   TArrayI clind3(20);
672   TArrayI clind4(20);
673   TArrayI clind5(20);
674
675   TArrayI mark0(20);
676   TArrayI mark1(20);
677   TArrayI mark2(20);
678   TArrayI mark3(20);
679   TArrayI mark4(20);
680   TArrayI mark5(20);
681
682
683   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
684   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
685   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
686   
687   for(Int_t ncl=0;ncl<nclusters;ncl++){
688     Int_t index = tr->GetClusterIndexSA(ncl); 
689     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
690     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
691     Int_t lay = (index & 0xf0000000) >> 28;
692     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
693     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
694     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
695     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
696     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
697     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}    
698   }
699   delete [] nnn;
700
701   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
702     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
703       Int_t mark = tr->GetClusterMark(nlay,ncl);
704       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
705       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
706       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
707       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
708       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
709       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
710
711     }
712   }
713
714   delete [] kkk;
715
716  
717   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
718   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
719     if(listlayer[i]->GetEntries()==0) end[i]=1;
720     else end[i]=listlayer[i]->GetEntries();
721   }
722
723   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
724   TClonesArray &tri = *listSA;
725   Int_t nlist=0;
726
727   if(end[0]==0) end[0]=1; //for tracks with cluster on layer 0 missing
728   for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
729     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l0); 
730     Double_t x1,y1,z1,sx1,sy1,sz1;
731     Double_t x2,y2,z2,sx2,sy2,sz2;
732     AliITSRecPoint* p1=0;
733     AliITSRecPoint* p2=0;
734     Int_t index1=clind0[l0];
735     Int_t index2=0;
736     Int_t mrk1 = mark0[l0];
737     Int_t mrk2 = 0;
738     Int_t lay1=0;
739     Int_t lay2=1;
740     Int_t module1=-1;
741     for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
742       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l1); 
743       index2=clind1[l1];
744       mrk2 = mark1[l1];
745       for(Int_t l2=0;l2<end[2];l2++){  //loop on layer 3
746         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l2);
747
748         if(cl0==0 && cl1!=0) {
749           p1=cl1;
750           p2=cl2;
751           index1=clind1[l1];mrk1=mark1[l1];
752           index2=clind2[l2];mrk2=mark2[l2];
753           lay1=1;
754           lay2=2;
755           module1 = p1->GetDetectorIndex()+firstmod[1]; 
756         }
757         if(cl0!=0 && cl1==0){
758           p1=cl0;
759           p2=cl2;index2=clind2[l2];mrk2=mark2[l2];
760           lay1=0;
761           lay2=2;
762           module1 = p1->GetDetectorIndex()+firstmod[0];
763         }
764         if(cl0!=0 && cl1!=0){
765           p1=cl0;
766           p2=cl1;
767           lay1=0;
768           lay2=1;
769           module1 = p1->GetDetectorIndex()+firstmod[0];
770         }
771         
772         Int_t cln1=mrk1;
773         Int_t cln2=mrk2;
774         AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(lay1,cln1);
775         AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(lay2,cln2);
776         x1 = arr1->GetX();
777         x2 = arr2->GetX();
778         y1 = arr1->GetY();
779         y2 = arr2->GetY();
780         z1 = arr1->GetZ();
781         z2 = arr2->GetZ();
782         sx1 = arr1->GetSx();
783         sx2 = arr2->GetSx();
784         sy1 = arr1->GetSy();
785         sy2 = arr2->GetSy();
786         sz1 = arr1->GetSz();
787         sz2 = arr2->GetSz();
788         
789         Int_t layer,ladder,detector;
790         AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
791         Float_t yclu1 = p1->GetY();
792         Float_t zclu1 = p1->GetZ();
793         Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);        
794         Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
795         Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
796
797         for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4   
798           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l3);
799           for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
800             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l4);
801             for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6  
802               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l5);
803               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
804                               
805               if(cl5!=0) {
806                 trac->AddClusterV2(5,(clind5[l5] & 0x0fffffff)>>0);
807                 trac->AddClusterMark(5,mark5[l5]);
808               }
809               if(cl4!=0){
810                 trac->AddClusterV2(4,(clind4[l4] & 0x0fffffff)>>0);
811                 trac->AddClusterMark(4,mark4[l4]);
812               }
813               if(cl3!=0){
814                 trac->AddClusterV2(3,(clind3[l3] & 0x0fffffff)>>0);
815                 trac->AddClusterMark(3,mark3[l3]);
816               }
817               if(cl2!=0){
818                 trac->AddClusterV2(2,(clind2[l2] & 0x0fffffff)>>0);
819                 trac->AddClusterMark(2,mark2[l2]);
820               }
821               if(cl1!=0){
822                 trac->AddClusterV2(1,(clind1[l1] & 0x0fffffff)>>0);
823                 trac->AddClusterMark(1,mark1[l1]);
824               }
825               if(cl0!=0){
826                 trac->AddClusterV2(0,(clind0[l0] & 0x0fffffff)>>0);
827                 trac->AddClusterMark(0,mark0[l0]);
828               }
829               //fit with Kalman filter using AliITStrackerMI::RefitAt()
830           
831
832               AliITStrackMI* ot = new AliITStrackSA(*trac);
833               
834               ot->ResetCovariance(10.);
835               ot->ResetClusters();
836               
837               if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
838                 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
839                 otrack2->ResetCovariance(10.); 
840                 otrack2->ResetClusters();
841                 //fit from layer 6 to layer 1
842                 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
843                   fListOfTracks->AddLast(otrack2);
844                   new (tri[nlist]) AliITStrackSA(*trac);
845                   nlist++;
846                 } else {
847                   delete otrack2;
848                 }
849                               
850               }       
851           
852               delete ot;
853               delete trac;
854             }//end loop layer 6
855           }//end loop layer 5
856         }//end loop layer 4        
857       }//end loop layer 3
858     }//end loop layer 2 
859   }//end loop layer 1
860
861   delete [] end;
862
863  
864
865   Int_t dim=fListOfTracks->GetEntries();
866   if(dim==0){
867     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
868       delete listlayer[i];
869     }
870     delete [] listlayer;
871     listSA->Delete();
872     delete listSA;
873     delete [] firstmod;
874     return 0;
875   }
876
877   Int_t lowchi2 = FindTrackLowChiSquare(fListOfTracks,dim);
878   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
879   AliITStrackSA* trsa = (AliITStrackSA*)listSA->At(lowchi2);
880  
881   if(otrack==0) {
882     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
883       delete listlayer[i];
884     }
885     delete [] listlayer; 
886     listSA->Delete();
887     delete listSA;
888     delete [] firstmod;
889     return 0;
890   }
891   Int_t * indexc = new Int_t[AliITSgeomTGeo::GetNLayers()];
892   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
893   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
894     indexc[nind] = otrack->GetClusterIndex(nind);
895   }      
896   AliITSRecPoint* cl0 = (AliITSRecPoint*)GetCluster(indexc[0]);
897   AliITSRecPoint* cl1 = (AliITSRecPoint*)GetCluster(indexc[1]);     
898   AliITSRecPoint* cl2 = (AliITSRecPoint*)GetCluster(indexc[2]);     
899   AliITSRecPoint* cl3 = (AliITSRecPoint*)GetCluster(indexc[3]);
900   AliITSRecPoint* cl4 = (AliITSRecPoint*)GetCluster(indexc[4]);
901   Int_t labl[3]={-1,-1,-1};
902   if(otrack->GetNumberOfClusters()==AliITSgeomTGeo::GetNLayers()){
903     AliITSRecPoint* cl5 = (AliITSRecPoint*)GetCluster(indexc[5]);
904     labl[0]=cl5->GetLabel(0);
905     labl[1]=cl5->GetLabel(1);
906     labl[2]=cl5->GetLabel(2);
907   }
908   delete [] indexc;
909   if(otrack->GetNumberOfClusters()==(AliITSgeomTGeo::GetNLayers()-1)){
910     labl[0]=-1;
911     labl[1]=-1;
912     labl[2]=-1;
913   }
914   Int_t minNPoints;
915   if(fSixPoints) minNPoints=6;
916   else minNPoints=5;
917   CookLabel(otrack,0.); //MI change - to see fake ratio
918   Int_t label =  Label(cl0->GetLabel(0),cl1->GetLabel(0), 
919                        cl2->GetLabel(0),cl3->GetLabel(0),
920                        cl4->GetLabel(0),labl[0],
921                        cl0->GetLabel(1),cl1->GetLabel(1),
922                        cl2->GetLabel(1),cl3->GetLabel(1),
923                        cl4->GetLabel(1),labl[1],
924                        cl0->GetLabel(2),cl1->GetLabel(2),
925                        cl2->GetLabel(2),cl3->GetLabel(2),
926                        cl4->GetLabel(2),labl[2],minNPoints);
927   
928   otrack->SetLabel(label);  
929   //remove clusters of found track
930
931   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
932     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
933       Int_t index = trsa->GetClusterMark(nlay,cln);
934       fCluLayer[nlay]->RemoveAt(index);
935       RemoveClusterCoord(nlay,index);
936       fCluLayer[nlay]->Compress();
937     }    
938   }
939   listSA->Delete();
940   delete listSA;
941
942   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
943     delete listlayer[i];
944   }
945   delete [] listlayer; 
946   delete [] firstmod;
947   return otrack;
948
949 }
950
951 //________________________________________________________________________
952
953 AliITStrackV2* AliITStrackerSA::FitShortTrack(AliITStrackSA* tr,Double_t *primaryVertex,Int_t innLay){
954   //fit of the found track (short tracks, <6 points, for cosmics). A.Dainese 31.07.07 
955
956   
957   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
958   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
959     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
960   }  
961
962   Int_t nclusters = tr->GetNumberOfClustersSA();
963   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
964   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
965     listlayer[i] = new TObjArray(0,0);
966   }
967
968   TArrayI clind0(20);
969   TArrayI clind1(20);
970   TArrayI clind2(20);
971   TArrayI clind3(20);
972   TArrayI clind4(20);
973   TArrayI clind5(20);
974
975   TArrayI mark0(20);
976   TArrayI mark1(20);
977   TArrayI mark2(20);
978   TArrayI mark3(20);
979   TArrayI mark4(20);
980   TArrayI mark5(20);
981
982
983   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
984   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
985   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
986   
987   for(Int_t ncl=0;ncl<nclusters;ncl++){
988     Int_t index = tr->GetClusterIndexSA(ncl); 
989     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
990     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
991     Int_t lay = (index & 0xf0000000) >> 28;
992     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
993     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
994     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
995     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
996     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
997     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}    
998   }
999   delete [] nnn;
1000
1001   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
1002     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
1003       Int_t mark = tr->GetClusterMark(nlay,ncl);
1004       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
1005       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
1006       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
1007       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
1008       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
1009       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
1010
1011     }
1012   }
1013
1014   delete [] kkk;
1015
1016  
1017   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
1018   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1019     if(listlayer[i]->GetEntries()==0) end[i]=1;
1020     else end[i]=listlayer[i]->GetEntries();
1021   }
1022
1023   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
1024   TClonesArray &tri = *listSA;
1025   Int_t nlist=0;
1026
1027   for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
1028     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l0); 
1029     for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
1030       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l1); 
1031       for(Int_t l2=0;l2<end[2];l2++){  //loop on layer 3
1032         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l2);
1033         for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4   
1034           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l3);
1035           for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
1036             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l4);
1037             for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6  
1038               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l5);
1039
1040
1041               Double_t x1,y1,z1,sx1,sy1,sz1;
1042               Double_t x2,y2,z2,sx2,sy2,sz2;
1043               AliITSRecPoint* p1=0;
1044               AliITSRecPoint* p2=0;
1045               Int_t index1=0,index2=0;
1046               Int_t mrk1=0,mrk2=0;
1047
1048               switch(innLay) {
1049               case 1: 
1050                 p1=cl1;
1051                 p2=cl2;
1052                 index1=clind1[l1];mrk1=mark1[l1];
1053                 index2=clind2[l2];mrk2=mark2[l2];
1054                 break;
1055               case 2:
1056                 p1=cl2;
1057                 p2=cl3;
1058                 index1=clind2[l2];mrk1=mark2[l2];
1059                 index2=clind3[l3];mrk2=mark3[l3];
1060                 break;
1061               case 3:
1062                 p1=cl3;
1063                 p2=cl4;
1064                 index1=clind3[l3];mrk1=mark3[l3];
1065                 index2=clind4[l4];mrk2=mark4[l4];
1066                 break;
1067               case 4:
1068                 p1=cl4;
1069                 p2=cl5;
1070                 index1=clind4[l4];mrk1=mark4[l4];
1071                 index2=clind5[l5];mrk2=mark5[l5];
1072                 break;
1073               }
1074               Int_t lay1=innLay;
1075               Int_t lay2=innLay+1;
1076               Int_t module1 = p1->GetDetectorIndex()+firstmod[innLay]; 
1077               
1078               Int_t cln1=mrk1;
1079               Int_t cln2=mrk2;
1080               AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(lay1,cln1);
1081               AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(lay2,cln2);
1082               x1 = arr1->GetX();
1083               x2 = arr2->GetX();
1084               y1 = arr1->GetY();
1085               y2 = arr2->GetY();
1086               z1 = arr1->GetZ();
1087               z2 = arr2->GetZ();
1088               sx1 = arr1->GetSx();
1089               sx2 = arr2->GetSx();
1090               sy1 = arr1->GetSy();
1091               sy2 = arr2->GetSy();
1092               sz1 = arr1->GetSz();
1093               sz2 = arr2->GetSz();
1094               
1095               Int_t layer,ladder,detector;
1096               AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
1097               Float_t yclu1 = p1->GetY();
1098               Float_t zclu1 = p1->GetZ();
1099               Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);        
1100               Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
1101               Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
1102               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
1103                               
1104               if(cl5!=0) {
1105                 trac->AddClusterV2(5,(clind5[l5] & 0x0fffffff)>>0);
1106                 trac->AddClusterMark(5,mark5[l5]);
1107               }
1108               if(cl4!=0){
1109                 trac->AddClusterV2(4,(clind4[l4] & 0x0fffffff)>>0);
1110                 trac->AddClusterMark(4,mark4[l4]);
1111               }
1112               if(cl3!=0){
1113                 trac->AddClusterV2(3,(clind3[l3] & 0x0fffffff)>>0);
1114                 trac->AddClusterMark(3,mark3[l3]);
1115               }
1116               if(cl2!=0){
1117                 trac->AddClusterV2(2,(clind2[l2] & 0x0fffffff)>>0);
1118                 trac->AddClusterMark(2,mark2[l2]);
1119               }
1120               if(cl1!=0){
1121                 trac->AddClusterV2(1,(clind1[l1] & 0x0fffffff)>>0);
1122                 trac->AddClusterMark(1,mark1[l1]);
1123               }
1124               if(cl0!=0){
1125                 trac->AddClusterV2(0,(clind0[l0] & 0x0fffffff)>>0);
1126                 trac->AddClusterMark(0,mark0[l0]);
1127               }
1128
1129
1130               //fit with Kalman filter using AliITStrackerMI::RefitAt()
1131               AliITStrackMI* ot = new AliITStrackSA(*trac);
1132               
1133               ot->ResetCovariance(10.);
1134               ot->ResetClusters();
1135               
1136               if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
1137                 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
1138                 otrack2->ResetCovariance(10.); 
1139                 otrack2->ResetClusters();
1140                 //fit from layer 6 to layer 1
1141                 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
1142                   fListOfTracks->AddLast(otrack2);
1143                   new (tri[nlist]) AliITStrackSA(*trac);
1144                   nlist++;
1145                 } else {
1146                   delete otrack2;
1147                 }
1148                               
1149               }       
1150           
1151               delete ot;
1152               delete trac;
1153             }//end loop layer 6
1154           }//end loop layer 5
1155         }//end loop layer 4        
1156       }//end loop layer 3
1157     }//end loop layer 2 
1158   }//end loop layer 1
1159
1160   delete [] end;
1161
1162  
1163
1164   Int_t dim=fListOfTracks->GetEntries();
1165   if(dim==0){
1166     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1167       delete listlayer[i];
1168     }
1169     delete [] listlayer;
1170     listSA->Delete();
1171     delete listSA;
1172     delete [] firstmod;
1173     return 0;
1174   }
1175
1176   Int_t lowchi2 = FindTrackLowChiSquare(fListOfTracks,dim);
1177   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
1178   AliITStrackSA* trsa = (AliITStrackSA*)listSA->At(lowchi2);
1179  
1180   if(otrack==0) {
1181     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1182       delete listlayer[i];
1183     }
1184     delete [] listlayer; 
1185     listSA->Delete();
1186     delete listSA;
1187     delete [] firstmod;
1188     return 0;
1189   }
1190   Int_t * indexc = new Int_t[AliITSgeomTGeo::GetNLayers()];
1191   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
1192   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
1193     indexc[nind] = otrack->GetClusterIndex(nind);
1194   }      
1195   Int_t labl[6][3];
1196   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
1197     if(i<otrack->GetNumberOfClusters()) {
1198       AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
1199       labl[i][0]=cl->GetLabel(0);
1200       labl[i][1]=cl->GetLabel(1);
1201       labl[i][2]=cl->GetLabel(2);
1202     } else {
1203       labl[i][0]=-1;
1204       labl[i][1]=-1;
1205       labl[i][2]=-1;
1206     }
1207   }
1208   delete [] indexc;
1209   CookLabel(otrack,0.); //MI change - to see fake ratio
1210  
1211   Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
1212   Int_t lflag=0;
1213   for(Int_t i=0;i<otrack->GetNumberOfClusters();i++)
1214     if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
1215   
1216   if(lflag<otrack->GetNumberOfClusters()) label = -label;
1217   otrack->SetLabel(label);  
1218
1219   //remove clusters of found track
1220   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
1221     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
1222       Int_t index = trsa->GetClusterMark(nlay,cln);
1223       fCluLayer[nlay]->RemoveAt(index);
1224       RemoveClusterCoord(nlay,index);
1225       fCluLayer[nlay]->Compress();
1226     }    
1227   }
1228   listSA->Delete();
1229   delete listSA;
1230
1231   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1232     delete listlayer[i];
1233   }
1234   delete [] listlayer; 
1235   delete [] firstmod;
1236   return otrack;
1237
1238 }
1239
1240
1241
1242 //________________________________________________________________________
1243
1244 AliITStrackV2* AliITStrackerSA::FitJumpyTrack(AliITStrackSA* tr,Double_t *primaryVertex) {
1245   //fit of the found track (jumpy tracks, <6 points, layers missing). A.Dainese 16.11.07 
1246
1247   
1248   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
1249   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1250     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
1251   }  
1252
1253   Int_t nclusters = tr->GetNumberOfClustersSA();
1254   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
1255   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1256     listlayer[i] = new TObjArray(0,0);
1257   }
1258
1259   TArrayI clind0(20);
1260   TArrayI clind1(20);
1261   TArrayI clind2(20);
1262   TArrayI clind3(20);
1263   TArrayI clind4(20);
1264   TArrayI clind5(20);
1265
1266   TArrayI mark0(20);
1267   TArrayI mark1(20);
1268   TArrayI mark2(20);
1269   TArrayI mark3(20);
1270   TArrayI mark4(20);
1271   TArrayI mark5(20);
1272
1273
1274   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
1275   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
1276   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
1277   
1278   for(Int_t ncl=0;ncl<nclusters;ncl++){
1279     Int_t index = tr->GetClusterIndexSA(ncl); 
1280     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
1281     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
1282     Int_t lay = (index & 0xf0000000) >> 28;
1283     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
1284     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
1285     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
1286     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
1287     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
1288     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}    
1289   }
1290   delete [] nnn;
1291
1292   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
1293     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
1294       Int_t mark = tr->GetClusterMark(nlay,ncl);
1295       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
1296       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
1297       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
1298       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
1299       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
1300       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
1301
1302     }
1303   }
1304
1305   delete [] kkk;
1306
1307   Int_t firstLay=-1,secondLay=-1;
1308   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
1309   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
1310     if(listlayer[i]->GetEntries()==0) {
1311       end[i]=1;
1312     } else {
1313       end[i]=listlayer[i]->GetEntries();
1314       if(firstLay==-1) {
1315         firstLay=i;
1316       } else if(secondLay==-1) {
1317         secondLay=i;
1318       }
1319     }
1320   }
1321
1322   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
1323   TClonesArray &tri = *listSA;
1324   Int_t nlist=0;
1325
1326   for(Int_t l0=0;l0<end[0];l0++){ //loop on layer 1
1327     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l0); 
1328     for(Int_t l1=0;l1<end[1];l1++){ //loop on layer 2
1329       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l1); 
1330       for(Int_t l2=0;l2<end[2];l2++){  //loop on layer 3
1331         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l2);
1332         for(Int_t l3=0;l3<end[3];l3++){ //loop on layer 4   
1333           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l3);
1334           for(Int_t l4=0;l4<end[4];l4++){ //loop on layer 5
1335             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l4);
1336             for(Int_t l5=0;l5<end[5];l5++){ //loop on layer 6  
1337               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l5);
1338
1339
1340               Double_t x1,y1,z1,sx1,sy1,sz1;
1341               Double_t x2,y2,z2,sx2,sy2,sz2;
1342               AliITSRecPoint* p1=0;
1343               AliITSRecPoint* p2=0;
1344               Int_t index1=0,index2=0;
1345               Int_t mrk1=0,mrk2=0;
1346
1347               switch(firstLay) {
1348               case 0:
1349                 p1=cl0;
1350                 index1=clind0[l0];mrk1=mark0[l0];
1351                 break;
1352               case 1:
1353                 p1=cl1;
1354                 index1=clind1[l1];mrk1=mark1[l1];
1355                 break;
1356               case 2:
1357                 p1=cl2;
1358                 index1=clind2[l2];mrk1=mark2[l2];
1359                 break;
1360               case 3:
1361                 p1=cl3;
1362                 index1=clind3[l3];mrk1=mark3[l3];
1363                 break;
1364               case 4:
1365                 p1=cl4;
1366                 index1=clind4[l4];mrk1=mark4[l4];
1367                 break;
1368               }
1369
1370               switch(secondLay) {
1371               case 1:
1372                 p2=cl1;
1373                 index2=clind1[l1];mrk2=mark1[l1];
1374                 break;
1375               case 2:
1376                 p2=cl2;
1377                 index2=clind2[l2];mrk2=mark2[l2];
1378                 break;
1379               case 3:
1380                 p2=cl3;
1381                 index2=clind3[l3];mrk2=mark3[l3];
1382                 break;
1383               case 4:
1384                 p2=cl4;
1385                 index2=clind4[l4];mrk2=mark4[l4];
1386                 break;
1387               case 5:
1388                 p2=cl5;
1389                 index2=clind5[l5];mrk2=mark5[l5];
1390                 break;
1391               }
1392
1393               Int_t module1 = p1->GetDetectorIndex()+firstmod[firstLay]; 
1394               
1395               Int_t cln1=mrk1;
1396               Int_t cln2=mrk2;
1397               AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(firstLay,cln1);
1398               AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(secondLay,cln2);
1399               x1 = arr1->GetX();
1400               x2 = arr2->GetX();
1401               y1 = arr1->GetY();
1402               y2 = arr2->GetY();
1403               z1 = arr1->GetZ();
1404               z2 = arr2->GetZ();
1405               sx1 = arr1->GetSx();
1406               sx2 = arr2->GetSx();
1407               sy1 = arr1->GetSy();
1408               sy2 = arr2->GetSy();
1409               sz1 = arr1->GetSz();
1410               sz2 = arr2->GetSz();
1411               
1412               Int_t layer,ladder,detector;
1413               AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
1414               Float_t yclu1 = p1->GetY();
1415               Float_t zclu1 = p1->GetZ();
1416               Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);        
1417               Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
1418               Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
1419               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
1420                               
1421               if(cl5!=0) {
1422                 trac->AddClusterV2(5,(clind5[l5] & 0x0fffffff)>>0);
1423                 trac->AddClusterMark(5,mark5[l5]);
1424               }
1425               if(cl4!=0){
1426                 trac->AddClusterV2(4,(clind4[l4] & 0x0fffffff)>>0);
1427                 trac->AddClusterMark(4,mark4[l4]);
1428               }
1429               if(cl3!=0){
1430                 trac->AddClusterV2(3,(clind3[l3] & 0x0fffffff)>>0);
1431                 trac->AddClusterMark(3,mark3[l3]);
1432               }
1433               if(cl2!=0){
1434                 trac->AddClusterV2(2,(clind2[l2] & 0x0fffffff)>>0);
1435                 trac->AddClusterMark(2,mark2[l2]);
1436               }
1437               if(cl1!=0){
1438                 trac->AddClusterV2(1,(clind1[l1] & 0x0fffffff)>>0);
1439                 trac->AddClusterMark(1,mark1[l1]);
1440               }
1441               if(cl0!=0){
1442                 trac->AddClusterV2(0,(clind0[l0] & 0x0fffffff)>>0);
1443                 trac->AddClusterMark(0,mark0[l0]);
1444               }
1445
1446
1447               //fit with Kalman filter using AliITStrackerMI::RefitAt()
1448               AliITStrackMI* ot = new AliITStrackSA(*trac);
1449               
1450               ot->ResetCovariance(10.);
1451               ot->ResetClusters();
1452               
1453               if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
1454                 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
1455                 otrack2->ResetCovariance(10.); 
1456                 otrack2->ResetClusters();
1457                 //fit from layer 6 to layer 1
1458                 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
1459                   fListOfTracks->AddLast(otrack2);
1460                   new (tri[nlist]) AliITStrackSA(*trac);
1461                   nlist++;
1462                 } else {
1463                   delete otrack2;
1464                 }
1465                               
1466               }       
1467           
1468               delete ot;
1469               delete trac;
1470             }//end loop layer 6
1471           }//end loop layer 5
1472         }//end loop layer 4        
1473       }//end loop layer 3
1474     }//end loop layer 2 
1475   }//end loop layer 1
1476
1477   delete [] end;
1478
1479  
1480
1481   Int_t dim=fListOfTracks->GetEntries();
1482   if(dim==0){
1483     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1484       delete listlayer[i];
1485     }
1486     delete [] listlayer;
1487     listSA->Delete();
1488     delete listSA;
1489     delete [] firstmod;
1490     return 0;
1491   }
1492
1493   Int_t lowchi2 = FindTrackLowChiSquare(fListOfTracks,dim);
1494   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
1495   AliITStrackSA* trsa = (AliITStrackSA*)listSA->At(lowchi2);
1496  
1497   if(otrack==0) {
1498     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1499       delete listlayer[i];
1500     }
1501     delete [] listlayer; 
1502     listSA->Delete();
1503     delete listSA;
1504     delete [] firstmod;
1505     return 0;
1506   }
1507   Int_t * indexc = new Int_t[AliITSgeomTGeo::GetNLayers()];
1508   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
1509   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
1510     indexc[nind] = otrack->GetClusterIndex(nind);
1511   }      
1512   Int_t labl[6][3];
1513   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
1514     if(i<otrack->GetNumberOfClusters()) {
1515       AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
1516       labl[i][0]=cl->GetLabel(0);
1517       labl[i][1]=cl->GetLabel(1);
1518       labl[i][2]=cl->GetLabel(2);
1519     } else {
1520       labl[i][0]=-1;
1521       labl[i][1]=-1;
1522       labl[i][2]=-1;
1523     }
1524   }
1525   delete [] indexc;
1526
1527   CookLabel(otrack,0.); //MI change - to see fake ratio
1528  
1529   Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
1530   Int_t lflag=0;
1531   for(Int_t i=0;i<otrack->GetNumberOfClusters();i++)
1532     if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
1533   
1534   if(lflag<otrack->GetNumberOfClusters()) label = -label;
1535   otrack->SetLabel(label);  
1536
1537   //remove clusters of found track
1538   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
1539     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
1540       Int_t index = trsa->GetClusterMark(nlay,cln);
1541       fCluLayer[nlay]->RemoveAt(index);
1542       RemoveClusterCoord(nlay,index);
1543       fCluLayer[nlay]->Compress();
1544     }    
1545   }
1546   listSA->Delete();
1547   delete listSA;
1548
1549   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1550     delete listlayer[i];
1551   }
1552   delete [] listlayer; 
1553   delete [] firstmod;
1554   return otrack;
1555
1556 }
1557
1558
1559
1560 //_______________________________________________________
1561 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
1562   //function used to to find the clusters associated to the track
1563
1564   if(AliITSReconstructor::GetRecoParam()->GetLayersToSkip(layer)) return 0;
1565
1566   Int_t nc=0;
1567   AliITSlayer &lay = fgLayers[layer];
1568   Double_t r=lay.GetR();
1569   if(pflag==1){      
1570     Float_t cx1,cx2,cy1,cy2;
1571     FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
1572     if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
1573        return 0;
1574     Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
1575     Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
1576     fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
1577   }
1578
1579  
1580   Int_t ncl = fCluLayer[layer]->GetEntries();
1581   for (Int_t index=0; index<ncl; index++) {
1582     AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
1583     if (!c) continue;
1584     if (c->GetQ()<=0) continue;
1585     
1586      AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
1587      Double_t phi = arr->GetPhi();
1588      if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
1589
1590      Double_t lambda = arr->GetLambda();
1591      if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
1592
1593      if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
1594      if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
1595      Int_t orind = arr->GetOrInd();
1596      trs->AddClusterSA(layer,orind);
1597      trs->AddClusterMark(layer,index);
1598        
1599      nc++;
1600      fLambdac=lambda;
1601      fPhiEstimate=phi;
1602
1603      fPointc[0]=arr->GetX();
1604      fPointc[1]=arr->GetY();
1605
1606   }
1607   return nc;
1608 }
1609
1610 //________________________________________________________________
1611 void AliITStrackerSA::UpdatePoints(){
1612   //update of points for the estimation of the curvature  
1613
1614   fPoint2[0]=fPoint3[0];
1615   fPoint2[1]=fPoint3[1];
1616   fPoint3[0]=fPointc[0];
1617   fPoint3[1]=fPointc[1];
1618
1619   
1620 }
1621
1622 //___________________________________________________________________
1623 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){
1624
1625    //given (x,y) of three recpoints (in global coordinates) 
1626    //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1627
1628    Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1629    if(den==0) return 0;
1630    a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1631    b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1632    c = -x1*x1-y1*y1-a*x1-b*y1;
1633    return 1;
1634  }
1635 //__________________________________________________________________________
1636  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){
1637  
1638  //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1639  //c2 is -rlayer*rlayer
1640
1641   if(a1==0) return 0;
1642  Double_t m = c2-c1; 
1643  Double_t aA = (b1*b1)/(a1*a1)+1;
1644  Double_t bB = (-2*m*b1/(a1*a1));
1645  Double_t cC = c2+(m*m)/(a1*a1);
1646  Double_t dD = bB*bB-4*aA*cC;
1647  if(dD<0) return 0;
1648  
1649  y1 = (-bB+TMath::Sqrt(dD))/(2*aA); 
1650  y2 = (-bB-TMath::Sqrt(dD))/(2*aA); 
1651  x1 = (c2-c1-b1*y1)/a1;
1652  x2 = (c2-c1-b1*y2)/a1;
1653
1654  return 1; 
1655 }
1656 //____________________________________________________________________
1657 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t 
1658 x2,Double_t y2,Double_t x3,Double_t y3){
1659
1660   //calculates the curvature of track  
1661   Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1662   if(den==0) return 0;
1663   Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1664   Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1665   Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1666   Double_t xc=-a/2.;
1667
1668   if((a*a+b*b-4*c)<0) return 0;
1669   Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1670   if(rad==0) return 0;
1671   
1672   if((x1>0 && y1>0 && x1<xc)) rad*=-1;
1673   if((x1<0 && y1>0 && x1<xc)) rad*=-1;
1674   //  if((x1<0 && y1<0 && x1<xc)) rad*=-1;
1675   // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1676   
1677   return 1/rad;
1678  
1679 }
1680
1681
1682 //____________________________________________________________________
1683 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1684
1685   //Returns the point closest to pp
1686
1687   Double_t diff1 = p1-pp;
1688   Double_t diff2 = p2-pp;
1689   
1690   if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1691   else fPhiEstimate=p2;  
1692   return fPhiEstimate;
1693   
1694 }
1695
1696
1697 //_________________________________________________________________
1698 Int_t AliITStrackerSA::FindTrackLowChiSquare(TObjArray* tracklist, Int_t dim) const {
1699   // returns track with lowers chi square  
1700   if(dim==1){
1701     //AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(0);
1702     //return trk;
1703     return 0;
1704   }
1705   if(dim==0) return 0;
1706   Double_t * chi2 = new Double_t[dim];
1707   Int_t * index = new Int_t[dim];
1708   for(Int_t i=0;i<dim;i++){
1709     AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(i);
1710     chi2[i]=trk->GetChi2();
1711     index[i]=i;
1712   }
1713
1714   Int_t w=0;Double_t value;
1715   Int_t lp;
1716   while(w<dim){
1717     for(Int_t j=w+1;j<dim;j++){
1718       if(chi2[w]<chi2[j]){
1719         value=chi2[w];
1720         chi2[w]=chi2[j];
1721         chi2[j]=value;
1722         lp=index[w];
1723         index[w]=index[j];
1724         index[j]=lp;
1725       }
1726     }
1727     w++;
1728   }
1729   Int_t tmp = index[dim-1];
1730   delete [] chi2;
1731   delete [] index;
1732   return tmp;
1733 }
1734
1735 //__________________________________________________________
1736 Int_t AliITStrackerSA::FindLabel(Int_t l0, Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5){
1737
1738   //function used to determine the track label
1739   
1740   Int_t lb[6] = {l0,l1,l2,l3,l4,l5};
1741   Int_t aa[6]={1,1,1,1,1,1};
1742   Int_t ff=0; 
1743   Int_t ll=0;
1744   Int_t k=0;Int_t w=0;Int_t num=6;
1745   for(Int_t i=5;i>=0;i--) if(lb[i]==-1) num=i;
1746   
1747   while(k<num){
1748   
1749     for(Int_t i=k+1;i<num;i++){
1750     
1751       if(lb[k]==lb[i] && aa[k]!=0){
1752       
1753         aa[k]+=1;
1754         aa[i]=0;
1755       }
1756     }
1757     k++;
1758   }
1759
1760   while(w<num){
1761  
1762     for(Int_t j=0;j<6;j++){
1763       if(aa[w]<aa[j]){
1764       ff=aa[w];
1765       aa[w]=aa[j];
1766       aa[j]=ff;
1767       ll=lb[w];
1768       lb[w]=lb[j];
1769       lb[j]=ll;
1770      }
1771     }
1772     w++;
1773   }
1774   
1775   if(num<1) return -1; 
1776   return lb[num-1];
1777 }
1778
1779 //_____________________________________________________________________________
1780 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,
1781 Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t minNPoints){
1782
1783  
1784   //function used to assign label to the found track. If track is fake, the label is negative
1785
1786   Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1787   Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1788   Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1789   Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1790   Int_t lflag=0;Int_t num=6;
1791   if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1792
1793   for(Int_t i=0;i<num;i++){
1794     if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1795   }
1796
1797   if(lflag>=minNPoints) return ll;
1798   else return -ll;
1799
1800   
1801 }
1802
1803 //_____________________________________________________________________________
1804 void AliITStrackerSA::SetWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1805   // Set sizes of the phi and lambda windows used for track finding
1806   fNloop = n;
1807   if(phi){ // user defined values
1808     fPhiWin = new Double_t[fNloop];
1809     fLambdaWin = new Double_t[fNloop];
1810     for(Int_t k=0;k<fNloop;k++){
1811       fPhiWin[k]=phi[k];
1812       fLambdaWin[k]=lam[k];
1813     }
1814   }
1815   else {  // default values
1816             
1817     Double_t phid[33]   = {0.002,0.003,0.004,0.0045,0.0047,
1818                            0.005,0.0053,0.0055,
1819                            0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1820                            0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1821                            0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1822     Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1823                             0.005,0.005,0.006,
1824                             0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1825                             0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1826                             0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1827     
1828     if(fNloop!=33){
1829       fNloop = 33;
1830     }
1831     
1832     
1833     fPhiWin = new Double_t[fNloop];
1834     fLambdaWin = new Double_t[fNloop];
1835
1836     Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
1837   
1838     for(Int_t k=0;k<fNloop;k++){
1839       fPhiWin[k]=phid[k]*factor;
1840       fLambdaWin[k]=lambdad[k]*factor;
1841     }
1842   
1843   }
1844
1845 }
1846 //_______________________________________________________________________
1847 void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
1848   //Returns values of phi (azimuthal) and lambda angles for a given cluster
1849 /*  
1850   Double_t rot[9];     fGeom->GetRotMatrix(module,rot);
1851   Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1852   Float_t tx,ty,tz;  fGeom->GetTrans(lay,lad,det,tx,ty,tz);     
1853
1854   Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1855   Double_t phi1=TMath::Pi()/2+alpha;
1856   if (lay==1) phi1+=TMath::Pi();
1857
1858   Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1859   Float_t r=tx*cp+ty*sp;
1860
1861   xyz= r*cp - cl->GetY()*sp;
1862   y= r*sp + cl->GetY()*cp;
1863   z=cl->GetZ();
1864 */
1865   Float_t xyz[3];
1866   cl->GetGlobalXYZ(xyz);
1867   x=xyz[0];
1868   y=xyz[1];
1869   z=xyz[2];
1870  
1871   phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
1872   lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1873 }
1874
1875 //________________________________________________________________________
1876 void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
1877
1878   //returns sigmax, y, z of cluster in global coordinates
1879 /*
1880   Double_t rot[9];     fGeom->GetRotMatrix(module,rot);
1881   Int_t lay,lad,det; 
1882   AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
1883  
1884   Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1885   Double_t phi=TMath::Pi()/2+alpha;
1886   if (lay==1) phi+=TMath::Pi();
1887
1888   Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
1889 */
1890   Float_t covm[6];
1891   cl->GetGlobalCov(covm);
1892   sx=TMath::Sqrt(covm[0]);
1893   sy=TMath::Sqrt(covm[3]);
1894   sz=TMath::Sqrt(covm[5]);
1895 /*
1896   sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1897   sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1898   sz = TMath::Sqrt(cl->GetSigmaZ2());
1899 */
1900 }
1901