ITS/AliITStrackerSA.cxx

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