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