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   if(fCluLayer){
 221     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 222       if(fCluLayer[i]){
 223         fCluLayer[i]->Delete();
 224         delete fCluLayer[i];
 225       }
 226     }
 227     delete [] fCluLayer;
 228   }
 229   if(fCluCoord){
 230     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 231       if(fCluCoord[i]){
 232         fCluCoord[i]->Delete();
 233         delete fCluCoord[i];
 234       }
 235     }
 236     delete [] fCluCoord;
 237   }
 238
 239 }
 240
 241 //____________________________________________________________________________
 242 Int_t AliITStrackerSA::Clusters2Tracks(AliESDEvent *event){
 243 // This method is used to find and fit the tracks. By default the corresponding
 244 // method in the parent class is invoked. In this way a combined tracking
 245 // TPC+ITS is performed. If the flag fITSStandAlone is true, the tracking
 246 // is done in the ITS only. In the standard reconstruction chain this option
 247 // can be set via AliReconstruction::SetOption("ITS","onlyITS")
 248   Int_t rc=0;
 249   if(!fITSStandAlone){
 250     rc=AliITStrackerMI::Clusters2Tracks(event);
 251   }
 252   else {
 253     AliDebug(1,"Stand Alone flag set: doing tracking in ITS alone\n");
 254   }
 255   if(!rc) rc=FindTracks(event);
 256   return rc;
 257 }
 258
 259 //____________________________________________________________________________
 260 void AliITStrackerSA::Init(){
 261   //  Reset all data members
 262     fPhiEstimate=0;
 263     for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
 264     fLambdac=0;
 265     fPhic=0;
 266     fCoef1=0;
 267     fCoef2=0;
 268     fCoef3=0;
 269     fPointc[0]=0;
 270     fPointc[1]=0;
 271     fVert = 0;
 272     fVertexer = 0;
 273     SetWindowSizes();
 274     fITSclusters = 0;
 275     SetSixPoints();
 276     SetOuterStartLayer(0);
 277     SetSAFlag(kFALSE);
 278     fListOfTracks=new TObjArray(0,0);
 279     fCluLayer = 0;
 280     fCluCoord = 0;
 281  }
 282 //_______________________________________________________________________
 283 void AliITStrackerSA::ResetForFinding(){
 284   //  Reset data members used in all loops during track finding
 285     fPhiEstimate=0;
 286     for(Int_t i=0;i<3;i++){fPoint1[i]=0;fPoint2[i]=0;fPoint3[i]=0;}
 287     fLambdac=0;
 288     fPhic=0;
 289     fCoef1=0;
 290     fCoef2=0;
 291     fCoef3=0;
 292     fPointc[0]=0;
 293     fPointc[1]=0;
 294     fListOfTracks->Delete();
 295 }
 296
 297
 298
 299 //______________________________________________________________________
 300 Int_t AliITStrackerSA::FindTracks(AliESDEvent* event){
 301
 302 // Track finder using the ESD object
 303
 304
 305   //controllare numero cluster sui layer1 e 2 (morti?)
 306   //non trova tracce...controllare..
 307
 308   if(!fITSclusters){
 309     Fatal("FindTracks","ITS cluster tree is not accessed - Abort!!!\n Please use method SetClusterTree to pass the pointer to the tree\n");
 310     return -1;
 311   }
 312
 313
 314   //Reads event and mark clusters of traks already found, with flag kITSin
 315    Int_t nentr=event->GetNumberOfTracks();
 316    while (nentr--) {
 317      AliESDtrack *track=event->GetTrack(nentr);
 318      if (track->GetStatus()&AliESDtrack::kITSin==AliESDtrack::kITSin){
 319        Int_t idx[12];
 320        Int_t ncl = track->GetITSclusters(idx);
 321        for(Int_t k=0;k<ncl;k++){
 322          AliITSRecPoint* cll = (AliITSRecPoint*)GetCluster(idx[k]);
 323          cll->SetBit(kSAflag);
 324        }
 325      }
 326    }
 327
 328    //Get primary vertex
 329    Double_t primaryVertex[3];
 330    event->GetVertex()->GetXYZ(primaryVertex);
 331    //Creates TClonesArray with clusters for each layer. The clusters already used
 332    //by AliITStrackerMI are not considered
 333
 334    Int_t nclusters[6]={0,0,0,0,0,0};
 335    Int_t dmar[6]={0,0,0,0,0,0};
 336    fCluLayer = new TClonesArray*[AliITSgeomTGeo::GetNLayers()];
 337    fCluCoord = new TClonesArray*[AliITSgeomTGeo::GetNLayers()];
 338
 339    for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 340      AliITSlayer &layer=fgLayers[i];
 341      for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
 342        AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
 343        if(cls->TestBit(kSAflag)==kTRUE) continue; //clusters used by TPC prol.
 344        if(cls->GetQ()==0) continue; //fake clusters dead zones
 345        nclusters[i]++;
 346      }
 347      dmar[i]=0;
 348      fCluLayer[i] = new TClonesArray("AliITSRecPoint",nclusters[i]);
 349      fCluCoord[i] = new TClonesArray("AliITSclusterTable",nclusters[i]);
 350    }
 351
 352
 353    Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
 354    for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 355      firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
 356    }
 357
 358    for(Int_t ilay=0;ilay<AliITSgeomTGeo::GetNLayers();ilay++){
 359      TClonesArray &clulay = *fCluLayer[ilay];
 360      TClonesArray &clucoo = *fCluCoord[ilay];
 361      AliITSlayer &layer=fgLayers[ilay];
 362      for(Int_t cli=0;cli<layer.GetNumberOfClusters();cli++){
 363        AliITSRecPoint* cls = (AliITSRecPoint*)layer.GetCluster(cli);
 364        if(cls->TestBit(kSAflag)==kTRUE) continue;
 365        if(cls->GetQ()==0) continue;
 366        Double_t phi=0;Double_t lambda=0;
 367        Float_t x=0;Float_t y=0;Float_t z=0;
 368        Float_t sx=0;Float_t sy=0;Float_t sz=0;
 369        GetCoorAngles(cls,phi,lambda,x,y,z,primaryVertex);
 370        GetCoorErrors(cls,sx,sy,sz);
 371        new (clulay[dmar[ilay]]) AliITSRecPoint(*cls);
 372        new (clucoo[dmar[ilay]]) AliITSclusterTable(x,y,z,sx,sy,sz,phi,lambda,cli);
 373        dmar[ilay]++;
 374      }
 375    }
 376
 377
 378
 379    Int_t ntrack=0;
 380    //loop on the different windows
 381    for(Int_t nloop=0;nloop<fNloop;nloop++){
 382      for(Int_t ncl=0;ncl<fCluLayer[0]->GetEntries();ncl++){ //loop starting from layer 0
 383
 384        ResetForFinding();
 385        Int_t pflag=0;
 386
 387        AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[0]->At(ncl);
 388
 389        if(!cl) continue;
 390
 391
 392        if (cl->GetQ()<=0) continue;
 393
 394        AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(0,ncl);
 395        fPhic = arr->GetPhi();
 396        fLambdac = arr->GetLambda();
 397        if (TMath::Abs(fLambdac)>0.26*TMath::Pi()) continue;
 398        fPhiEstimate = fPhic;
 399        AliITStrackSA* trs = new AliITStrackSA();
 400        fPoint1[0]=primaryVertex[0];
 401        fPoint1[1]=primaryVertex[1];
 402
 403
 404        fPoint2[0]=arr->GetX();
 405        fPoint2[1]=arr->GetY();
 406        Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()];//counter for clusters on each layer
 407        for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){ nn[i]=0;}
 408        nn[0] = SearchClusters(0,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
 409
 410        nn[1] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
 411       if(nn[1]>0){
 412         pflag=1;
 413         fPoint3[0] = fPointc[0];
 414         fPoint3[1] = fPointc[1];
 415       }
 416       nn[2] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
 417       if(nn[1]==0 && nn[2]==0) pflag=0;
 418       if(nn[2]!=0 && nn[1]!=0){ pflag=1; UpdatePoints();}
 419       if(nn[2]!=0 && nn[1]==0){
 420         pflag=1;
 421         fPoint3[0]=fPointc[0];
 422         fPoint3[1]=fPointc[1];
 423       }
 424
 425       nn[3] = SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
 426       pflag=1;
 427       if(nn[3]!=0) UpdatePoints();
 428       nn[4] = SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
 429       pflag=1;
 430       if(nn[4]!=0) UpdatePoints();
 431       nn[5] = SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],trs,primaryVertex[2],pflag);
 432
 433
 434       Int_t layOK=0;
 435       Int_t numberofpoints;
 436       if(fSixPoints) numberofpoints=6;  //check of the candidate track
 437       else numberofpoints=5;           //if track is good (with the required number
 438       for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers();nnp++){    //of points) it is written on file
 439         if(nn[nnp]!=0) layOK+=1;
 440       }
 441       if(layOK>=numberofpoints){
 442
 443         AliITStrackV2* tr2 = FitTrack(trs,primaryVertex);
 444
 445         if(tr2==0) continue;
 446
 447         AliESDtrack outtrack;
 448         outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
 449         event->AddTrack(&outtrack);
 450         ntrack++;
 451       }
 452
 453
 454       delete trs;
 455       delete[] nn;
 456
 457      }//end loop on clusters of layer1
 458
 459      //end loop2
 460    }
 461
 462    //if 5/6 points are required, second loop starting
 463    //from second layer (SPD2), to find tracks with point of
 464    //layer 1 missing
 465    if(!fSixPoints) {
 466      //printf("looking from SPD2\n");
 467      //   counter for clusters on each layer
 468      Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()-1];
 469      for(Int_t nloop=0;nloop<fNloop;nloop++){
 470        Int_t ncl2=fCluLayer[1]->GetEntries();
 471        while(ncl2--){ //loop starting from layer 2
 472          ResetForFinding();
 473          Int_t pflag=0;
 474          AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[1]->At(ncl2);
 475
 476          if(!cl) continue;
 477          AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(1,ncl2);
 478          fPhic = arr->GetPhi();
 479          fLambdac = arr->GetLambda();
 480          fPhiEstimate = fPhic;
 481
 482          AliITStrackSA* trs = new AliITStrackSA();
 483          fPoint1[0]=primaryVertex[0];
 484          fPoint1[1]=primaryVertex[1];
 485
 486          fPoint2[0]=arr->GetX();
 487          fPoint2[1]=arr->GetY();
 488          for(Int_t kk=0;kk<AliITSgeomTGeo::GetNLayers()-1;kk++)nn[kk] = 0;
 489          nn[0] = SearchClusters(1,fPhiWin[nloop],fLambdaWin[nloop],
 490                                    trs,primaryVertex[2],pflag);
 491          nn[1] = SearchClusters(2,fPhiWin[nloop],fLambdaWin[nloop],
 492                                    trs,primaryVertex[2],pflag);
 493          if(nn[1]!=0){
 494            pflag=1;
 495            fPoint3[0]=fPointc[0];
 496            fPoint3[1]=fPointc[1];
 497          }
 498          nn[2]= SearchClusters(3,fPhiWin[nloop],fLambdaWin[nloop],
 499                                    trs,primaryVertex[2],pflag);
 500          if(nn[2]!=0){
 501            pflag=1;
 502            UpdatePoints();
 503          }
 504          nn[3]= SearchClusters(4,fPhiWin[nloop],fLambdaWin[nloop],
 505                                    trs,primaryVertex[2],pflag);
 506          if(nn[3]!=0){
 507            pflag=1;
 508            UpdatePoints();
 509          }
 510          nn[4]=SearchClusters(5,fPhiWin[nloop],fLambdaWin[nloop],
 511                                    trs,primaryVertex[2],pflag);
 512
 513          Int_t fl=0;
 514          for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-1;nnp++){
 515            if(nn[nnp]!=0) fl+=1;
 516          }
 517          if(fl>=5){  // 5/6
 518            AliITStrackV2* tr2 = FitTrack(trs,primaryVertex);
 519            if(tr2==0){
 520              continue;
 521            }
 522
 523            AliESDtrack outtrack;
 524            outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
 525            event->AddTrack(&outtrack);
 526            ntrack++;
 527
 528          }
 529
 530         delete trs;
 531       }//end loop on clusters of layer2
 532      }
 533
 534     delete [] nn;
 535   } //end opt="5/6"
 536
 537
 538
 539    // search for tracks starting from SPD2, SDD1, SDD2, SSD2
 540    // for cosmics (A. Dainese 31.07.07)
 541    if(fOuterStartLayer>0) {
 542      for(Int_t innLay=1; innLay<=fOuterStartLayer; innLay++) {
 543        printf("Searching from layer %d outward\n",innLay);
 544        //   counter for clusters on each layer
 545        Int_t * nn = new Int_t[AliITSgeomTGeo::GetNLayers()-innLay];
 546        for(Int_t nloop=0;nloop<fNloop;nloop++){
 547          Int_t nclInnLay=fCluLayer[innLay]->GetEntries();
 548          while(nclInnLay--){ //loop starting from layer innLay
 549            ResetForFinding();
 550            Int_t pflag=0;
 551            AliITSRecPoint* cl = (AliITSRecPoint*)fCluLayer[innLay]->At(nclInnLay);
 552
 553            if(!cl) continue;
 554            AliITSclusterTable* arr = (AliITSclusterTable*)GetClusterCoord(innLay,nclInnLay);
 555            fPhic = arr->GetPhi();
 556            fLambdac = arr->GetLambda();
 557            fPhiEstimate = fPhic;
 558
 559            AliITStrackSA* trs = new AliITStrackSA();
 560            fPoint1[0]=primaryVertex[0];
 561            fPoint1[1]=primaryVertex[1];
 562            fPoint2[0]=arr->GetX();
 563            fPoint2[1]=arr->GetY();
 564
 565            Int_t kk;
 566            for(kk=0;kk<AliITSgeomTGeo::GetNLayers()-innLay;kk++) nn[kk] = 0;
 567
 568            kk=0;
 569            nn[kk] = SearchClusters(innLay,fPhiWin[nloop],fLambdaWin[nloop],
 570                                   trs,primaryVertex[2],pflag);
 571            for(Int_t nextLay=innLay+1; nextLay<AliITSgeomTGeo::GetNLayers(); nextLay++) {
 572              kk++;
 573              nn[kk] = SearchClusters(nextLay,fPhiWin[nloop],fLambdaWin[nloop],
 574                                      trs,primaryVertex[2],pflag);
 575              if(nn[kk]!=0){
 576                pflag=1;
 577                if(kk==1) {
 578                  fPoint3[0]=fPointc[0];
 579                  fPoint3[1]=fPointc[1];
 580                } else {
 581                  UpdatePoints();
 582                }
 583              }
 584            }
 585
 586            Int_t fl=0;
 587            for(Int_t nnp=0;nnp<AliITSgeomTGeo::GetNLayers()-innLay;nnp++){
 588              if(nn[nnp]!=0) fl+=1;
 589            }
 590            if(fl>=AliITSgeomTGeo::GetNLayers()-innLay){
 591              AliITStrackV2* tr2 = FitShortTrack(trs,primaryVertex,innLay);
 592              if(tr2==0){
 593                continue;
 594              }
 595
 596              AliESDtrack outtrack;
 597              outtrack.UpdateTrackParams(tr2,AliESDtrack::kITSin);
 598              event->AddTrack(&outtrack);
 599              ntrack++;
 600
 601            }
 602
 603            delete trs;
 604          }//end loop on clusters of innLay
 605        } //end loop on window sizes
 606
 607        delete [] nn;
 608      } //end loop on innLay
 609    } //end if(fOuterStartLayer>0)
 610
 611
 612   delete [] firstmod;
 613   Info("FindTracks","Number of found tracks: %d",event->GetNumberOfTracks());
 614   return 0;
 615
 616 }
 617
 618
 619
 620
 621
 622 //________________________________________________________________________
 623
 624 AliITStrackV2* AliITStrackerSA::FitTrack(AliITStrackSA* tr,Double_t *primaryVertex){
 625   //fit of the found track
 626
 627
 628   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
 629   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 630     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
 631   }
 632
 633   Int_t nclusters = tr->GetNumberOfClustersSA();
 634   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
 635   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 636     listlayer[i] = new TObjArray(0,0);
 637   }
 638
 639   TArrayI clind0(20);
 640   TArrayI clind1(20);
 641   TArrayI clind2(20);
 642   TArrayI clind3(20);
 643   TArrayI clind4(20);
 644   TArrayI clind5(20);
 645
 646   TArrayI mark0(20);
 647   TArrayI mark1(20);
 648   TArrayI mark2(20);
 649   TArrayI mark3(20);
 650   TArrayI mark4(20);
 651   TArrayI mark5(20);
 652
 653
 654   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
 655   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
 656   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
 657
 658   for(Int_t ncl=0;ncl<nclusters;ncl++){
 659     Int_t index = tr->GetClusterIndexSA(ncl);
 660     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
 661     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
 662     Int_t lay = (index & 0xf0000000) >> 28;
 663     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
 664     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
 665     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
 666     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
 667     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
 668     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}
 669   }
 670   delete [] nnn;
 671
 672   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
 673     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
 674       Int_t mark = tr->GetClusterMark(nlay,ncl);
 675       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
 676       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
 677       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
 678       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
 679       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
 680       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
 681
 682     }
 683   }
 684
 685   delete [] kkk;
 686
 687
 688   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
 689   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 690     if(listlayer[i]->GetEntries()==0) end[i]=1;
 691     else end[i]=listlayer[i]->GetEntries();
 692   }
 693
 694   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
 695   TClonesArray &tri = *listSA;
 696   Int_t nlist=0;
 697
 698   if(end[0]==0) end[0]=1; //for tracks with cluster on layer 0 missing
 699   for(Int_t l1=0;l1<end[0];l1++){ //loop on layer 1
 700     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l1);
 701     Double_t x1,y1,z1,sx1,sy1,sz1;
 702     Double_t x2,y2,z2,sx2,sy2,sz2;
 703     AliITSRecPoint* p1=0;
 704     AliITSRecPoint* p2=0;
 705     Int_t index1=clind0[l1];
 706     Int_t index2=0;
 707     Int_t mrk1 = mark0[l1];
 708     Int_t mrk2 = 0;
 709     Int_t lay1=0;
 710     Int_t lay2=1;
 711     Int_t module1=-1;
 712     for(Int_t l2=0;l2<end[1];l2++){ //loop on layer 2
 713       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l2);
 714       index2=clind1[l2];
 715       mrk2 = mark1[l2];
 716       for(Int_t l3=0;l3<end[2];l3++){  //loop on layer 3
 717         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l3);
 718
 719         if(cl0==0 && cl1!=0) {
 720           p1=cl1;
 721           p2=cl2;
 722           index1=clind1[l2];mrk1=mark1[l2];
 723           index2=clind2[l3];mrk2=mark2[l3];
 724           lay1=1;
 725           lay2=2;
 726           module1 = p1->GetDetectorIndex()+firstmod[1];
 727         }
 728         if(cl0!=0 && cl1==0){
 729           p1=cl0;
 730           p2=cl2;index2=clind2[l3];mrk2=mark2[l3];
 731           lay1=0;
 732           lay2=2;
 733           module1 = p1->GetDetectorIndex()+firstmod[0];
 734         }
 735         if(cl0!=0 && cl1!=0){
 736           p1=cl0;
 737           p2=cl1;
 738           lay1=0;
 739           lay2=1;
 740           module1 = p1->GetDetectorIndex()+firstmod[0];
 741         }
 742
 743         Int_t cln1=mrk1;
 744         Int_t cln2=mrk2;
 745         AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(lay1,cln1);
 746         AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(lay2,cln2);
 747         x1 = arr1->GetX();
 748         x2 = arr2->GetX();
 749         y1 = arr1->GetY();
 750         y2 = arr2->GetY();
 751         z1 = arr1->GetZ();
 752         z2 = arr2->GetZ();
 753         sx1 = arr1->GetSx();
 754         sx2 = arr2->GetSx();
 755         sy1 = arr1->GetSy();
 756         sy2 = arr2->GetSy();
 757         sz1 = arr1->GetSz();
 758         sz2 = arr2->GetSz();
 759
 760         Int_t layer,ladder,detector;
 761         AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
 762         Float_t yclu1 = p1->GetY();
 763         Float_t zclu1 = p1->GetZ();
 764         Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
 765         Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
 766         Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
 767
 768         for(Int_t l4=0;l4<end[3];l4++){ //loop on layer 4
 769           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l4);
 770           for(Int_t l5=0;l5<end[4];l5++){ //loop on layer 5
 771             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l5);
 772             for(Int_t l6=0;l6<end[5];l6++){ //loop on layer 6
 773               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l6);
 774               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
 775
 776               if(cl5!=0) {
 777                 trac->AddClusterV2(5,(clind5[l6] & 0x0fffffff)>>0);
 778                 trac->AddClusterMark(5,mark5[l6]);
 779               }
 780               if(cl4!=0){
 781                 trac->AddClusterV2(4,(clind4[l5] & 0x0fffffff)>>0);
 782                 trac->AddClusterMark(4,mark4[l5]);
 783               }
 784               if(cl3!=0){
 785                 trac->AddClusterV2(3,(clind3[l4] & 0x0fffffff)>>0);
 786                 trac->AddClusterMark(3,mark3[l4]);
 787               }
 788               if(cl2!=0){
 789                 trac->AddClusterV2(2,(clind2[l3] & 0x0fffffff)>>0);
 790                 trac->AddClusterMark(2,mark2[l3]);
 791               }
 792               if(cl1!=0){
 793                 trac->AddClusterV2(1,(clind1[l2] & 0x0fffffff)>>0);
 794                 trac->AddClusterMark(1,mark1[l2]);
 795               }
 796               if(cl0!=0){
 797                 trac->AddClusterV2(0,(clind0[l1] & 0x0fffffff)>>0);
 798                 trac->AddClusterMark(0,mark0[l1]);
 799               }
 800               //fit with Kalman filter using AliITStrackerMI::RefitAt()
 801
 802
 803               AliITStrackMI* ot = new AliITStrackSA(*trac);
 804
 805               ot->ResetCovariance(10.);
 806               ot->ResetClusters();
 807
 808               if(RefitAt(krInsideITSscreen,ot,trac)){ //fit from layer 1 to layer 6
 809                 AliITStrackMI *otrack2 = new AliITStrackMI(*ot);
 810                 otrack2->ResetCovariance(10.);
 811                 otrack2->ResetClusters();
 812                 //fit from layer 6 to layer 1
 813                 if(RefitAt(krInsideSPD1,otrack2,ot)) {
 814                   fListOfTracks->AddLast(otrack2);
 815                   new (tri[nlist]) AliITStrackSA(*trac);
 816                   nlist++;
 817                 } else {
 818                   delete otrack2;
 819                 }
 820
 821               }
 822
 823               delete ot;
 824               delete trac;
 825             }//end loop layer 6
 826           }//end loop layer 5
 827         }//end loop layer 4
 828       }//end loop layer 3
 829     }//end loop layer 2
 830   }//end loop layer 1
 831
 832   delete [] end;
 833
 834
 835
 836   Int_t dim=fListOfTracks->GetEntries();
 837   if(dim==0){
 838     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 839       delete listlayer[i];
 840     }
 841     delete [] listlayer;
 842     listSA->Delete();
 843     delete listSA;
 844     delete [] firstmod;
 845     return 0;
 846   }
 847
 848   Int_t lowchi2 = FindTrackLowChiSquare(fListOfTracks,dim);
 849   AliITStrackV2* otrack =(AliITStrackV2*)fListOfTracks->At(lowchi2);
 850   AliITStrackSA* trsa = (AliITStrackSA*)listSA->At(lowchi2);
 851
 852   if(otrack==0) {
 853     for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 854       delete listlayer[i];
 855     }
 856     delete [] listlayer;
 857     listSA->Delete();
 858     delete listSA;
 859     delete [] firstmod;
 860     return 0;
 861   }
 862   Int_t * indexc = new Int_t[AliITSgeomTGeo::GetNLayers()];
 863   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) indexc[i]=0;
 864   for(Int_t nind=0;nind<otrack->GetNumberOfClusters();nind++){
 865     indexc[nind] = otrack->GetClusterIndex(nind);
 866   }
 867   AliITSRecPoint* cl0 = (AliITSRecPoint*)GetCluster(indexc[0]);
 868   AliITSRecPoint* cl1 = (AliITSRecPoint*)GetCluster(indexc[1]);
 869   AliITSRecPoint* cl2 = (AliITSRecPoint*)GetCluster(indexc[2]);
 870   AliITSRecPoint* cl3 = (AliITSRecPoint*)GetCluster(indexc[3]);
 871   AliITSRecPoint* cl4 = (AliITSRecPoint*)GetCluster(indexc[4]);
 872   Int_t labl[3]={-1,-1,-1};
 873   if(otrack->GetNumberOfClusters()==AliITSgeomTGeo::GetNLayers()){
 874     AliITSRecPoint* cl5 = (AliITSRecPoint*)GetCluster(indexc[5]);
 875     labl[0]=cl5->GetLabel(0);
 876     labl[1]=cl5->GetLabel(1);
 877     labl[2]=cl5->GetLabel(2);
 878   }
 879   delete [] indexc;
 880   if(otrack->GetNumberOfClusters()==(AliITSgeomTGeo::GetNLayers()-1)){
 881     labl[0]=-1;
 882     labl[1]=-1;
 883     labl[2]=-1;
 884   }
 885   Int_t numberofpoints;
 886   if(fSixPoints) numberofpoints=6;
 887   else numberofpoints=5;
 888   CookLabel(otrack,0.); //MI change - to see fake ratio
 889   Int_t label =  Label(cl0->GetLabel(0),cl1->GetLabel(0),
 890                        cl2->GetLabel(0),cl3->GetLabel(0),
 891                        cl4->GetLabel(0),labl[0],
 892                        cl0->GetLabel(1),cl1->GetLabel(1),
 893                        cl2->GetLabel(1),cl3->GetLabel(1),
 894                        cl4->GetLabel(1),labl[1],
 895                        cl0->GetLabel(2),cl1->GetLabel(2),
 896                        cl2->GetLabel(2),cl3->GetLabel(2),
 897                        cl4->GetLabel(2),labl[2],numberofpoints);
 898
 899   otrack->SetLabel(label);
 900   //remove clusters of found track
 901
 902   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
 903     for(Int_t cln=0;cln<trsa->GetNumberOfMarked(nlay);cln++){
 904       Int_t index = trsa->GetClusterMark(nlay,cln);
 905       fCluLayer[nlay]->RemoveAt(index);
 906       RemoveClusterCoord(nlay,index);
 907       fCluLayer[nlay]->Compress();
 908     }
 909   }
 910   listSA->Delete();
 911   delete listSA;
 912
 913   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 914     delete listlayer[i];
 915   }
 916   delete [] listlayer;
 917   delete [] firstmod;
 918   return otrack;
 919
 920 }
 921
 922 //________________________________________________________________________
 923
 924 AliITStrackV2* AliITStrackerSA::FitShortTrack(AliITStrackSA* tr,Double_t *primaryVertex,Int_t innLay){
 925   //fit of the found track (short tracks, <6 points, for cosmics). A.Dainese 31.07.07
 926
 927
 928   Int_t * firstmod = new Int_t[AliITSgeomTGeo::GetNLayers()];
 929   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 930     firstmod[i]=AliITSgeomTGeo::GetModuleIndex(i+1,1,1);
 931   }
 932
 933   Int_t nclusters = tr->GetNumberOfClustersSA();
 934   TObjArray** listlayer = new TObjArray*[AliITSgeomTGeo::GetNLayers()];
 935   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 936     listlayer[i] = new TObjArray(0,0);
 937   }
 938
 939   TArrayI clind0(20);
 940   TArrayI clind1(20);
 941   TArrayI clind2(20);
 942   TArrayI clind3(20);
 943   TArrayI clind4(20);
 944   TArrayI clind5(20);
 945
 946   TArrayI mark0(20);
 947   TArrayI mark1(20);
 948   TArrayI mark2(20);
 949   TArrayI mark3(20);
 950   TArrayI mark4(20);
 951   TArrayI mark5(20);
 952
 953
 954   Int_t * nnn = new Int_t[AliITSgeomTGeo::GetNLayers()];
 955   Int_t * kkk = new Int_t[AliITSgeomTGeo::GetNLayers()];
 956   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++) {nnn[i]=0;kkk[i]=0;}
 957
 958   for(Int_t ncl=0;ncl<nclusters;ncl++){
 959     Int_t index = tr->GetClusterIndexSA(ncl);
 960     AliITSRecPoint* cl = (AliITSRecPoint*)GetCluster(index);
 961     if(cl->TestBit(kSAflag)==kTRUE) cl->ResetBit(kSAflag);
 962     Int_t lay = (index & 0xf0000000) >> 28;
 963     if(lay==0) { listlayer[0]->AddLast(cl); clind0[nnn[0]]=index;nnn[0]++;}
 964     if(lay==1) { listlayer[1]->AddLast(cl); clind1[nnn[1]]=index;nnn[1]++;}
 965     if(lay==2) { listlayer[2]->AddLast(cl); clind2[nnn[2]]=index;nnn[2]++;}
 966     if(lay==3) { listlayer[3]->AddLast(cl); clind3[nnn[3]]=index;nnn[3]++;}
 967     if(lay==4) { listlayer[4]->AddLast(cl); clind4[nnn[4]]=index;nnn[4]++;}
 968     if(lay==5) { listlayer[5]->AddLast(cl); clind5[nnn[5]]=index;nnn[5]++;}
 969   }
 970   delete [] nnn;
 971
 972   for(Int_t nlay=0;nlay<AliITSgeomTGeo::GetNLayers();nlay++){
 973     for(Int_t ncl=0;ncl<tr->GetNumberOfMarked(nlay);ncl++){
 974       Int_t mark = tr->GetClusterMark(nlay,ncl);
 975       if(nlay==0) { mark0[kkk[0]]=mark;kkk[0]++;}
 976       if(nlay==1) { mark1[kkk[1]]=mark;kkk[1]++;}
 977       if(nlay==2) { mark2[kkk[2]]=mark;kkk[2]++;}
 978       if(nlay==3) { mark3[kkk[3]]=mark;kkk[3]++;}
 979       if(nlay==4) { mark4[kkk[4]]=mark;kkk[4]++;}
 980       if(nlay==5) { mark5[kkk[5]]=mark;kkk[5]++;}
 981
 982     }
 983   }
 984
 985   delete [] kkk;
 986
 987
 988   Int_t * end = new Int_t[AliITSgeomTGeo::GetNLayers()];
 989   for(Int_t i=0;i<AliITSgeomTGeo::GetNLayers();i++){
 990     if(listlayer[i]->GetEntries()==0) end[i]=1;
 991     else end[i]=listlayer[i]->GetEntries();
 992   }
 993
 994   TClonesArray* listSA = new TClonesArray("AliITStrackSA");
 995   TClonesArray &tri = *listSA;
 996   Int_t nlist=0;
 997
 998   //  if(end[0]==0) end[0]=1; //for tracks with cluster on layer 0 missing
 999   for(Int_t l1=0;l1<end[0];l1++){ //loop on layer 1
1000     AliITSRecPoint* cl0 = (AliITSRecPoint*)listlayer[0]->At(l1);
1001     for(Int_t l2=0;l2<end[1];l2++){ //loop on layer 2
1002       AliITSRecPoint* cl1 = (AliITSRecPoint*)listlayer[1]->At(l2);
1003       for(Int_t l3=0;l3<end[2];l3++){  //loop on layer 3
1004         AliITSRecPoint* cl2 = (AliITSRecPoint*)listlayer[2]->At(l3);
1005         for(Int_t l4=0;l4<end[3];l4++){ //loop on layer 4
1006           AliITSRecPoint* cl3 = (AliITSRecPoint*)listlayer[3]->At(l4);
1007           for(Int_t l5=0;l5<end[4];l5++){ //loop on layer 5
1008             AliITSRecPoint* cl4 = (AliITSRecPoint*)listlayer[4]->At(l5);
1009             for(Int_t l6=0;l6<end[5];l6++){ //loop on layer 6
1010               AliITSRecPoint* cl5 = (AliITSRecPoint*)listlayer[5]->At(l6);
1011
1012
1013               Double_t x1,y1,z1,sx1,sy1,sz1;
1014               Double_t x2,y2,z2,sx2,sy2,sz2;
1015               AliITSRecPoint* p1=0;
1016               AliITSRecPoint* p2=0;
1017               Int_t index1=0,index2=0;
1018               Int_t mrk1=0,mrk2=0;
1019               if(innLay==1) {
1020                 p1=cl1;
1021                 p2=cl2;
1022                 index1=clind1[l2];mrk1=mark1[l2];
1023                 index2=clind2[l3];mrk2=mark2[l3];
1024               } else if(innLay==2) {
1025                 p1=cl2;
1026                 p2=cl3;
1027                 index1=clind2[l3];mrk1=mark2[l3];
1028                 index2=clind3[l4];mrk2=mark3[l4];
1029               } else if(innLay==3) {
1030                 p1=cl3;
1031                 p2=cl4;
1032                 index1=clind3[l4];mrk1=mark3[l4];
1033                 index2=clind4[l5];mrk2=mark4[l5];
1034               } else if(innLay==4) {
1035                 p1=cl4;
1036                 p2=cl5;
1037                 index1=clind4[l5];mrk1=mark4[l5];
1038                 index2=clind5[l6];mrk2=mark5[l6];
1039               }
1040               Int_t lay1=innLay;
1041               Int_t lay2=innLay+1;
1042               Int_t module1 = p1->GetDetectorIndex()+firstmod[innLay];
1043
1044               Int_t cln1=mrk1;
1045               Int_t cln2=mrk2;
1046               AliITSclusterTable* arr1 = (AliITSclusterTable*)GetClusterCoord(lay1,cln1);
1047               AliITSclusterTable* arr2 = (AliITSclusterTable*)GetClusterCoord(lay2,cln2);
1048               x1 = arr1->GetX();
1049               x2 = arr2->GetX();
1050               y1 = arr1->GetY();
1051               y2 = arr2->GetY();
1052               z1 = arr1->GetZ();
1053               z2 = arr2->GetZ();
1054               sx1 = arr1->GetSx();
1055               sx2 = arr2->GetSx();
1056               sy1 = arr1->GetSy();
1057               sy2 = arr2->GetSy();
1058               sz1 = arr1->GetSz();
1059               sz2 = arr2->GetSz();
1060
1061               Int_t layer,ladder,detector;
1062               AliITSgeomTGeo::GetModuleId(module1,layer,ladder,detector);
1063               Float_t yclu1 = p1->GetY();
1064               Float_t zclu1 = p1->GetZ();
1065               Double_t cv=Curvature(primaryVertex[0],primaryVertex[1],x1,y1,x2,y2);
1066               Double_t tgl2 = (z2-z1)/TMath::Sqrt((x2-x1)*(x2-x1)+(y2-y1)*(y2-y1));
1067               Double_t phi2 = TMath::ATan2((y2-y1),(x2-x1));
1068               //printf("phi2: %f   tgl2: %f\n",phi2,tgl2);
1069               AliITStrackSA* trac = new AliITStrackSA(layer,ladder,detector,yclu1,zclu1,phi2,tgl2,cv,1);
1070
1071               //Float_t clxyz[3];
1072               if(cl5!=0) {
1073                 trac->AddClusterV2(5,(clind5[l6] & 0x0fffffff)>>0);
1074                 trac->AddClusterMark(5,mark5[l6]);
1075                 //cl5->GetGlobalXYZ(clxyz);printf("gcl %f %f %f\n",clxyz[0],clxyz[1],clxyz[2]);
1076               }
1077               if(cl4!=0){
1078                 trac->AddClusterV2(4,(clind4[l5] & 0x0fffffff)>>0);
1079                 trac->AddClusterMark(4,mark4[l5]);
1080                 //cl4->GetGlobalXYZ(clxyz);printf("gcl %f %f %f\n",clxyz[0],clxyz[1],clxyz[2]);
1081               }
1082               if(cl3!=0){
1083                 trac->AddClusterV2(3,(clind3[l4] & 0x0fffffff)>>0);
1084                 trac->AddClusterMark(3,mark3[l4]);
1085                 //cl3->GetGlobalXYZ(clxyz);printf("gcl %f %f %f\n",clxyz[0],clxyz[1],clxyz[2]);
1086               }
1087               if(cl2!=0){
1088                 trac->AddClusterV2(2,(clind2[l3] & 0x0fffffff)>>0);
1089                 trac->AddClusterMark(2,mark2[l3]);
1090                 //cl2->GetGlobalXYZ(clxyz);printf("gcl %f %f %f\n",clxyz[0],clxyz[1],clxyz[2]);
1091               }
1092               if(cl1!=0){
1093                 trac->AddClusterV2(1,(clind1[l2] & 0x0fffffff)>>0);
1094                 trac->AddClusterMark(1,mark1[l2]);
1095                 //cl1->GetGlobalXYZ(clxyz);printf("gcl %f %f %f\n",clxyz[0],clxyz[1],clxyz[2]);
1096               }
1097               if(cl0!=0){
1098                 trac->AddClusterV2(0,(clind0[l1] & 0x0fffffff)>>0);
1099                 trac->AddClusterMark(0,mark0[l1]);
1100                 //cl0->GetGlobalXYZ(clxyz);printf("gcl %f %f %f\n",clxyz[0],clxyz[1],clxyz[2]);
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               //printf("%d  \n",trac->GetNumberOfClusters());for(Int_t ijk=0;ijk<trac->GetNumberOfClusters();ijk++) printf("    %d",trac->GetClusterIndex(ijk)); printf("\n");
1111               if(RefitAt(krInsideITSscreen,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(krInsideSPD1,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 lowes 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
1385   delete [] chi2;
1386   delete [] index;
1387   return index[dim-1];
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   if(lb[5]==-1) num=5;
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   if(num==6)  return lb[5];
1430   else return lb[4];
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
1554