Possibility of selecting larger windows in SA tracker, in order to use the tracker...
[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
388    Int_t ntrack=0;
389    //loop on the different windows
390    for(Int_t nloop=0;nloop<fNloop;nloop++){
391      for(Int_t ncl=0;ncl<fCluLayer[0]->GetEntries();ncl++){ //loop starting from layer 0
392        
393        ResetForFinding();
394        Int_t pflag=0;
395        
396        AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[0]->At(ncl);
397
398        if(!cl) continue;
399
400
401        if (cl->GetQ()<=0) continue;
402        
403        AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(0,ncl); 
404        fPhic = arr->GetPhi();
405        fLambdac = arr->GetLambda();
406        if (TMath::Abs(fLambdac)>0.26*TMath::Pi()) continue;
407        fPhiEstimate = fPhic;
408        AliITStrackSA* trs = new AliITStrackSA(); 
409        fPoint1[0]=primaryVertex[0];
410        fPoint1[1]=primaryVertex[1];
411
412
413        fPoint2[0]=arr->GetX();
414        fPoint2[1]=arr->GetY();
415        Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()];//counter for clusters on each layer
416        for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){ nn[i]=0;}
417        nn[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
418        
419        nn[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
420       if(nn[1]>0){
421         pflag=1;
422         fPoint3[0] = fPointc[0];
423         fPoint3[1] = fPointc[1];
424       }
425       nn[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
426       if(nn[1]==0 && nn[2]==0) pflag=0;
427       if(nn[2]!=0 && nn[1]!=0){ pflag=1; UpdatePoints();}
428       if(nn[2]!=0 && nn[1]==0){
429         pflag=1;
430         fPoint3[0]=fPointc[0];
431         fPoint3[1]=fPointc[1];
432       }
433
434       nn[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
435       pflag=1;
436       if(nn[3]!=0) UpdatePoints();
437       nn[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag); 
438       pflag=1;
439       if(nn[4]!=0) UpdatePoints();
440       nn[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag); 
441           
442
443       Int_t layOK=0;
444       Int_t numberofpoints;
445       if(fSixPoints) numberofpoints=6;  //check of the candidate track
446       else numberofpoints=5;           //if track is good (with the required number        
447       for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers();nnp++){    //of points) it is written on file
448         if(nn[nnp]!=0) layOK+=1;
449       }
450       if(layOK>=numberofpoints){
451
452         AliITStrackV2* tr2 = FitTrack(trs,primaryVertex);
453         
454         if(tr2==0) continue;
455
456         AliESDtrack outtrack;
457         outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
458         event->AddTrack(&outtrack);
459         ntrack++;
460       }
461         
462
463       delete trs;
464       delete[] nn;
465       
466      }//end loop on clusters of layer1
467     
468      //end loop2
469    }
470    
471    //if 5/6 points are required, second loop starting 
472    //from second layer (SPD2), to find tracks with point of 
473    //layer 1 missing   
474    if(!fSixPoints) {
475      //printf("looking from SPD2\n");
476      //   counter for clusters on each layer  
477      Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()-1];      
478      for(Int_t nloop=0;nloop<fNloop;nloop++){
479        Int_t ncl2=fCluLayer[1]->GetEntries();
480        while(ncl2--){ //loop starting from layer 2
481          ResetForFinding();
482          Int_t pflag=0;
483          AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[1]->At(ncl2);
484          
485          if(!cl) continue;
486          AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(1,ncl2);
487          fPhic = arr->GetPhi();
488          fLambdac = arr->GetLambda();
489          fPhiEstimate = fPhic;
490          
491          AliITStrackSA* trs = new AliITStrackSA(); 
492          fPoint1[0]=primaryVertex[0];
493          fPoint1[1]=primaryVertex[1];
494          
495          fPoint2[0]=arr->GetX();
496          fPoint2[1]=arr->GetY();
497          for(Int_t kk=0;kk<AliITSgeomTGeo::GetNLayers()-1;kk++)nn[kk] = 0;
498          nn[0] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],
499                                    trs,primaryVertex[2],pflag);
500          nn[1] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],
501                                    trs,primaryVertex[2],pflag);
502          if(nn[1]!=0){
503            pflag=1;
504            fPoint3[0]=fPointc[0];
505            fPoint3[1]=fPointc[1];
506          }
507          nn[2]= SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],
508                                    trs,primaryVertex[2],pflag);
509          if(nn[2]!=0){
510            pflag=1;
511            UpdatePoints();
512          }
513          nn[3]= SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],
514                                    trs,primaryVertex[2],pflag);
515          if(nn[3]!=0){
516            pflag=1;
517            UpdatePoints();
518          }
519          nn[4]=SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],
520                                    trs,primaryVertex[2],pflag);
521
522          Int_t fl=0;
523          for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-1;nnp++){
524            if(nn[nnp]!=0) fl+=1;
525          }
526          if(fl>=5){  // 5/6       
527            AliITStrackV2* tr2 = FitTrack(trs,primaryVertex);
528            if(tr2==0){
529              continue;
530            }
531            
532            AliESDtrack outtrack;
533            outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
534            event->AddTrack(&outtrack);
535            ntrack++;
536
537          }   
538               
539         delete trs;
540       }//end loop on clusters of layer2
541      }
542    
543     delete [] nn;
544   } //end opt="5/6"  
545
546
547
548    // search for tracks starting from SPD2, SDD1, SDD2, SSD2
549    // for cosmics (A. Dainese 31.07.07)
550    if(fOuterStartLayer>0) {
551      for(Int_t innLay=1; innLay<=fOuterStartLayer; innLay++) {
552        printf("Searching from layer %d outward\n",innLay);
553        //   counter for clusters on each layer  
554        Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()-innLay];      
555        for(Int_t nloop=0;nloop<fNloop;nloop++){
556          Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
557          while(nclInnLay--){ //loop starting from layer innLay
558            ResetForFinding();
559            Int_t pflag=0;
560            AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
561          
562            if(!cl) continue;
563            AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
564            fPhic = arr->GetPhi();
565            fLambdac = arr->GetLambda();
566            fPhiEstimate = fPhic;
567          
568            AliITStrackSA* trs = new AliITStrackSA(); 
569            fPoint1[0]=primaryVertex[0];
570            fPoint1[1]=primaryVertex[1];
571            fPoint2[0]=arr->GetX();
572            fPoint2[1]=arr->GetY();
573
574            Int_t kk;
575            for(kk=0;kk<AliITSgeomTGeo::GetNLayers()-innLay;kk++) nn[kk] = 0;
576
577            kk=0;
578            nn[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
579                                   trs,primaryVertex[2],pflag);
580            for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
581              kk++;
582              nn[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
583                                      trs,primaryVertex[2],pflag);
584              if(nn[kk]!=0){
585                pflag=1;
586                if(kk==1) {
587                  fPoint3[0]=fPointc[0];
588                  fPoint3[1]=fPointc[1];
589                } else {
590                  UpdatePoints();
591                }
592              }
593            }
594
595            Int_t fl=0;
596            for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
597              if(nn[nnp]!=0) fl+=1;
598            }
599            if(fl>=AliITSgeomTGeo::GetNLayers()-innLay){ 
600              AliITStrackV2* tr2 = FitShortTrack(trs,primaryVertex,innLay);
601              if(tr2==0){
602                continue;
603              }
604              
605              AliESDtrack outtrack;
606              outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
607              event->AddTrack(&outtrack);
608              ntrack++;
609              
610            }   
611               
612            delete trs;
613          }//end loop on clusters of innLay
614        } //end loop on window sizes
615    
616        delete [] nn;
617      } //end loop on innLay
618    } //end if(fOuterStartLayer>0)
619
620    
621   delete [] firstmod;
622   Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
623   return 0;
624
625 }
626  
627
628
629
630
631 //________________________________________________________________________
632
633 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex){
634   //fit of the found track
635
636   
637   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
638   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
639     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
640   }  
641
642   Int_t nclusters = tr->GetNumberOfClustersSA();
643   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
644   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
645     listlayer[i] = new TObjArray(0,0);
646   }
647
648   TArrayI clind0(20);
649   TArrayI clind1(20);
650   TArrayI clind2(20);
651   TArrayI clind3(20);
652   TArrayI clind4(20);
653   TArrayI clind5(20);
654
655   TArrayI mark0(20);
656   TArrayI mark1(20);
657   TArrayI mark2(20);
658   TArrayI mark3(20);
659   TArrayI mark4(20);
660   TArrayI mark5(20);
661
662
663   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
664   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
665   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
666   
667   for(Int_t ncl=0;ncl<nclusters;ncl++){
668     Int_t index = tr->GetClusterIndexSA(ncl); 
669     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
670     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
671     Int_t lay = (index & 0xf0000000) >> 28;
672     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
673     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
674     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
675     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
676     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
677     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}    
678   }
679   delete [] nnn;
680
681   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
682     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
683       Int_t mark = tr->GetClusterMark(nlay,ncl);
684       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
685       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
686       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
687       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
688       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
689       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
690
691     }
692   }
693
694   delete [] kkk;
695
696  
697   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
698   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
699     if(listlayer[i]->GetEntries()==0) end[i]=1;
700     else end[i]=listlayer[i]->GetEntries();
701   }
702
703   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
704   TClonesArray &tri = *listSA;
705   Int_t nlist=0;
706
707   if(end[0]==0) end[0]=1; //for tracks with cluster on layer 0 missing
708   for(Int_t l1=0;l1<end[0];l1++){ //loop on layer 1
709     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l1); 
710     Double_t x1,y1,z1,sx1,sy1,sz1;
711     Double_t x2,y2,z2,sx2,sy2,sz2;
712     AliITSRecPoint* p1=0;
713     AliITSRecPoint* p2=0;
714     Int_t index1=clind0[l1];
715     Int_t index2=0;
716     Int_t mrk1 = mark0[l1];
717     Int_t mrk2 = 0;
718     Int_t lay1=0;
719     Int_t lay2=1;
720     Int_t module1=-1;
721     for(Int_t l2=0;l2<end[1];l2++){ //loop on layer 2
722       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l2); 
723       index2=clind1[l2];
724       mrk2 = mark1[l2];
725       for(Int_t l3=0;l3<end[2];l3++){  //loop on layer 3
726         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l3);
727
728         if(cl0==0 && cl1!=0) {
729           p1=cl1;
730           p2=cl2;
731           index1=clind1[l2];mrk1=mark1[l2];
732           index2=clind2[l3];mrk2=mark2[l3];
733           lay1=1;
734           lay2=2;
735           module1 = p1->GetDetectorIndex()+firstmod[1]; 
736         }
737         if(cl0!=0 && cl1==0){
738           p1=cl0;
739           p2=cl2;index2=clind2[l3];mrk2=mark2[l3];
740           lay1=0;
741           lay2=2;
742           module1 = p1->GetDetectorIndex()+firstmod[0];
743         }
744         if(cl0!=0 && cl1!=0){
745           p1=cl0;
746           p2=cl1;
747           lay1=0;
748           lay2=1;
749           module1 = p1->GetDetectorIndex()+firstmod[0];
750         }
751         
752         Int_t cln1=mrk1;
753         Int_t cln2=mrk2;
754         AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(lay1,cln1);
755         AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(lay2,cln2);
756         x1 = arr1->GetX();
757         x2 = arr2->GetX();
758         y1 = arr1->GetY();
759         y2 = arr2->GetY();
760         z1 = arr1->GetZ();
761         z2 = arr2->GetZ();
762         sx1 = arr1->GetSx();
763         sx2 = arr2->GetSx();
764         sy1 = arr1->GetSy();
765         sy2 = arr2->GetSy();
766         sz1 = arr1->GetSz();
767         sz2 = arr2->GetSz();
768         
769         Int_t layer,ladder,detector;
770         AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
771         Float_t yclu1 = p1->GetY();
772         Float_t zclu1 = p1->GetZ();
773         Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);        
774         Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
775         Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
776
777         for(Int_t l4=0;l4<end[3];l4++){ //loop on layer 4   
778           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l4);
779           for(Int_t l5=0;l5<end[4];l5++){ //loop on layer 5
780             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l5);
781             for(Int_t l6=0;l6<end[5];l6++){ //loop on layer 6  
782               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l6);
783               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
784                               
785               if(cl5!=0) {
786                 trac->AddClusterV2(5,(clind5[l6] & 0x0fffffff)>>0);
787                 trac->AddClusterMark(5,mark5[l6]);
788               }
789               if(cl4!=0){
790                 trac->AddClusterV2(4,(clind4[l5] & 0x0fffffff)>>0);
791                 trac->AddClusterMark(4,mark4[l5]);
792               }
793               if(cl3!=0){
794                 trac->AddClusterV2(3,(clind3[l4] & 0x0fffffff)>>0);
795                 trac->AddClusterMark(3,mark3[l4]);
796               }
797               if(cl2!=0){
798                 trac->AddClusterV2(2,(clind2[l3] & 0x0fffffff)>>0);
799                 trac->AddClusterMark(2,mark2[l3]);
800               }
801               if(cl1!=0){
802                 trac->AddClusterV2(1,(clind1[l2] & 0x0fffffff)>>0);
803                 trac->AddClusterMark(1,mark1[l2]);
804               }
805               if(cl0!=0){
806                 trac->AddClusterV2(0,(clind0[l1] & 0x0fffffff)>>0);
807                 trac->AddClusterMark(0,mark0[l1]);
808               }
809               //fit with Kalman filter using AliITStrackerMI::RefitAt()
810           
811
812               AliITStrackMI* ot = new AliITStrackSA(*trac);
813               
814               ot->ResetCovariance(10.);
815               ot->ResetClusters();
816               
817               if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
818                 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
819                 otrack2->ResetCovariance(10.); 
820                 otrack2->ResetClusters();
821                 //fit from layer 6 to layer 1
822                 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
823                   fListOfTracks->AddLast(otrack2);
824                   new (tri[nlist]) AliITStrackSA(*trac);
825                   nlist++;
826                 } else {
827                   delete otrack2;
828                 }
829                               
830               }       
831           
832               delete ot;
833               delete trac;
834             }//end loop layer 6
835           }//end loop layer 5
836         }//end loop layer 4        
837       }//end loop layer 3
838     }//end loop layer 2 
839   }//end loop layer 1
840
841   delete [] end;
842
843  
844
845   Int_t dim=fListOfTracks->GetEntries();
846   if(dim==0){
847     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
848       delete listlayer[i];
849     }
850     delete [] listlayer;
851     listSA->Delete();
852     delete listSA;
853     delete [] firstmod;
854     return 0;
855   }
856
857   Int_t lowchi2 = FindTrackLowChiSquare(fListOfTracks,dim);
858   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
859   AliITStrackSA* trsa = (AliITStrackSA*)listSA->At(lowchi2);
860  
861   if(otrack==0) {
862     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
863       delete listlayer[i];
864     }
865     delete [] listlayer; 
866     listSA->Delete();
867     delete listSA;
868     delete [] firstmod;
869     return 0;
870   }
871   Int_t * indexc = new Int_t[AliITSgeomTGeo::GetNLayers()];
872   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
873   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
874     indexc[nind] = otrack->GetClusterIndex(nind);
875   }      
876   AliITSRecPoint* cl0 = (AliITSRecPoint*)GetCluster(indexc[0]);
877   AliITSRecPoint* cl1 = (AliITSRecPoint*)GetCluster(indexc[1]);     
878   AliITSRecPoint* cl2 = (AliITSRecPoint*)GetCluster(indexc[2]);     
879   AliITSRecPoint* cl3 = (AliITSRecPoint*)GetCluster(indexc[3]);
880   AliITSRecPoint* cl4 = (AliITSRecPoint*)GetCluster(indexc[4]);
881   Int_t labl[3]={-1,-1,-1};
882   if(otrack->GetNumberOfClusters()==AliITSgeomTGeo::GetNLayers()){
883     AliITSRecPoint* cl5 = (AliITSRecPoint*)GetCluster(indexc[5]);
884     labl[0]=cl5->GetLabel(0);
885     labl[1]=cl5->GetLabel(1);
886     labl[2]=cl5->GetLabel(2);
887   }
888   delete [] indexc;
889   if(otrack->GetNumberOfClusters()==(AliITSgeomTGeo::GetNLayers()-1)){
890     labl[0]=-1;
891     labl[1]=-1;
892     labl[2]=-1;
893   }
894   Int_t numberofpoints;
895   if(fSixPoints) numberofpoints=6;
896   else numberofpoints=5;
897   CookLabel(otrack,0.); //MI change - to see fake ratio
898   Int_t label =  Label(cl0->GetLabel(0),cl1->GetLabel(0), 
899                        cl2->GetLabel(0),cl3->GetLabel(0),
900                        cl4->GetLabel(0),labl[0],
901                        cl0->GetLabel(1),cl1->GetLabel(1),
902                        cl2->GetLabel(1),cl3->GetLabel(1),
903                        cl4->GetLabel(1),labl[1],
904                        cl0->GetLabel(2),cl1->GetLabel(2),
905                        cl2->GetLabel(2),cl3->GetLabel(2),
906                        cl4->GetLabel(2),labl[2],numberofpoints);
907   
908   otrack->SetLabel(label);  
909   //remove clusters of found track
910
911   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
912     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
913       Int_t index = trsa->GetClusterMark(nlay,cln);
914       fCluLayer[nlay]->RemoveAt(index);
915       RemoveClusterCoord(nlay,index);
916       fCluLayer[nlay]->Compress();
917     }    
918   }
919   listSA->Delete();
920   delete listSA;
921
922   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
923     delete listlayer[i];
924   }
925   delete [] listlayer; 
926   delete [] firstmod;
927   return otrack;
928
929 }
930
931 //________________________________________________________________________
932
933 AliITStrackV2* AliITStrackerSA::FitShortTrack(AliITStrackSA* tr,Double_t *primaryVertex,Int_t innLay){
934   //fit of the found track (short tracks, <6 points, for cosmics). A.Dainese 31.07.07 
935
936   
937   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
938   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
939     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
940   }  
941
942   Int_t nclusters = tr->GetNumberOfClustersSA();
943   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
944   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
945     listlayer[i] = new TObjArray(0,0);
946   }
947
948   TArrayI clind0(20);
949   TArrayI clind1(20);
950   TArrayI clind2(20);
951   TArrayI clind3(20);
952   TArrayI clind4(20);
953   TArrayI clind5(20);
954
955   TArrayI mark0(20);
956   TArrayI mark1(20);
957   TArrayI mark2(20);
958   TArrayI mark3(20);
959   TArrayI mark4(20);
960   TArrayI mark5(20);
961
962
963   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
964   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
965   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
966   
967   for(Int_t ncl=0;ncl<nclusters;ncl++){
968     Int_t index = tr->GetClusterIndexSA(ncl); 
969     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
970     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
971     Int_t lay = (index & 0xf0000000) >> 28;
972     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
973     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
974     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
975     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
976     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
977     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}    
978   }
979   delete [] nnn;
980
981   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
982     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
983       Int_t mark = tr->GetClusterMark(nlay,ncl);
984       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
985       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
986       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
987       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
988       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
989       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
990
991     }
992   }
993
994   delete [] kkk;
995
996  
997   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
998   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
999     if(listlayer[i]->GetEntries()==0) end[i]=1;
1000     else end[i]=listlayer[i]->GetEntries();
1001   }
1002
1003   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
1004   TClonesArray &tri = *listSA;
1005   Int_t nlist=0;
1006
1007   //  if(end[0]==0) end[0]=1; //for tracks with cluster on layer 0 missing
1008   for(Int_t l1=0;l1<end[0];l1++){ //loop on layer 1
1009     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l1); 
1010     for(Int_t l2=0;l2<end[1];l2++){ //loop on layer 2
1011       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l2); 
1012       for(Int_t l3=0;l3<end[2];l3++){  //loop on layer 3
1013         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l3);
1014         for(Int_t l4=0;l4<end[3];l4++){ //loop on layer 4   
1015           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l4);
1016           for(Int_t l5=0;l5<end[4];l5++){ //loop on layer 5
1017             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l5);
1018             for(Int_t l6=0;l6<end[5];l6++){ //loop on layer 6  
1019               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l6);
1020
1021
1022               Double_t x1,y1,z1,sx1,sy1,sz1;
1023               Double_t x2,y2,z2,sx2,sy2,sz2;
1024               AliITSRecPoint* p1=0;
1025               AliITSRecPoint* p2=0;
1026               Int_t index1=0,index2=0;
1027               Int_t mrk1=0,mrk2=0;
1028               if(innLay==1) {
1029                 p1=cl1;
1030                 p2=cl2;
1031                 index1=clind1[l2];mrk1=mark1[l2];
1032                 index2=clind2[l3];mrk2=mark2[l3];
1033               } else if(innLay==2) {
1034                 p1=cl2;
1035                 p2=cl3;
1036                 index1=clind2[l3];mrk1=mark2[l3];
1037                 index2=clind3[l4];mrk2=mark3[l4];
1038               } else if(innLay==3) {
1039                 p1=cl3;
1040                 p2=cl4;
1041                 index1=clind3[l4];mrk1=mark3[l4];
1042                 index2=clind4[l5];mrk2=mark4[l5];
1043               } else if(innLay==4) {
1044                 p1=cl4;
1045                 p2=cl5;
1046                 index1=clind4[l5];mrk1=mark4[l5];
1047                 index2=clind5[l6];mrk2=mark5[l6];
1048               }
1049               Int_t lay1=innLay;
1050               Int_t lay2=innLay+1;
1051               Int_t module1 = p1->GetDetectorIndex()+firstmod[innLay]; 
1052               
1053               Int_t cln1=mrk1;
1054               Int_t cln2=mrk2;
1055               AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(lay1,cln1);
1056               AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(lay2,cln2);
1057               x1 = arr1->GetX();
1058               x2 = arr2->GetX();
1059               y1 = arr1->GetY();
1060               y2 = arr2->GetY();
1061               z1 = arr1->GetZ();
1062               z2 = arr2->GetZ();
1063               sx1 = arr1->GetSx();
1064               sx2 = arr2->GetSx();
1065               sy1 = arr1->GetSy();
1066               sy2 = arr2->GetSy();
1067               sz1 = arr1->GetSz();
1068               sz2 = arr2->GetSz();
1069               
1070               Int_t layer,ladder,detector;
1071               AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
1072               Float_t yclu1 = p1->GetY();
1073               Float_t zclu1 = p1->GetZ();
1074               Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);        
1075               Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
1076               Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
1077               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
1078                               
1079               if(cl5!=0) {
1080                 trac->AddClusterV2(5,(clind5[l6] & 0x0fffffff)>>0);
1081                 trac->AddClusterMark(5,mark5[l6]);
1082               }
1083               if(cl4!=0){
1084                 trac->AddClusterV2(4,(clind4[l5] & 0x0fffffff)>>0);
1085                 trac->AddClusterMark(4,mark4[l5]);
1086               }
1087               if(cl3!=0){
1088                 trac->AddClusterV2(3,(clind3[l4] & 0x0fffffff)>>0);
1089                 trac->AddClusterMark(3,mark3[l4]);
1090               }
1091               if(cl2!=0){
1092                 trac->AddClusterV2(2,(clind2[l3] & 0x0fffffff)>>0);
1093                 trac->AddClusterMark(2,mark2[l3]);
1094               }
1095               if(cl1!=0){
1096                 trac->AddClusterV2(1,(clind1[l2] & 0x0fffffff)>>0);
1097                 trac->AddClusterMark(1,mark1[l2]);
1098               }
1099               if(cl0!=0){
1100                 trac->AddClusterV2(0,(clind0[l1] & 0x0fffffff)>>0);
1101                 trac->AddClusterMark(0,mark0[l1]);
1102               }
1103
1104
1105               //fit with Kalman filter using AliITStrackerMI::RefitAt()
1106               AliITStrackMI* ot = new AliITStrackSA(*trac);
1107               
1108               ot->ResetCovariance(10.);
1109               ot->ResetClusters();
1110               
1111               if(RefitAt(AliITSRecoParam::GetrInsideITSscreen(),ot,trac)){ //fit from layer 1 to layer 6
1112                 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
1113                 otrack2->ResetCovariance(10.); 
1114                 otrack2->ResetClusters();
1115                 //fit from layer 6 to layer 1
1116                 if(RefitAt(AliITSRecoParam::GetrInsideSPD1(),otrack2,ot)) {
1117                   fListOfTracks->AddLast(otrack2);
1118                   new (tri[nlist]) AliITStrackSA(*trac);
1119                   nlist++;
1120                 } else {
1121                   delete otrack2;
1122                 }
1123                               
1124               }       
1125           
1126               delete ot;
1127               delete trac;
1128             }//end loop layer 6
1129           }//end loop layer 5
1130         }//end loop layer 4        
1131       }//end loop layer 3
1132     }//end loop layer 2 
1133   }//end loop layer 1
1134
1135   delete [] end;
1136
1137  
1138
1139   Int_t dim=fListOfTracks->GetEntries();
1140   if(dim==0){
1141     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1142       delete listlayer[i];
1143     }
1144     delete [] listlayer;
1145     listSA->Delete();
1146     delete listSA;
1147     delete [] firstmod;
1148     return 0;
1149   }
1150
1151   Int_t lowchi2 = FindTrackLowChiSquare(fListOfTracks,dim);
1152   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
1153   AliITStrackSA* trsa = (AliITStrackSA*)listSA->At(lowchi2);
1154  
1155   if(otrack==0) {
1156     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1157       delete listlayer[i];
1158     }
1159     delete [] listlayer; 
1160     listSA->Delete();
1161     delete listSA;
1162     delete [] firstmod;
1163     return 0;
1164   }
1165   Int_t * indexc = new Int_t[AliITSgeomTGeo::GetNLayers()];
1166   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
1167   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
1168     indexc[nind] = otrack->GetClusterIndex(nind);
1169   }      
1170   Int_t labl[6][3];
1171   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {
1172     if(i<otrack->GetNumberOfClusters()) {
1173       AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(indexc[i]);
1174       labl[i][0]=cl->GetLabel(0);
1175       labl[i][1]=cl->GetLabel(1);
1176       labl[i][2]=cl->GetLabel(2);
1177     } else {
1178       labl[i][0]=-1;
1179       labl[i][1]=-1;
1180       labl[i][2]=-1;
1181     }
1182   }
1183   delete [] indexc;
1184   Int_t numberofpoints=AliITSgeomTGeo::GetNLayers()-innLay;
1185   CookLabel(otrack,0.); //MI change - to see fake ratio
1186  
1187   Int_t label=FindLabel(labl[0][0],labl[1][0],labl[2][0],labl[3][0],labl[4][0],labl[5][0]);
1188   Int_t lflag=0;
1189   for(Int_t i=0;i<numberofpoints;i++)
1190     if(labl[i][0]==label || labl[i][1]==label || labl[i][2]==label) lflag++;
1191   
1192   if(lflag<numberofpoints) label = -label;
1193   otrack->SetLabel(label);  
1194
1195   //remove clusters of found track
1196   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
1197     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
1198       Int_t index = trsa->GetClusterMark(nlay,cln);
1199       fCluLayer[nlay]->RemoveAt(index);
1200       RemoveClusterCoord(nlay,index);
1201       fCluLayer[nlay]->Compress();
1202     }    
1203   }
1204   listSA->Delete();
1205   delete listSA;
1206
1207   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
1208     delete listlayer[i];
1209   }
1210   delete [] listlayer; 
1211   delete [] firstmod;
1212   return otrack;
1213
1214 }
1215
1216
1217
1218 //_______________________________________________________
1219 Int_t AliITStrackerSA::SearchClusters(Int_t layer,Double_t phiwindow,Double_t lambdawindow, AliITStrackSA* trs,Double_t /*zvertex*/,Int_t pflag){
1220   //function used to to find the clusters associated to the track
1221   Int_t nc=0;
1222   AliITSlayer &lay = fgLayers[layer];
1223   Double_t r=lay.GetR();
1224    if(pflag==1){      
1225     Float_t cx1,cx2,cy1,cy2;
1226     FindEquation(fPoint1[0],fPoint1[1],fPoint2[0],fPoint2[1],fPoint3[0],fPoint3[1],fCoef1,fCoef2,fCoef3);
1227     if (FindIntersection(fCoef1,fCoef2,fCoef3,-r*r,cx1,cy1,cx2,cy2)==0)
1228        return 0;
1229     Double_t fi1=TMath::ATan2(cy1-fPoint1[1],cx1-fPoint1[0]);
1230     Double_t fi2=TMath::ATan2(cy2-fPoint1[1],cx2-fPoint1[0]);
1231     fPhiEstimate=ChoosePoint(fi1,fi2,fPhic);
1232   }
1233
1234  
1235   Int_t ncl = fCluLayer[layer]->GetEntries();
1236   for (Int_t index=0; index<ncl; index++) {
1237     AliITSRecPoint *c = (AliITSRecPoint*)fCluLayer[layer]->At(index);
1238     if (!c) continue;
1239     if (c->GetQ()<=0) continue;
1240     
1241      AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(layer,index);
1242      Double_t phi = arr->GetPhi();
1243      if (TMath::Abs(phi-fPhiEstimate)>phiwindow) continue;
1244
1245      Double_t lambda = arr->GetLambda();
1246      if (TMath::Abs(lambda-fLambdac)>lambdawindow) continue;
1247
1248      if(trs->GetNumberOfClustersSA()==trs->GetMaxNumberOfClusters()) return 0;
1249      if(trs->GetNumberOfMarked(layer)==trs->GetMaxNMarkedPerLayer()) return 0;
1250      Int_t orind = arr->GetOrInd();
1251      trs->AddClusterSA(layer,orind);
1252      trs->AddClusterMark(layer,index);
1253        
1254      nc++;
1255      fLambdac=lambda;
1256      fPhiEstimate=phi;
1257
1258      fPointc[0]=arr->GetX();
1259      fPointc[1]=arr->GetY();
1260
1261   }
1262   return nc;
1263 }
1264
1265 //________________________________________________________________
1266 void AliITStrackerSA::UpdatePoints(){
1267   //update of points for the estimation of the curvature  
1268
1269   fPoint2[0]=fPoint3[0];
1270   fPoint2[1]=fPoint3[1];
1271   fPoint3[0]=fPointc[0];
1272   fPoint3[1]=fPointc[1];
1273
1274   
1275 }
1276
1277 //___________________________________________________________________
1278 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){
1279
1280    //given (x,y) of three recpoints (in global coordinates) 
1281    //returns the parameters a,b,c of circonference x*x + y*y +a*x + b*y +c
1282
1283    Float_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1284    if(den==0) return 0;
1285    a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1286    b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1287    c = -x1*x1-y1*y1-a*x1-b*y1;
1288    return 1;
1289  }
1290 //__________________________________________________________________________
1291  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){
1292  
1293  //Finds the intersection between the circonference of the track and the circonference centered in (0,0) represented by one layer
1294  //c2 is -rlayer*rlayer
1295
1296   if(a1==0) return 0;
1297  Double_t m = c2-c1; 
1298  Double_t aA = (b1*b1)/(a1*a1)+1;
1299  Double_t bB = (-2*m*b1/(a1*a1));
1300  Double_t cC = c2+(m*m)/(a1*a1);
1301  Double_t dD = bB*bB-4*aA*cC;
1302  if(dD<0) return 0;
1303  
1304  y1 = (-bB+TMath::Sqrt(dD))/(2*aA); 
1305  y2 = (-bB-TMath::Sqrt(dD))/(2*aA); 
1306  x1 = (c2-c1-b1*y1)/a1;
1307  x2 = (c2-c1-b1*y2)/a1;
1308
1309  return 1; 
1310 }
1311 //____________________________________________________________________
1312 Double_t AliITStrackerSA::Curvature(Double_t x1,Double_t y1,Double_t 
1313 x2,Double_t y2,Double_t x3,Double_t y3){
1314
1315   //calculates the curvature of track  
1316   Double_t den = (x3-x1)*(y2-y1)-(x2-x1)*(y3-y1);
1317   if(den==0) return 0;
1318   Double_t a = ((y3-y1)*(x2*x2+y2*y2-x1*x1-y1*y1)-(y2-y1)*(x3*x3+y3*y3-x1*x1-y1*y1))/den;
1319   Double_t b = -(x2*x2-x1*x1+y2*y2-y1*y1+a*(x2-x1))/(y2-y1);
1320   Double_t c = -x1*x1-y1*y1-a*x1-b*y1;
1321   Double_t xc=-a/2.;
1322
1323   if((a*a+b*b-4*c)<0) return 0;
1324   Double_t rad = TMath::Sqrt(a*a+b*b-4*c)/2.;
1325   if(rad==0) return 0;
1326   
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   // if((x1>0 && y1<0 && x1<xc)) rad*=-1;
1331   
1332   return 1/rad;
1333  
1334 }
1335
1336
1337 //____________________________________________________________________
1338 Double_t AliITStrackerSA::ChoosePoint(Double_t p1, Double_t p2, Double_t pp){
1339
1340   //Returns the point closest to pp
1341
1342   Double_t diff1 = p1-pp;
1343   Double_t diff2 = p2-pp;
1344   
1345   if(TMath::Abs(diff1)<TMath::Abs(diff2)) fPhiEstimate=p1;
1346   else fPhiEstimate=p2;  
1347   return fPhiEstimate;
1348   
1349 }
1350
1351
1352 //_________________________________________________________________
1353 Int_t AliITStrackerSA::FindTrackLowChiSquare(TObjArray* tracklist, Int_t dim) const {
1354   // returns track with lowers chi square  
1355   if(dim==1){
1356     //AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(0);
1357     //return trk;
1358     return 0;
1359   }
1360   if(dim==0) return 0;
1361   Double_t * chi2 = new Double_t[dim];
1362   Int_t * index = new Int_t[dim];
1363   for(Int_t i=0;i<dim;i++){
1364     AliITStrackV2* trk = (AliITStrackV2*)tracklist->At(i);
1365     chi2[i]=trk->GetChi2();
1366     index[i]=i;
1367   }
1368
1369   Int_t w=0;Double_t value;
1370   Int_t lp;
1371   while(w<dim){
1372     for(Int_t j=w+1;j<dim;j++){
1373       if(chi2[w]<chi2[j]){
1374         value=chi2[w];
1375         chi2[w]=chi2[j];
1376         chi2[j]=value;
1377         lp=index[w];
1378         index[w]=index[j];
1379         index[j]=lp;
1380       }
1381     }
1382     w++;
1383   }
1384   Int_t tmp = index[dim-1];
1385   delete [] chi2;
1386   delete [] index;
1387   return tmp;
1388 }
1389
1390 //__________________________________________________________
1391 Int_t AliITStrackerSA::FindLabel(Int_t l1, Int_t l2, Int_t l3, Int_t l4, Int_t l5, Int_t l6){
1392
1393   //function used to determine the track label
1394   
1395   Int_t lb[6] = {l1,l2,l3,l4,l5,l6};
1396   Int_t aa[6]={1,1,1,1,1,1};
1397   Int_t ff=0; 
1398   Int_t ll=0;
1399   Int_t k=0;Int_t w=0;Int_t num=6;
1400   for(Int_t i=5;i>=0;i--) if(lb[i]==-1) num=i;
1401   
1402   while(k<num){
1403   
1404     for(Int_t i=k+1;i<num;i++){
1405     
1406       if(lb[k]==lb[i] && aa[k]!=0){
1407       
1408         aa[k]+=1;
1409         aa[i]=0;
1410       }
1411     }
1412     k++;
1413   }
1414
1415   while(w<num){
1416  
1417     for(Int_t j=0;j<6;j++){
1418       if(aa[w]<aa[j]){
1419       ff=aa[w];
1420       aa[w]=aa[j];
1421       aa[j]=ff;
1422       ll=lb[w];
1423       lb[w]=lb[j];
1424       lb[j]=ll;
1425      }
1426     }
1427     w++;
1428   }
1429   
1430   if(num<1) return -1; 
1431   return lb[num-1];
1432 }
1433
1434 //_____________________________________________________________________________
1435 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,
1436 Int_t gl12, Int_t gl13, Int_t gl14,Int_t gl15, Int_t gl16, Int_t gl17, Int_t gl18, Int_t numberofpoints){
1437
1438  
1439   //function used to assign label to the found track. If track is fake, the label is negative
1440
1441   Int_t lb0[6] = {gl1,gl2,gl3,gl4,gl5,gl6};
1442   Int_t lb1[6] = {gl7,gl8,gl9,gl10,gl11,gl12};
1443   Int_t lb2[6] = {gl13,gl14,gl15,gl16,gl17,gl18};
1444   Int_t ll=FindLabel(lb0[0],lb0[1],lb0[2],lb0[3],lb0[4],lb0[5]);
1445   Int_t lflag=0;Int_t num=6;
1446   if(lb0[5]==-1 && lb1[5]==-1 && lb2[5]==-1) num=5;
1447
1448   for(Int_t i=0;i<num;i++){
1449     if(lb0[i]==ll || lb1[i]==ll || lb2[i]==ll) lflag+=1;
1450   }
1451
1452   if(lflag>=numberofpoints) return ll;
1453   else return -ll;
1454
1455   
1456 }
1457
1458 //_____________________________________________________________________________
1459 void AliITStrackerSA::SetWindowSizes(Int_t n, Double_t *phi, Double_t *lam){
1460   // Set sizes of the phi and lambda windows used for track finding
1461   fNloop = n;
1462   if(phi){ // user defined values
1463     fPhiWin = new Double_t[fNloop];
1464     fLambdaWin = new Double_t[fNloop];
1465     for(Int_t k=0;k<fNloop;k++){
1466       fPhiWin[k]=phi[k];
1467       fLambdaWin[k]=lam[k];
1468     }
1469   }
1470   else {  // default values
1471             
1472     Double_t phid[33]   = {0.002,0.003,0.004,0.0045,0.0047,
1473                            0.005,0.0053,0.0055,
1474                            0.006,0.0063,0.0065,0.007,0.0073,0.0075,0.0077,
1475                            0.008,0.0083,0.0085,0.0087,0.009,0.0095,0.0097,
1476                            0.01,0.0105,0.011,0.0115,0.012,0.0125,0.013,0.0135,0.0140,0.0145};
1477     Double_t lambdad[33] = {0.003,0.004,0.005,0.005,0.005,
1478                             0.005,0.005,0.006,
1479                             0.006,0.006,0.006,0.007,0.007,0.007,0.007,
1480                             0.007,0.007,0.007,0.007,0.007,0.007,0.007,
1481                             0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008,0.008};
1482     
1483     if(fNloop!=33){
1484       fNloop = 33;
1485     }
1486     
1487     
1488     fPhiWin = new Double_t[fNloop];
1489     fLambdaWin = new Double_t[fNloop];
1490
1491
1492
1493     Double_t factor=AliITSReconstructor::GetRecoParam()->GetFactorSAWindowSizes(); // possibility to enlarge windows for cosmics reco with large misalignments (A.Dainese)
1494     for(Int_t k=0;k<fNloop;k++) {
1495       fPhiWin[k]=phid[k]*factor;
1496       fLambdaWin[k]=lambdad[k]*factor;
1497     }
1498   
1499   }
1500
1501 }
1502 //_______________________________________________________________________
1503 void AliITStrackerSA::GetCoorAngles(AliITSRecPoint* cl,Double_t &phi,Double_t &lambda, Float_t &x, Float_t &y,Float_t &z,Double_t* vertex){
1504   //Returns values of phi (azimuthal) and lambda angles for a given cluster
1505 /*  
1506   Double_t rot[9];     fGeom->GetRotMatrix(module,rot);
1507   Int_t lay,lad,det; fGeom->GetModuleId(module,lay,lad,det);
1508   Float_t tx,ty,tz;  fGeom->GetTrans(lay,lad,det,tx,ty,tz);     
1509
1510   Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1511   Double_t phi1=TMath::Pi()/2+alpha;
1512   if (lay==1) phi1+=TMath::Pi();
1513
1514   Float_t cp=TMath::Cos(phi1), sp=TMath::Sin(phi1);
1515   Float_t r=tx*cp+ty*sp;
1516
1517   xyz= r*cp - cl->GetY()*sp;
1518   y= r*sp + cl->GetY()*cp;
1519   z=cl->GetZ();
1520 */
1521   Float_t xyz[3];
1522   cl->GetGlobalXYZ(xyz);
1523   x=xyz[0];
1524   y=xyz[1];
1525   z=xyz[2];
1526  
1527   phi=TMath::ATan2(y-vertex[1],x-vertex[0]);
1528   lambda=TMath::ATan2(z-vertex[2],TMath::Sqrt((x-vertex[0])*(x-vertex[0])+(y-vertex[1])*(y-vertex[1])));
1529 }
1530
1531 //________________________________________________________________________
1532 void AliITStrackerSA::GetCoorErrors(AliITSRecPoint* cl,Float_t &sx,Float_t &sy, Float_t &sz){
1533
1534   //returns sigmax, y, z of cluster in global coordinates
1535 /*
1536   Double_t rot[9];     fGeom->GetRotMatrix(module,rot);
1537   Int_t lay,lad,det; 
1538   AliITSgeomTGeo::GetModuleId(module,lay,lad,det);
1539  
1540   Double_t alpha=TMath::ATan2(rot[1],rot[0])+TMath::Pi();
1541   Double_t phi=TMath::Pi()/2+alpha;
1542   if (lay==1) phi+=TMath::Pi();
1543
1544   Float_t cp=TMath::Cos(phi), sp=TMath::Sin(phi);
1545 */
1546   Float_t covm[6];
1547   cl->GetGlobalCov(covm);
1548   sx=TMath::Sqrt(covm[0]);
1549   sy=TMath::Sqrt(covm[3]);
1550   sz=TMath::Sqrt(covm[5]);
1551 /*
1552   sx = TMath::Sqrt(sp*sp*cl->GetSigmaY2());
1553   sy = TMath::Sqrt(cp*cp*cl->GetSigmaY2());
1554   sz = TMath::Sqrt(cl->GetSigmaZ2());
1555 */
1556 }
1557