HLT/ITS/tracking/AliITStrackerHLT.cxx

   1 /**************************************************************************
   2  * Copyright(c) 2007-2009, 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 /* $Id: AliITStrackerHLT.cxx 32466 2009-05-20 07:51:56Z hristov $ */
  16
  17 //-------------------------------------------------------------------------
  18 //               Implementation of the ITS tracker class
  19 //    It reads AliITSRecPoint clusters and creates AliHLTITSTrack tracks
  20 //                   and fills with them the ESD
  21 //          Origin: Marian Ivanov, CERN, Marian.Ivanov@cern.ch
  22 //          Current support and development:
  23 //                     Andrea Dainese, andrea.dainese@lnl.infn.it
  24 //     dE/dx analysis by: Boris Batyunya, JINR, Boris.Batiounia@cern.ch
  25 //     Params moved to AliITSRecoParam by: Andrea Dainese, INFN
  26 //     Material budget from TGeo by: Ludovic Gaudichet & Andrea Dainese, INFN
  27 //-------------------------------------------------------------------------
  28
  29 //#include <TMatrixD.h>
  30 #include <TTree.h>
  31 #include <TDatabasePDG.h>
  32 #include <TString.h>
  33 #include <TRandom.h>
  34 #include <TTreeStream.h>
  35
  36
  37 #include "AliLog.h"
  38 #include "AliITSCalibrationSPD.h"
  39 #include "AliITSCalibrationSDD.h"
  40 #include "AliITSCalibrationSSD.h"
  41 #include "AliCDBEntry.h"
  42 #include "AliCDBManager.h"
  43 #include "AliAlignObj.h"
  44 #include "AliTrackPointArray.h"
  45 #include "AliESDVertex.h"
  46 #include "AliESDEvent.h"
  47 #include "AliESDtrack.h"
  48 #include "AliV0.h"
  49 #include "AliHelix.h"
  50 #include "AliITSChannelStatus.h"
  51 #include "AliITSDetTypeRec.h"
  52 #include "AliITSRecPoint.h"
  53 #include "AliITSgeomTGeo.h"
  54 #include "AliITSReconstructor.h"
  55 #include "AliITSClusterParam.h"
  56 #include "AliITSsegmentation.h"
  57 #include "AliITSCalibration.h"
  58 #include "AliITSV0Finder.h"
  59 #include "AliITStrackerHLT.h"
  60 //#include "AliHLTTPCCATrackParam.h"
  61 //#include "AliHLTVertexer.h"
  62
  63
  64 ClassImp(AliITStrackerHLT)
  65
  66
  67
  68 Bool_t AliITStrackerHLT::TransportToX( AliExternalTrackParam *t, double x ) const
  69 {
  70   return t->PropagateTo( x, t->GetBz() );
  71 }
  72
  73 Bool_t AliITStrackerHLT::TransportToPhiX( AliExternalTrackParam *t, double phi, double x ) const
  74 {
  75   return t->Propagate( phi, x, t->GetBz() );
  76 }
  77
  78
  79 //------------------------------------------------------------------------
  80 Int_t AliITStrackerHLT::UpdateMI(AliHLTITSTrack* track, const AliITSRecPoint* cl,Double_t /*chi2*/,Int_t index) const
  81 {
  82   //
  83   // Update ITS track
  84   //
  85
  86   if (cl->GetQ()<=0) return 0;  // ingore the "virtual" clusters
  87
  88   Int_t layer = (index & 0xf0000000) >> 28;
  89
  90   // Take into account the mis-alignment (bring track to cluster plane)
  91
  92   Double_t xTrOrig=track->GetX();
  93   if (!TransportToX( track, xTrOrig + cl->GetX() ) ) return 0;
  94
  95   Double_t err2Y, err2Z;
  96
  97   GetClusterErrors2( layer, cl, track, err2Y, err2Z );
  98
  99   Double_t p[2]={ cl->GetY(), cl->GetZ()};
 100   Double_t cov[3]={err2Y, 0., err2Z};
 101
 102   Int_t updated = 1;
 103   //if( layer!=2 && layer!=3 )
 104   updated = track->AliExternalTrackParam::Update(p,cov);
 105
 106   int n = track->GetNumberOfClusters();
 107   track->SetClusterIndex(n,index);
 108   track->SetNumberOfClusters(n+1);
 109
 110   return updated;
 111 }
 112
 113
 114 AliHLTITSLayer AliITStrackerHLT::fgLayers[AliITSgeomTGeo::kNLayers]; // ITS layers
 115
 116 AliITStrackerHLT::AliITStrackerHLT()
 117   :AliTracker(),
 118    fEsd(0),
 119    fUseTGeo(3),
 120    fxOverX0Pipe(-1.),
 121    fxTimesRhoPipe(-1.),
 122    fDebugStreamer(0),
 123    fITSChannelStatus(0),
 124    fkDetTypeRec(0),
 125    fTracks()
 126 {
 127   //Default constructor
 128   Int_t i;
 129   for(i=0;i<4;i++) fSPDdetzcentre[i]=0.;
 130   for(i=0;i<2;i++) {fxOverX0Shield[i]=-1.;fxTimesRhoShield[i]=-1.;}
 131   for(i=0;i<6;i++) {fxOverX0Layer[i]=-1.;fxTimesRhoLayer[i]=-1.;}
 132 }
 133 //------------------------------------------------------------------------
 134 AliITStrackerHLT::AliITStrackerHLT(const Char_t *geom)
 135 : AliTracker(),
 136   fEsd(0),
 137   fUseTGeo(3),
 138   fxOverX0Pipe(-1.),
 139   fxTimesRhoPipe(-1.),
 140   fDebugStreamer(0),
 141   fITSChannelStatus(0),
 142   fkDetTypeRec(0),
 143   fTracks()
 144 {
 145   //--------------------------------------------------------------------
 146   //This is the AliITStrackerHLT constructor
 147   //--------------------------------------------------------------------
 148   if (geom) {
 149     AliWarning("\"geom\" is actually a dummy argument !");
 150   }
 151
 152   for (Int_t i=1; i<AliITSgeomTGeo::GetNLayers()+1; i++) {
 153     Int_t nlad=AliITSgeomTGeo::GetNLadders(i);
 154     Int_t ndet=AliITSgeomTGeo::GetNDetectors(i);
 155
 156     Double_t xyz[3], &x=xyz[0], &y=xyz[1], &z=xyz[2];
 157     AliITSgeomTGeo::GetOrigTranslation(i,1,1,xyz);
 158     Double_t poff=TMath::ATan2(y,x);
 159     Double_t zoff=z;
 160     Double_t r=TMath::Sqrt(x*x + y*y);
 161
 162     AliITSgeomTGeo::GetOrigTranslation(i,1,2,xyz);
 163     r += TMath::Sqrt(x*x + y*y);
 164     AliITSgeomTGeo::GetOrigTranslation(i,2,1,xyz);
 165     r += TMath::Sqrt(x*x + y*y);
 166     AliITSgeomTGeo::GetOrigTranslation(i,2,2,xyz);
 167     r += TMath::Sqrt(x*x + y*y);
 168     r*=0.25;
 169
 170     new (fgLayers+i-1) AliHLTITSLayer(r,poff,zoff,nlad,ndet);
 171
 172     for (Int_t j=1; j<nlad+1; j++) {
 173       for (Int_t k=1; k<ndet+1; k++) { //Fill this layer with detectors
 174         TGeoHMatrix m; AliITSgeomTGeo::GetOrigMatrix(i,j,k,m);
 175         const TGeoHMatrix *tm=AliITSgeomTGeo::GetTracking2LocalMatrix(i,j,k);
 176         m.Multiply(tm);
 177         Double_t txyz[3]={0.};
 178         xyz[0]=0.;xyz[1]=0.;xyz[2]=0.;
 179         m.LocalToMaster(txyz,xyz);
 180         r=TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]);
 181         Double_t phi=TMath::ATan2(xyz[1],xyz[0]);
 182
 183         if (phi<0) phi+=TMath::TwoPi();
 184         else if (phi>=TMath::TwoPi()) phi-=TMath::TwoPi();
 185
 186         AliHLTITSDetector &det=fgLayers[i-1].GetDetector((j-1)*ndet + k-1);
 187         new(&det) AliHLTITSDetector(r,phi);
 188         // compute the real radius (with misalignment)
 189         TGeoHMatrix mmisal(*(AliITSgeomTGeo::GetMatrix(i,j,k)));
 190         mmisal.Multiply(tm);
 191         xyz[0]=0.;xyz[1]=0.;xyz[2]=0.;
 192         mmisal.LocalToMaster(txyz,xyz);
 193         Double_t rmisal=TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]*xyz[1]);
 194         det.SetRmisal(rmisal);
 195
 196       } // end loop on detectors
 197     } // end loop on ladders
 198   } // end loop on layers
 199
 200
 201
 202   Double_t xyzVtx[]={AliITSReconstructor::GetRecoParam()->GetXVdef(),
 203                      AliITSReconstructor::GetRecoParam()->GetYVdef(),
 204                      AliITSReconstructor::GetRecoParam()->GetZVdef()};
 205   Double_t ersVtx[]={AliITSReconstructor::GetRecoParam()->GetSigmaXVdef(),
 206                      AliITSReconstructor::GetRecoParam()->GetSigmaYVdef(),
 207                      AliITSReconstructor::GetRecoParam()->GetSigmaZVdef()};
 208   SetVertex(xyzVtx,ersVtx);
 209
 210   // store positions of centre of SPD modules (in z)
 211   Double_t tr[3];
 212   AliITSgeomTGeo::GetTranslation(1,1,1,tr);
 213   fSPDdetzcentre[0] = tr[2];
 214   AliITSgeomTGeo::GetTranslation(1,1,2,tr);
 215   fSPDdetzcentre[1] = tr[2];
 216   AliITSgeomTGeo::GetTranslation(1,1,3,tr);
 217   fSPDdetzcentre[2] = tr[2];
 218   AliITSgeomTGeo::GetTranslation(1,1,4,tr);
 219   fSPDdetzcentre[3] = tr[2];
 220
 221   fUseTGeo = AliITSReconstructor::GetRecoParam()->GetUseTGeoInTracker();
 222   if(AliITSReconstructor::GetRecoParam()->GetExtendedEtaAcceptance() && fUseTGeo!=1 && fUseTGeo!=3) {
 223     AliWarning("fUseTGeo changed to 3 because fExtendedEtaAcceptance is kTRUE");
 224     fUseTGeo = 3;
 225   }
 226
 227   for(Int_t i=0;i<2;i++) {fxOverX0Shield[i]=-1.;fxTimesRhoShield[i]=-1.;}
 228   for(Int_t i=0;i<6;i++) {fxOverX0Layer[i]=-1.;fxTimesRhoLayer[i]=-1.;}
 229
 230   fDebugStreamer = new TTreeSRedirector("ITSdebug.root");
 231
 232 }
 233 //------------------------------------------------------------------------
 234 AliITStrackerHLT::AliITStrackerHLT(const AliITStrackerHLT &tracker)
 235 :AliTracker(tracker),
 236  fEsd(tracker.fEsd),
 237  fUseTGeo(tracker.fUseTGeo),
 238  fxOverX0Pipe(tracker.fxOverX0Pipe),
 239  fxTimesRhoPipe(tracker.fxTimesRhoPipe),
 240  fDebugStreamer(tracker.fDebugStreamer),
 241  fITSChannelStatus(tracker.fITSChannelStatus),
 242  fkDetTypeRec(tracker.fkDetTypeRec),
 243  fTracks()
 244 {
 245   //Copy constructor
 246   Int_t i;
 247   for(i=0;i<4;i++) {
 248     fSPDdetzcentre[i]=tracker.fSPDdetzcentre[i];
 249   }
 250   for(i=0;i<6;i++) {
 251     fxOverX0Layer[i]=tracker.fxOverX0Layer[i];
 252     fxTimesRhoLayer[i]=tracker.fxTimesRhoLayer[i];
 253   }
 254   for(i=0;i<2;i++) {
 255     fxOverX0Shield[i]=tracker.fxOverX0Shield[i];
 256     fxTimesRhoShield[i]=tracker.fxTimesRhoShield[i];
 257   }
 258 }
 259 //------------------------------------------------------------------------
 260 AliITStrackerHLT & AliITStrackerHLT::operator=(const AliITStrackerHLT &tracker){
 261   //Assignment operator
 262   this->~AliITStrackerHLT();
 263   new(this) AliITStrackerHLT(tracker);
 264   return *this;
 265 }
 266 //------------------------------------------------------------------------
 267 AliITStrackerHLT::~AliITStrackerHLT()
 268 {
 269   //
 270   //destructor
 271   //
 272   if (fDebugStreamer) {
 273     //fDebugStreamer->Close();
 274     delete fDebugStreamer;
 275   }
 276   if(fITSChannelStatus) delete fITSChannelStatus;
 277 }
 278
 279
 280 //------------------------------------------------------------------------
 281 void AliITStrackerHLT::ReadBadFromDetTypeRec() {
 282   //--------------------------------------------------------------------
 283   //This function read ITS bad detectors, chips, channels from AliITSDetTypeRec
 284   //i.e. from OCDB
 285   //--------------------------------------------------------------------
 286
 287   if(!AliITSReconstructor::GetRecoParam()->GetUseBadZonesFromOCDB()) return;
 288
 289   Info("ReadBadFromDetTypeRec","Reading info about bad ITS detectors and channels");
 290
 291   if(!fkDetTypeRec) Error("ReadBadFromDetTypeRec","AliITSDetTypeRec nof found!\n");
 292
 293   // ITS channels map
 294   if(fITSChannelStatus) delete fITSChannelStatus;
 295   fITSChannelStatus = new AliITSChannelStatus(fkDetTypeRec);
 296
 297   // ITS detectors and chips
 298   Int_t i=0,j=0,k=0,ndet=0;
 299   for (i=1; i<AliITSgeomTGeo::GetNLayers()+1; i++) {
 300     Int_t nBadDetsPerLayer=0;
 301     ndet=AliITSgeomTGeo::GetNDetectors(i);
 302     for (j=1; j<AliITSgeomTGeo::GetNLadders(i)+1; j++) {
 303       for (k=1; k<ndet+1; k++) {
 304         AliHLTITSDetector &det=fgLayers[i-1].GetDetector((j-1)*ndet + k-1);
 305         det.ReadBadDetectorAndChips(i-1,(j-1)*ndet + k-1,fkDetTypeRec);
 306         if(det.IsBad()) {nBadDetsPerLayer++;}
 307       } // end loop on detectors
 308     } // end loop on ladders
 309     Info("ReadBadFromDetTypeRec",Form("Layer %d: %d bad out of %d",i-1,nBadDetsPerLayer,ndet*AliITSgeomTGeo::GetNLadders(i)));
 310   } // end loop on layers
 311
 312   return;
 313 }
 314
 315
 316 //------------------------------------------------------------------------
 317 void AliITStrackerHLT::LoadClusters( std::vector<AliITSRecPoint> clusters)
 318 {
 319   //--------------------------------------------------------------------
 320   //This function loads ITS clusters
 321   //--------------------------------------------------------------------
 322
 323   //SignDeltas(clusters,GetZ());
 324
 325   for( unsigned int icl=0; icl<clusters.size(); icl++ ){
 326
 327     AliITSRecPoint &cl = clusters[icl];
 328     Int_t i=cl.GetLayer();
 329     //Int_t ndet=fgLayers[i].GetNdetectors();
 330     //Int_t detector=cl.GetDetectorIndex();
 331     if (!cl.Misalign()) AliWarning("Can't misalign this cluster !");
 332     fgLayers[i].InsertCluster(new AliITSRecPoint(cl));
 333   }
 334
 335   for( int i=0; i<AliITSgeomTGeo::GetNLayers(); i++ ){
 336     for( int detector = 0; detector<fgLayers[i].GetNdetectors(); detector++ ){
 337
 338       // add dead zone "virtual" cluster in SPD, if there is a cluster within
 339       // zwindow cm from the dead zone
 340       if (i<2 && AliITSReconstructor::GetRecoParam()->GetAddVirtualClustersInDeadZone()) {
 341         for (Float_t xdead = 0; xdead < AliITSRecoParam::GetSPDdetxlength(); xdead += (i+1.)*AliITSReconstructor::GetRecoParam()->GetXPassDeadZoneHits()) {
 342           Int_t lab[4]   = {0,0,0,detector};
 343           Int_t info[3]  = {0,0,i};
 344           Float_t q      = 0.; // this identifies virtual clusters
 345           Float_t hit[5] = {xdead,
 346                             0.,
 347                             AliITSReconstructor::GetRecoParam()->GetSigmaXDeadZoneHit2(),
 348                             AliITSReconstructor::GetRecoParam()->GetSigmaZDeadZoneHit2(),
 349                             q};
 350           Bool_t local   = kTRUE;
 351           Double_t zwindow = AliITSReconstructor::GetRecoParam()->GetZWindowDeadZone();
 352           hit[1] = fSPDdetzcentre[0]+0.5*AliITSRecoParam::GetSPDdetzlength();
 353           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 354             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 355           hit[1] = fSPDdetzcentre[1]-0.5*AliITSRecoParam::GetSPDdetzlength();
 356           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 357             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 358           hit[1] = fSPDdetzcentre[1]+0.5*AliITSRecoParam::GetSPDdetzlength();
 359           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 360             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 361           hit[1] = fSPDdetzcentre[2]-0.5*AliITSRecoParam::GetSPDdetzlength();
 362           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 363             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 364           hit[1] = fSPDdetzcentre[2]+0.5*AliITSRecoParam::GetSPDdetzlength();
 365           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 366             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 367           hit[1] = fSPDdetzcentre[3]-0.5*AliITSRecoParam::GetSPDdetzlength();
 368           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 369             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 370         }
 371       } // "virtual" clusters in SPD
 372     }
 373     fgLayers[i].ResetRoad(); //road defined by the cluster density
 374     fgLayers[i].SortClusters();
 375   }
 376 }
 377
 378
 379 //------------------------------------------------------------------------
 380 Int_t AliITStrackerHLT::LoadClusters(TTree *cTree) {
 381   //--------------------------------------------------------------------
 382   //This function loads ITS clusters
 383   //--------------------------------------------------------------------
 384   TBranch *branch=cTree->GetBranch("ITSRecPoints");
 385   if (!branch) {
 386     Error("LoadClusters"," can't get the branch !\n");
 387     return 1;
 388   }
 389
 390   static TClonesArray dummy("AliITSRecPoint",10000), *clusters=&dummy;
 391   branch->SetAddress(&clusters);
 392
 393   Int_t i=0,j=0,ndet=0;
 394   Int_t detector=0;
 395   for (i=0; i<AliITSgeomTGeo::GetNLayers(); i++) {
 396     ndet=fgLayers[i].GetNdetectors();
 397     Int_t jmax = j + fgLayers[i].GetNladders()*ndet;
 398     for (; j<jmax; j++) {
 399       if (!cTree->GetEvent(j)) continue;
 400       Int_t ncl=clusters->GetEntriesFast();
 401       SignDeltas(clusters,GetZ());
 402
 403       while (ncl--) {
 404         AliITSRecPoint *c=(AliITSRecPoint*)clusters->UncheckedAt(ncl);
 405         detector=c->GetDetectorIndex();
 406
 407         if (!c->Misalign()) AliWarning("Can't misalign this cluster !");
 408
 409         fgLayers[i].InsertCluster(new AliITSRecPoint(*c));
 410       }
 411       clusters->Delete();
 412       // add dead zone "virtual" cluster in SPD, if there is a cluster within
 413       // zwindow cm from the dead zone
 414       if (i<2 && AliITSReconstructor::GetRecoParam()->GetAddVirtualClustersInDeadZone()) {
 415         for (Float_t xdead = 0; xdead < AliITSRecoParam::GetSPDdetxlength(); xdead += (i+1.)*AliITSReconstructor::GetRecoParam()->GetXPassDeadZoneHits()) {
 416           Int_t lab[4]   = {0,0,0,detector};
 417           Int_t info[3]  = {0,0,i};
 418           Float_t q      = 0.; // this identifies virtual clusters
 419           Float_t hit[5] = {xdead,
 420                             0.,
 421                             AliITSReconstructor::GetRecoParam()->GetSigmaXDeadZoneHit2(),
 422                             AliITSReconstructor::GetRecoParam()->GetSigmaZDeadZoneHit2(),
 423                             q};
 424           Bool_t local   = kTRUE;
 425           Double_t zwindow = AliITSReconstructor::GetRecoParam()->GetZWindowDeadZone();
 426           hit[1] = fSPDdetzcentre[0]+0.5*AliITSRecoParam::GetSPDdetzlength();
 427           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 428             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 429           hit[1] = fSPDdetzcentre[1]-0.5*AliITSRecoParam::GetSPDdetzlength();
 430           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 431             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 432           hit[1] = fSPDdetzcentre[1]+0.5*AliITSRecoParam::GetSPDdetzlength();
 433           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 434             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 435           hit[1] = fSPDdetzcentre[2]-0.5*AliITSRecoParam::GetSPDdetzlength();
 436           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 437             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 438           hit[1] = fSPDdetzcentre[2]+0.5*AliITSRecoParam::GetSPDdetzlength();
 439           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 440             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 441           hit[1] = fSPDdetzcentre[3]-0.5*AliITSRecoParam::GetSPDdetzlength();
 442           if (TMath::Abs(fgLayers[i].GetDetector(detector).GetZmax()-hit[1])<zwindow)
 443             fgLayers[i].InsertCluster(new AliITSRecPoint(lab,hit,info,local));
 444         }
 445       } // "virtual" clusters in SPD
 446
 447     }
 448     //
 449     fgLayers[i].ResetRoad(); //road defined by the cluster density
 450     fgLayers[i].SortClusters();
 451   }
 452
 453   dummy.Clear();
 454
 455   return 0;
 456 }
 457 //------------------------------------------------------------------------
 458 void AliITStrackerHLT::UnloadClusters() {
 459   //--------------------------------------------------------------------
 460   //This function unloads ITS clusters
 461   //--------------------------------------------------------------------
 462   for (Int_t i=0; i<AliITSgeomTGeo::GetNLayers(); i++) fgLayers[i].ResetClusters();
 463 }
 464
 465
 466
 467 //------------------------------------------------------------------------
 468 Int_t AliITStrackerHLT::CorrectForTPCtoITSDeadZoneMaterial(AliHLTITSTrack *t) {
 469   //--------------------------------------------------------------------
 470   // Correction for the material between the TPC and the ITS
 471   //--------------------------------------------------------------------
 472   if (t->GetX() > AliITSRecoParam::Getriw()) {   // inward direction
 473       if (!t->PropagateToTGeo(AliITSRecoParam::Getriw(),1)) return 0;// TPC inner wall
 474       if (!t->PropagateToTGeo(AliITSRecoParam::Getrcd(),1)) return 0;// TPC central drum
 475       if (!t->PropagateToTGeo(AliITSRecoParam::Getrs(),1))  return 0;// ITS screen
 476   } else if (t->GetX() < AliITSRecoParam::Getrs()) {  // outward direction
 477       if (!t->PropagateToTGeo(AliITSRecoParam::Getrs(),1))        return 0;// ITS screen
 478       if (!t->PropagateToTGeo(AliITSRecoParam::Getrcd(),1))       return 0;// TPC central drum
 479       if (!t->PropagateToTGeo(AliITSRecoParam::Getriw()+0.001,1)) return 0;// TPC inner wall
 480   } else {
 481     printf("CorrectForTPCtoITSDeadZoneMaterial: Track is already in the dead zone !\n");
 482     return 0;
 483   }
 484
 485   return 1;
 486 }
 487
 488 #include "TStopwatch.h"
 489
 490 void AliITStrackerHLT::Reconstruct( std::vector<AliExternalTrackParam> tracksTPC )
 491 {
 492   //--------------------------------------------------------------------
 493   // This functions reconstructs ITS tracks
 494   // The clusters must be already loaded !
 495   //--------------------------------------------------------------------
 496
 497   Double_t pimass = TDatabasePDG::Instance()->GetParticle(211)->Mass();
 498   fTracks.clear();
 499
 500   for( unsigned int itr=0; itr<tracksTPC.size(); itr++ ){
 501
 502     AliHLTITSTrack tMI( tracksTPC[itr] );
 503     AliHLTITSTrack *t = &tMI;
 504     t->SetTPCtrackId( itr );
 505     t->SetMass(pimass);
 506     t->SetExpQ(0);
 507
 508     if (!CorrectForTPCtoITSDeadZoneMaterial(t))  continue;
 509
 510     //Int_t tpcLabel=t->GetLabel(); //save the TPC track label
 511
 512     FollowProlongationTree(t);
 513     int nclu=0;
 514     for(Int_t i=0; i<6; i++) {
 515       if( t->GetClusterIndex(i)>=0 ) nclu++;
 516     }
 517     //cout<<"N assigned ITS clusters = "<<nclu<<std::endl;
 518     if( nclu>0 ){
 519       t->SetLabel(-1);//tpcLabel);
 520       t->SetFakeRatio(1.);
 521       CookLabel(t,0.); //For comparison only
 522       //cout<<"label = "<<t->GetLabel()<<" / "<<tpcLabel<<endl;
 523       TransportToX(t, 0 );
 524       //cout<<"\n fill track : parameters at "<<t->GetX()<<": "<< TMath::Sqrt(TMath::Abs(t->GetSigmaY2()))<<" "<< TMath::Sqrt(TMath::Abs(t->GetSigmaY2()))<<endl;
 525         //t->Print();
 526       fTracks.push_back( *t );
 527     }
 528   }
 529
 530   //AliHLTVertexer vertexer;
 531   //vertexer.SetESD( event );
 532   //vertexer.FindPrimaryVertex();
 533   //vertexer.FindV0s();
 534 }
 535
 536
 537
 538 //------------------------------------------------------------------------
 539 Int_t AliITStrackerHLT::Clusters2Tracks(AliESDEvent *event) {
 540   //--------------------------------------------------------------------
 541   // This functions reconstructs ITS tracks
 542   // The clusters must be already loaded !
 543   //--------------------------------------------------------------------
 544   std::cout<<"\n\n ITS starts...\n"<<std::endl;
 545   TStopwatch timer;
 546
 547   fEsd = event;         // store pointer to the esd
 548   std::vector<AliExternalTrackParam> tracksTPC;
 549
 550   for( int itr=0; itr<event->GetNumberOfTracks(); itr++ ){
 551
 552     AliESDtrack *esdTrack = event->GetTrack(itr);
 553
 554     if ((esdTrack->GetStatus()&AliESDtrack::kTPCin)==0) continue;
 555     if (esdTrack->GetStatus()&AliESDtrack::kTPCout) continue;
 556     if (esdTrack->GetStatus()&AliESDtrack::kITSin) continue;
 557     if (esdTrack->GetKinkIndex(0)>0) continue;   //kink daughter
 558
 559     AliHLTITSTrack t(*esdTrack);
 560     t.SetTPCtrackId( itr );
 561     tracksTPC.push_back( t );
 562   }
 563
 564   Reconstruct( tracksTPC );
 565
 566   for( unsigned int itr=0; itr<fTracks.size(); itr++ ){
 567     AliHLTITSTrack &t = fTracks[itr];
 568     UpdateESDtrack(event->GetTrack(t.TPCtrackId()), &t, AliESDtrack::kITSin);
 569   }
 570
 571   //AliHLTVertexer vertexer;
 572   //vertexer.SetESD( event );
 573   //vertexer.FindPrimaryVertex();
 574   //vertexer.FindV0s();
 575
 576   timer.Stop();
 577   static double totalTime = 0;
 578   static int nEvnts = 0;
 579   totalTime+=timer.CpuTime();
 580   nEvnts++;
 581   std::cout<<"\n\n ITS tracker time = "<<totalTime/nEvnts<<" [s/ev]  for "<<nEvnts<<" events\n\n "<<std::endl;
 582   return 0;
 583 }
 584
 585
 586 //------------------------------------------------------------------------
 587 Int_t AliITStrackerHLT::PropagateBack(AliESDEvent * /*event*/) {
 588   return 0;
 589 }
 590
 591 //------------------------------------------------------------------------
 592 Int_t AliITStrackerHLT::RefitInward(AliESDEvent * /*event*/ ) {
 593   return 0;
 594 }
 595 //------------------------------------------------------------------------
 596 AliCluster *AliITStrackerHLT::GetCluster(Int_t index) const {
 597   //--------------------------------------------------------------------
 598   //       Return pointer to a given cluster
 599   //--------------------------------------------------------------------
 600   Int_t l=(index & 0xf0000000) >> 28;
 601   Int_t c=(index & 0x0fffffff) >> 00;
 602   return fgLayers[l].GetCluster(c);
 603 }
 604 //------------------------------------------------------------------------
 605 Bool_t AliITStrackerHLT::GetTrackPoint(Int_t index, AliTrackPoint& p) const {
 606   //--------------------------------------------------------------------
 607   // Get track space point with index i
 608   //--------------------------------------------------------------------
 609
 610   Int_t l=(index & 0xf0000000) >> 28;
 611   Int_t c=(index & 0x0fffffff) >> 00;
 612   AliITSRecPoint *cl = fgLayers[l].GetCluster(c);
 613   Int_t idet = cl->GetDetectorIndex();
 614
 615   Float_t xyz[3];
 616   Float_t cov[6];
 617   cl->GetGlobalXYZ(xyz);
 618   cl->GetGlobalCov(cov);
 619   p.SetXYZ(xyz, cov);
 620   p.SetCharge(cl->GetQ());
 621   p.SetDriftTime(cl->GetDriftTime());
 622   AliGeomManager::ELayerID iLayer = AliGeomManager::kInvalidLayer;
 623   switch (l) {
 624   case 0:
 625     iLayer = AliGeomManager::kSPD1;
 626     break;
 627   case 1:
 628     iLayer = AliGeomManager::kSPD2;
 629     break;
 630   case 2:
 631     iLayer = AliGeomManager::kSDD1;
 632     break;
 633   case 3:
 634     iLayer = AliGeomManager::kSDD2;
 635     break;
 636   case 4:
 637     iLayer = AliGeomManager::kSSD1;
 638     break;
 639   case 5:
 640     iLayer = AliGeomManager::kSSD2;
 641     break;
 642   default:
 643     AliWarning(Form("Wrong layer index in ITS (%d) !",l));
 644     break;
 645   };
 646   UShort_t volid = AliGeomManager::LayerToVolUID(iLayer,idet);
 647   p.SetVolumeID((UShort_t)volid);
 648   return kTRUE;
 649 }
 650 //------------------------------------------------------------------------
 651 Bool_t AliITStrackerHLT::GetTrackPointTrackingError(Int_t index,
 652                         AliTrackPoint& p, const AliESDtrack *t) {
 653   //--------------------------------------------------------------------
 654   // Get track space point with index i
 655   // (assign error estimated during the tracking)
 656   //--------------------------------------------------------------------
 657
 658   Int_t l=(index & 0xf0000000) >> 28;
 659   Int_t c=(index & 0x0fffffff) >> 00;
 660   const AliITSRecPoint *cl = fgLayers[l].GetCluster(c);
 661   Int_t idet = cl->GetDetectorIndex();
 662
 663   const AliHLTITSDetector &det=fgLayers[l].GetDetector(idet);
 664
 665   // tgphi and tglambda of the track in tracking frame with alpha=det.GetPhi
 666   Float_t detxy[2];
 667   detxy[0] = det.GetR()*TMath::Cos(det.GetPhi());
 668   detxy[1] = det.GetR()*TMath::Sin(det.GetPhi());
 669   Double_t alpha = t->GetAlpha();
 670   Double_t xdetintrackframe = detxy[0]*TMath::Cos(alpha)+detxy[1]*TMath::Sin(alpha);
 671   Float_t phi = TMath::ASin(t->GetSnpAt(xdetintrackframe,GetBz()));
 672   phi += alpha-det.GetPhi();
 673   Float_t tgphi = TMath::Tan(phi);
 674
 675   Float_t tgl = t->GetTgl(); // tgl about const along track
 676   Float_t expQ = TMath::Max(0.8*t->GetTPCsignal(),30.);
 677
 678   Float_t errlocalx,errlocalz;
 679   Bool_t addMisalErr=kFALSE;
 680   AliITSClusterParam::GetError(l,cl,tgl,tgphi,expQ,errlocalx,errlocalz,addMisalErr);
 681
 682   Float_t xyz[3];
 683   Float_t cov[6];
 684   cl->GetGlobalXYZ(xyz);
 685   //  cl->GetGlobalCov(cov);
 686   Float_t pos[3] = {0.,0.,0.};
 687   AliCluster tmpcl((UShort_t)cl->GetVolumeId(),pos[0],pos[1],pos[2],errlocalx*errlocalx,errlocalz*errlocalz,0);
 688   tmpcl.GetGlobalCov(cov);
 689
 690   p.SetXYZ(xyz, cov);
 691   p.SetCharge(cl->GetQ());
 692   p.SetDriftTime(cl->GetDriftTime());
 693
 694   AliGeomManager::ELayerID iLayer = AliGeomManager::kInvalidLayer;
 695   switch (l) {
 696   case 0:
 697     iLayer = AliGeomManager::kSPD1;
 698     break;
 699   case 1:
 700     iLayer = AliGeomManager::kSPD2;
 701     break;
 702   case 2:
 703     iLayer = AliGeomManager::kSDD1;
 704     break;
 705   case 3:
 706     iLayer = AliGeomManager::kSDD2;
 707     break;
 708   case 4:
 709     iLayer = AliGeomManager::kSSD1;
 710     break;
 711   case 5:
 712     iLayer = AliGeomManager::kSSD2;
 713     break;
 714   default:
 715     AliWarning(Form("Wrong layer index in ITS (%d) !",l));
 716     break;
 717   };
 718   UShort_t volid = AliGeomManager::LayerToVolUID(iLayer,idet);
 719
 720   p.SetVolumeID((UShort_t)volid);
 721   return kTRUE;
 722 }
 723
 724
 725 //------------------------------------------------------------------------
 726 void AliITStrackerHLT::FollowProlongationTree(AliHLTITSTrack * track )
 727 {
 728     //cout<<endl;
 729   for (Int_t ilayer=5; ilayer>=0; ilayer--) {
 730     //cout<<"\nLayer "<<ilayer<<endl;
 731
 732     AliHLTITSLayer &layer=fgLayers[ilayer];
 733
 734     //cout<<" shield material.. "<<endl;
 735
 736     // material between SSD and SDD, SDD and SPD
 737     if (ilayer==3 && !CorrectForShieldMaterial(track,"SDD","inward")) continue;
 738     if (ilayer==1 && !CorrectForShieldMaterial(track,"SPD","inward")) continue;
 739
 740     int idet;
 741
 742     {
 743       // propagate to the layer radius
 744
 745       Double_t r = layer.GetR(), phi,z;
 746       Double_t xToGo;
 747       //cout<<" propagate to layer R= "<<r<<" .."<<endl;
 748       if( !track->GetLocalXat(r,xToGo) ) continue;
 749       if( !TransportToX(track, xToGo) ) continue;
 750
 751       // detector number
 752
 753       if (!track->GetPhiZat(r,phi,z)) continue;
 754       idet=layer.FindDetectorIndex(phi,z);
 755       //cout<<" detector number = "<<idet<<endl;
 756    }
 757
 758
 759     //cout<<" correct for the layer material .. "<<endl;
 760
 761     // correct for the layer material
 762     {
 763       Double_t trackGlobXYZ1[3];
 764       if (!track->GetXYZ(trackGlobXYZ1)) continue;
 765       CorrectForLayerMaterial(track,ilayer,trackGlobXYZ1,"inward");
 766     }
 767
 768     // track outside layer acceptance in z
 769
 770     if( idet<0 ) continue;
 771
 772     // propagate to the intersection with the detector plane
 773
 774     //cout<<" propagate to the intersection with the detector .. "<<endl;
 775
 776     const AliHLTITSDetector &det=layer.GetDetector( idet );
 777     if (!TransportToPhiX( track, det.GetPhi(), det.GetR() ) ) continue;
 778
 779     // DEFINITION OF SEARCH ROAD AND CLUSTERS SELECTION
 780
 781     // road in global (rphi,z) [i.e. in tracking ref. system]
 782
 783     Double_t zmin,zmax,ymin,ymax;
 784
 785     //cout<<" ComputeRoad .. "<<endl;
 786     if (!ComputeRoad(track,ilayer,idet,zmin,zmax,ymin,ymax)) continue;
 787
 788     //cout<<" Road: y["<<ymin<<","<<ymax<<"], z["<<zmin<<","<<zmax<<"] "<<endl;
 789
 790     // select clusters in road
 791
 792     //cout<<" SelectClusters .. "<<endl;
 793     layer.SelectClusters(zmin,zmax,ymin,ymax);
 794
 795     // Define criteria for track-cluster association
 796
 797     Double_t msz = track->GetSigmaZ2() +
 798       AliITSReconstructor::GetRecoParam()->GetNSigmaZLayerForRoadZ()*
 799       AliITSReconstructor::GetRecoParam()->GetNSigmaZLayerForRoadZ()*
 800       AliITSReconstructor::GetRecoParam()->GetSigmaZ2(ilayer);
 801
 802     Double_t msy = track->GetSigmaY2() +
 803       AliITSReconstructor::GetRecoParam()->GetNSigmaYLayerForRoadY()*
 804       AliITSReconstructor::GetRecoParam()->GetNSigmaYLayerForRoadY()*
 805       AliITSReconstructor::GetRecoParam()->GetSigmaY2(ilayer);
 806
 807      msz *= AliITSReconstructor::GetRecoParam()->GetNSigma2RoadZNonC();
 808      msy *= AliITSReconstructor::GetRecoParam()->GetNSigma2RoadYNonC();
 809
 810      msz = 1./msz; // 1/RoadZ^2
 811      msy = 1./msy; // 1/RoadY^2
 812
 813      //
 814      // LOOP OVER ALL POSSIBLE TRACK PROLONGATIONS ON THIS LAYER
 815      //
 816
 817      const AliITSRecPoint *cl=0;
 818      Int_t clidx=-1;
 819      Double_t chi2trkcl=AliITSReconstructor::GetRecoParam()->GetMaxChi2(); // init with big value
 820      Bool_t deadzoneSPD=kFALSE;
 821
 822      // check if the road contains a dead zone
 823      Bool_t noClusters = !layer.GetNextCluster(clidx,kTRUE);
 824
 825      Double_t dz=0.5*(zmax-zmin);
 826      Double_t dy=0.5*(ymax-ymin);
 827
 828      Int_t dead = CheckDeadZone(track,ilayer,idet,dz,dy,noClusters);
 829
 830      // create a prolongation without clusters (check also if there are no clusters in the road)
 831
 832      if (dead==1) { // dead zone at z=0,+-7cm in SPD
 833        deadzoneSPD=kTRUE;
 834      }
 835
 836      clidx=-1;
 837
 838      //cout<<" loop over clusters in the road .. "<<endl;
 839
 840      // loop over clusters in the road
 841      const AliITSRecPoint *bestCluster=0;
 842      double bestChi2 = 1.e10;
 843      AliHLTITSTrack bestTrack( *track );
 844      int bestIdx = -1;
 845      while ((cl=layer.GetNextCluster(clidx))!=0) {
 846        //cout<<" cluster: "<<cl->GetX()<<" "<<cl->GetY()<<" "<<cl->GetZ()<<endl;
 847        AliHLTITSTrack t(*track);
 848        if (cl->GetQ()==0 && deadzoneSPD==kTRUE) continue;
 849
 850        Int_t idetc=cl->GetDetectorIndex();
 851
 852        //cout<<" cluster detector: "<<idetc<<endl;
 853
 854        if ( idet !=idetc ) { // new cluster's detector
 855          const AliHLTITSDetector &detc=layer.GetDetector(idetc);
 856          if (!TransportToPhiX( track, detc.GetPhi(),detc.GetR()) ) continue;
 857          t = *track;
 858          idet = idetc;
 859        }
 860
 861        // take into account misalignment (bring track to real detector plane)
 862
 863        if (!TransportToX( &t, t.GetX() + cl->GetX() ) ) continue;
 864        double chi2 = ( (t.GetZ()-cl->GetZ())*(t.GetZ()-cl->GetZ())*msz +
 865                        (t.GetY()-cl->GetY())*(t.GetY()-cl->GetY())*msy   );
 866        //cout<<" chi2="<<chi2<<endl;
 867        if ( chi2 < bestChi2 ){
 868          bestChi2 = chi2;
 869          bestCluster = cl;
 870          bestTrack = t;
 871          bestIdx = clidx;
 872          continue;
 873        }
 874      }
 875
 876      //cout<<" best chi2= "<<bestChi2<<endl;
 877
 878      if( bestCluster && bestChi2 <=1. ){
 879
 880        //cout<<" cluster found "<<endl;
 881        *track = bestTrack;
 882
 883
 884       // calculate track-clusters chi2
 885        chi2trkcl = GetPredictedChi2MI(track,bestCluster,ilayer);
 886        //cout<<" track-clusters chi2 = "<<chi2trkcl<<endl;
 887        //cout<<" max chi2 = "<<AliITSReconstructor::GetRecoParam()->GetMaxChi2s(ilayer)<<endl;
 888
 889        // chi2 cut
 890        if (chi2trkcl < AliITSReconstructor::GetRecoParam()->GetMaxChi2s(ilayer)) {
 891          if (bestCluster->GetQ()==0) deadzoneSPD=kTRUE; // only 1 prolongation with virtual cluster
 892          //cout<<"set index.."<<endl;
 893          //cout<<"set index ok"<<endl;
 894          if (bestCluster->GetQ()!=0) { // real cluster
 895            //cout<<" UpdateMI ... "<<endl;
 896            if (!UpdateMI(track,bestCluster,chi2trkcl,(ilayer<<28)+bestIdx)) {
 897              continue;
 898            }
 899          }
 900        }
 901      }
 902      //cout<<" goto next layer "<<endl;
 903
 904   }
 905 }
 906
 907
 908 //------------------------------------------------------------------------
 909 AliHLTITSLayer & AliITStrackerHLT::GetLayer(Int_t layer) const
 910 {
 911   //--------------------------------------------------------------------
 912   //
 913   //
 914   return fgLayers[layer];
 915 }
 916
 917
 918
 919 //------------------------------------------------------------------------
 920 void AliITStrackerHLT::CookLabel(AliHLTITSTrack *track,Float_t wrong) const
 921 {
 922   // get MC label for the track
 923
 924   Int_t mcLabel = -1;
 925
 926   vector<int> labels;
 927   Int_t nClusters = track->GetNumberOfClusters();
 928
 929   for (Int_t i=0; i<nClusters; i++){
 930     Int_t cindex = track->GetClusterIndex(i);
 931     //Int_t l=(cindex & 0xf0000000) >> 28;
 932     AliITSRecPoint *cl = (AliITSRecPoint*)GetCluster(cindex);
 933     if ( cl->GetLabel(0) >= 0 ) labels.push_back(cl->GetLabel(0)) ;
 934     if ( cl->GetLabel(1) >= 0 ) labels.push_back(cl->GetLabel(1)) ;
 935     if ( cl->GetLabel(2) >= 0 ) labels.push_back(cl->GetLabel(2)) ;
 936   }
 937   std::sort( labels.begin(), labels.end() );
 938   labels.push_back( -1 ); // put -1 to the end
 939   int labelMax = -1, labelCur = -1, nLabelsMax = 0, nLabelsCurr = 0;
 940   for ( unsigned int iLab = 0; iLab < labels.size(); iLab++ ) {
 941     if ( labels[iLab] != labelCur ) {
 942       if ( labelCur >= 0 && nLabelsMax< nLabelsCurr ) {
 943         nLabelsMax = nLabelsCurr;
 944         labelMax = labelCur;
 945       }
 946       labelCur = labels[iLab];
 947       nLabelsCurr = 0;
 948     }
 949     nLabelsCurr++;
 950   }
 951
 952   if( labelMax>=0 && nLabelsMax <= wrong * nClusters ) labelMax = -labelMax;
 953
 954   mcLabel = labelMax;
 955
 956   track->SetLabel( mcLabel );
 957 }
 958
 959
 960
 961 void AliITStrackerHLT::GetClusterErrors2( Int_t layer, const AliITSRecPoint *cluster, AliHLTITSTrack* track, double &err2Y, double &err2Z ) const
 962 {
 963   Float_t erry,errz;
 964   Float_t theta = track->GetTgl();
 965   Float_t phi   = track->GetSnp();
 966   phi = TMath::Abs(phi)*TMath::Sqrt(1./((1.-phi)*(1.+phi)));
 967   AliITSClusterParam::GetError(layer,cluster,theta,phi,track->GetExpQ(),erry,errz);
 968   err2Y = erry*erry;
 969   err2Z = errz*errz;
 970 }
 971
 972
 973 //------------------------------------------------------------------------
 974 Double_t AliITStrackerHLT::GetPredictedChi2MI(AliHLTITSTrack* track, const AliITSRecPoint *cluster,Int_t layer)
 975 {
 976   //
 977   // Compute predicted chi2
 978   //
 979
 980   AliHLTITSTrack t(*track);
 981   if (!t.Propagate( t.GetX()+cluster->GetX())) return 1000.;
 982
 983   Double_t  err2Y,err2Z;
 984   GetClusterErrors2( layer, cluster, &t, err2Y, err2Z );
 985
 986   Double_t chi2 = t.GetPredictedChi2(cluster->GetY(),cluster->GetZ(),err2Y,err2Z);
 987
 988   Float_t ny,nz;
 989   Float_t theta = track->GetTgl();
 990   Float_t phi   = track->GetSnp();
 991   phi = TMath::Abs(phi)*TMath::Sqrt(1./((1.-phi)*(1.+phi)));
 992   AliITSClusterParam::GetNTeor(layer,cluster,theta,phi,ny,nz);
 993   Double_t delta = cluster->GetNy()+cluster->GetNz()-nz-ny;
 994   if (delta>1){
 995     chi2+=0.5*TMath::Min(delta/2,2.);
 996     chi2+=2.*cluster->GetDeltaProbability();
 997   }
 998   return chi2;
 999 }
1000
1001
1002
1003
1004 //------------------------------------------------------------------------
1005 void AliITStrackerHLT::SignDeltas(const TObjArray *clusterArray, Float_t vz)
1006 {
1007   //
1008   // Clusters from delta electrons?
1009   //
1010   Int_t entries = clusterArray->GetEntriesFast();
1011   if (entries<4) return;
1012   AliITSRecPoint* cluster = (AliITSRecPoint*)clusterArray->At(0);
1013   Int_t layer = cluster->GetLayer();
1014   if (layer>1) return;
1015   Int_t index[10000];
1016   Int_t ncandidates=0;
1017   Float_t r = (layer>0)? 7:4;
1018   //
1019   for (Int_t i=0;i<entries;i++){
1020     AliITSRecPoint* cl0 = (AliITSRecPoint*)clusterArray->At(i);
1021     Float_t nz = 1+TMath::Abs((cl0->GetZ()-vz)/r);
1022     if (cl0->GetNy()+cl0->GetNz()<=5+2*layer+nz) continue;
1023     index[ncandidates] = i;  //candidate to belong to delta electron track
1024     ncandidates++;
1025     if (cl0->GetNy()+cl0->GetNz()>9+2*layer+nz) {
1026       cl0->SetDeltaProbability(1);
1027     }
1028   }
1029   //
1030   //
1031   //
1032   for (Int_t i=0;i<ncandidates;i++){
1033     AliITSRecPoint* cl0 = (AliITSRecPoint*)clusterArray->At(index[i]);
1034     if (cl0->GetDeltaProbability()>0.8) continue;
1035     //
1036     Int_t ncl = 0;
1037     Float_t y[100],z[100],sumy,sumz,sumy2, sumyz, sumw;
1038     sumy=sumz=sumy2=sumyz=sumw=0.0;
1039     for (Int_t j=0;j<ncandidates;j++){
1040       if (i==j) continue;
1041       AliITSRecPoint* cl1 = (AliITSRecPoint*)clusterArray->At(index[j]);
1042       //
1043       Float_t dz = cl0->GetZ()-cl1->GetZ();
1044       Float_t dy = cl0->GetY()-cl1->GetY();
1045       if (TMath::Sqrt(dz*dz+dy*dy)<0.2){
1046         Float_t weight = cl1->GetNy()+cl1->GetNz()-2;
1047         y[ncl] = cl1->GetY();
1048         z[ncl] = cl1->GetZ();
1049         sumy+= y[ncl]*weight;
1050         sumz+= z[ncl]*weight;
1051         sumy2+=y[ncl]*y[ncl]*weight;
1052         sumyz+=y[ncl]*z[ncl]*weight;
1053         sumw+=weight;
1054         ncl++;
1055       }
1056     }
1057     if (ncl<4) continue;
1058     Float_t det = sumw*sumy2  - sumy*sumy;
1059     Float_t delta=1000;
1060     if (TMath::Abs(det)>0.01){
1061       Float_t z0  = (sumy2*sumz - sumy*sumyz)/det;
1062       Float_t k   = (sumyz*sumw - sumy*sumz)/det;
1063       delta = TMath::Abs(cl0->GetZ()-(z0+k*cl0->GetY()));
1064     }
1065     else{
1066       Float_t z0  = sumyz/sumy;
1067       delta = TMath::Abs(cl0->GetZ()-z0);
1068     }
1069     if ( delta<0.05) {
1070       cl0->SetDeltaProbability(1-20.*delta);
1071     }
1072   }
1073 }
1074 //------------------------------------------------------------------------
1075 void AliITStrackerHLT::UpdateESDtrack(AliESDtrack *tESD, AliHLTITSTrack* track, ULong_t flags) const
1076 {
1077   //
1078   // Update ESD track
1079   //
1080   tESD->UpdateTrackParams(track,flags);
1081   AliHLTITSTrack * oldtrack = (AliHLTITSTrack*)(tESD->GetITStrack());
1082   if (oldtrack) delete oldtrack;
1083   tESD->SetITStrack(new AliHLTITSTrack(*track));
1084 }
1085
1086
1087
1088
1089 //------------------------------------------------------------------------
1090 void AliITStrackerHLT::BuildMaterialLUT(TString material) {
1091   //--------------------------------------------------------------------
1092   // Fill a look-up table with mean material
1093   //--------------------------------------------------------------------
1094
1095   Int_t n=1000;
1096   Double_t mparam[7];
1097   Double_t point1[3],point2[3];
1098   Double_t phi,cosphi,sinphi,z;
1099   // 0-5 layers, 6 pipe, 7-8 shields
1100   Double_t rmin[9]={ 3.5, 5.5,13.0,22.0,35.0,41.0, 2.0, 8.0,25.0};
1101   Double_t rmax[9]={ 5.5, 8.0,17.0,26.0,41.0,47.0, 3.0,10.5,30.0};
1102
1103   Int_t ifirst=0,ilast=0;
1104   if(material.Contains("Pipe")) {
1105     ifirst=6; ilast=6;
1106   } else if(material.Contains("Shields")) {
1107     ifirst=7; ilast=8;
1108   } else if(material.Contains("Layers")) {
1109     ifirst=0; ilast=5;
1110   } else {
1111     Error("BuildMaterialLUT","Wrong layer name\n");
1112   }
1113
1114   for(Int_t imat=ifirst; imat<=ilast; imat++) {
1115     Double_t param[5]={0.,0.,0.,0.,0.};
1116     for (Int_t i=0; i<n; i++) {
1117       phi = 2.*TMath::Pi()*gRandom->Rndm();
1118       cosphi = TMath::Cos(phi); sinphi = TMath::Sin(phi);
1119       z = 14.*(-1.+2.*gRandom->Rndm()); // SPD barrel
1120       point1[0] = rmin[imat]*cosphi;
1121       point1[1] = rmin[imat]*sinphi;
1122       point1[2] = z;
1123       point2[0] = rmax[imat]*cosphi;
1124       point2[1] = rmax[imat]*sinphi;
1125       point2[2] = z;
1126       AliTracker::MeanMaterialBudget(point1,point2,mparam);
1127       for(Int_t j=0;j<5;j++) param[j]+=mparam[j];
1128     }
1129     for(Int_t j=0;j<5;j++) param[j]/=(Float_t)n;
1130     if(imat<=5) {
1131       fxOverX0Layer[imat] = param[1];
1132       fxTimesRhoLayer[imat] = param[0]*param[4];
1133     } else if(imat==6) {
1134       fxOverX0Pipe = param[1];
1135       fxTimesRhoPipe = param[0]*param[4];
1136     } else if(imat==7) {
1137       fxOverX0Shield[0] = param[1];
1138       fxTimesRhoShield[0] = param[0]*param[4];
1139     } else if(imat==8) {
1140       fxOverX0Shield[1] = param[1];
1141       fxTimesRhoShield[1] = param[0]*param[4];
1142     }
1143   }
1144   /*
1145   printf("%s\n",material.Data());
1146   printf("%f  %f\n",fxOverX0Pipe,fxTimesRhoPipe);
1147   printf("%f  %f\n",fxOverX0Shield[0],fxTimesRhoShield[0]);
1148   printf("%f  %f\n",fxOverX0Shield[1],fxTimesRhoShield[1]);
1149   printf("%f  %f\n",fxOverX0Layer[0],fxTimesRhoLayer[0]);
1150   printf("%f  %f\n",fxOverX0Layer[1],fxTimesRhoLayer[1]);
1151   printf("%f  %f\n",fxOverX0Layer[2],fxTimesRhoLayer[2]);
1152   printf("%f  %f\n",fxOverX0Layer[3],fxTimesRhoLayer[3]);
1153   printf("%f  %f\n",fxOverX0Layer[4],fxTimesRhoLayer[4]);
1154   printf("%f  %f\n",fxOverX0Layer[5],fxTimesRhoLayer[5]);
1155   */
1156   return;
1157 }
1158
1159 //------------------------------------------------------------------------
1160 Int_t AliITStrackerHLT::CorrectForPipeMaterial(AliHLTITSTrack *t,
1161                                               TString direction) {
1162   //-------------------------------------------------------------------
1163   // Propagate beyond beam pipe and correct for material
1164   // (material budget in different ways according to fUseTGeo value)
1165   // Add time if going outward (PropagateTo or PropagateToTGeo)
1166   //-------------------------------------------------------------------
1167
1168   // Define budget mode:
1169   // 0: material from AliITSRecoParam (hard coded)
1170   // 1: material from TGeo in one step (on the fly)
1171   // 2: material from lut
1172   // 3: material from TGeo in one step (same for all hypotheses)
1173   Int_t mode;
1174   switch(fUseTGeo) {
1175   case 0:
1176     mode=0;
1177     break;
1178   case 1:
1179     mode=1;
1180     break;
1181   case 2:
1182     mode=2;
1183     break;
1184   case 3:
1185     mode=3;
1186     break;
1187   case 4:
1188     mode=3;
1189     break;
1190   default:
1191     mode=0;
1192     break;
1193   }
1194
1195   Float_t  dir = (direction.Contains("inward") ? 1. : -1.);
1196   Double_t rToGo=(dir>0 ? AliITSRecoParam::GetrInsidePipe() : AliITSRecoParam::GetrOutsidePipe());
1197   Double_t xToGo;
1198   if (!t->GetLocalXat(rToGo,xToGo)) return 0;
1199
1200   Double_t xOverX0,x0,lengthTimesMeanDensity;
1201
1202   switch(mode) {
1203   case 0:
1204     xOverX0 = AliITSRecoParam::GetdPipe();
1205     x0 = AliITSRecoParam::GetX0Be();
1206     lengthTimesMeanDensity = xOverX0*x0;
1207     lengthTimesMeanDensity *= dir;
1208     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1209     break;
1210   case 1:
1211     if (!t->PropagateToTGeo(xToGo,1)) return 0;
1212     break;
1213   case 2:
1214     if(fxOverX0Pipe<0) BuildMaterialLUT("Pipe");
1215     xOverX0 = fxOverX0Pipe;
1216     lengthTimesMeanDensity = fxTimesRhoPipe;
1217     lengthTimesMeanDensity *= dir;
1218     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1219     break;
1220   case 3:
1221     double xOverX0PipeTrks, xTimesRhoPipeTrks;
1222     if (!t->PropagateToTGeo(xToGo,1,xOverX0,lengthTimesMeanDensity)) return 0;
1223     Double_t angle=TMath::Sqrt((1.+t->GetTgl()*t->GetTgl())/
1224                                ((1.-t->GetSnp())*(1.+t->GetSnp())));
1225     xOverX0PipeTrks = TMath::Abs(xOverX0)/angle;
1226     xTimesRhoPipeTrks = TMath::Abs(lengthTimesMeanDensity)/angle;
1227     xOverX0 = xOverX0PipeTrks;
1228     lengthTimesMeanDensity = xTimesRhoPipeTrks;
1229     lengthTimesMeanDensity *= dir;
1230     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1231     break;
1232   }
1233
1234   return 1;
1235 }
1236 //------------------------------------------------------------------------
1237 Int_t AliITStrackerHLT::CorrectForShieldMaterial(AliHLTITSTrack *t,
1238                                                 TString shield,
1239                                                 TString direction) {
1240   //-------------------------------------------------------------------
1241   // Propagate beyond SPD or SDD shield and correct for material
1242   // (material budget in different ways according to fUseTGeo value)
1243   // Add time if going outward (PropagateTo or PropagateToTGeo)
1244   //-------------------------------------------------------------------
1245
1246   // Define budget mode:
1247   // 0: material from AliITSRecoParam (hard coded)
1248   // 1: material from TGeo in steps of X cm (on the fly)
1249   //    X = AliITSRecoParam::GetStepSizeTGeo()
1250   // 2: material from lut
1251   // 3: material from TGeo in one step (same for all hypotheses)
1252   Int_t mode;
1253   switch(fUseTGeo) {
1254   case 0:
1255     mode=0;
1256     break;
1257   case 1:
1258     mode=1;
1259     break;
1260   case 2:
1261     mode=2;
1262     break;
1263   case 3:
1264     mode=3;
1265     break;
1266   case 4:
1267     mode=3;
1268     break;
1269   default:
1270     mode=0;
1271     break;
1272   }
1273
1274
1275   Float_t  dir = (direction.Contains("inward") ? 1. : -1.);
1276   Double_t rToGo;
1277   Int_t    shieldindex=0;
1278   if (shield.Contains("SDD")) { // SDDouter
1279     rToGo=(dir>0 ? AliITSRecoParam::GetrInsideShield(1) : AliITSRecoParam::GetrOutsideShield(1));
1280     shieldindex=1;
1281   } else if (shield.Contains("SPD")) {        // SPDouter
1282     rToGo=(dir>0 ? AliITSRecoParam::GetrInsideShield(0) : AliITSRecoParam::GetrOutsideShield(0));
1283     shieldindex=0;
1284   } else {
1285     Error("CorrectForShieldMaterial"," Wrong shield name\n");
1286     return 0;
1287   }
1288   Double_t xToGo;
1289   if (!t->GetLocalXat(rToGo,xToGo)) return 0;
1290
1291   Double_t xOverX0,x0,lengthTimesMeanDensity;
1292   Int_t nsteps=1;
1293
1294   switch(mode) {
1295   case 0:
1296     xOverX0 = AliITSRecoParam::Getdshield(shieldindex);
1297     x0 = AliITSRecoParam::GetX0shield(shieldindex);
1298     lengthTimesMeanDensity = xOverX0*x0;
1299     lengthTimesMeanDensity *= dir;
1300     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1301     break;
1302   case 1:
1303     nsteps= (Int_t)(TMath::Abs(t->GetX()-xToGo)/AliITSReconstructor::GetRecoParam()->GetStepSizeTGeo())+1;
1304     if (!t->PropagateToTGeo(xToGo,nsteps)) return 0; // cross the material and apply correction
1305     break;
1306   case 2:
1307     if(fxOverX0Shield[shieldindex]<0) BuildMaterialLUT("Shields");
1308     xOverX0 = fxOverX0Shield[shieldindex];
1309     lengthTimesMeanDensity = fxTimesRhoShield[shieldindex];
1310     lengthTimesMeanDensity *= dir;
1311     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1312     break;
1313   case 3:
1314     if (!t->PropagateToTGeo(xToGo,1,xOverX0,lengthTimesMeanDensity)) return 0;
1315     Double_t angle=TMath::Sqrt((1.+t->GetTgl()*t->GetTgl())/
1316                                ((1.-t->GetSnp())*(1.+t->GetSnp())));
1317     double xOverX0ShieldTrks = TMath::Abs(xOverX0)/angle;
1318     double xTimesRhoShieldTrks = TMath::Abs(lengthTimesMeanDensity)/angle;
1319     xOverX0 = xOverX0ShieldTrks;
1320     lengthTimesMeanDensity = xTimesRhoShieldTrks;
1321     lengthTimesMeanDensity *= dir;
1322     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1323     break;
1324   }
1325
1326   return 1;
1327 }
1328 //------------------------------------------------------------------------
1329 Int_t AliITStrackerHLT::CorrectForLayerMaterial(AliHLTITSTrack *t,
1330                                                Int_t layerindex,
1331                                                Double_t oldGlobXYZ[3],
1332                                                TString direction) {
1333   //-------------------------------------------------------------------
1334   // Propagate beyond layer and correct for material
1335   // (material budget in different ways according to fUseTGeo value)
1336   // Add time if going outward (PropagateTo or PropagateToTGeo)
1337   //-------------------------------------------------------------------
1338
1339   // Define budget mode:
1340   // 0: material from AliITSRecoParam (hard coded)
1341   // 1: material from TGeo in stepsof X cm (on the fly)
1342   //    X = AliITSRecoParam::GetStepSizeTGeo()
1343   // 2: material from lut
1344   // 3: material from TGeo in one step (same for all hypotheses)
1345   Int_t mode;
1346   switch(fUseTGeo) {
1347   case 0:
1348     mode=0;
1349     break;
1350   case 1:
1351     mode=1;
1352     break;
1353   case 2:
1354     mode=2;
1355     break;
1356   case 3:
1357     mode=3;
1358     break;
1359   case 4:
1360     mode=3;
1361     break;
1362   default:
1363     mode=0;
1364     break;
1365   }
1366
1367
1368   Float_t  dir = (direction.Contains("inward") ? 1. : -1.);
1369
1370   Double_t r=fgLayers[layerindex].GetR();
1371   Double_t deltar=(layerindex<2 ? 0.10*r : 0.05*r);
1372
1373   Double_t rToGo=TMath::Sqrt(t->GetX()*t->GetX()+t->GetY()*t->GetY())-deltar*dir;
1374   Double_t xToGo;
1375   if (!t->GetLocalXat(rToGo,xToGo)) return 0;
1376
1377
1378   Double_t xOverX0=0.0,x0=0.0,lengthTimesMeanDensity=0.0;
1379   Int_t nsteps=1;
1380
1381   // back before material (no correction)
1382   Double_t rOld,xOld;
1383   rOld=TMath::Sqrt(oldGlobXYZ[0]*oldGlobXYZ[0]+oldGlobXYZ[1]*oldGlobXYZ[1]);
1384   if (!t->GetLocalXat(rOld,xOld)) return 0;
1385   if (!TransportToX( t, xOld)) return 0;
1386
1387   switch(mode) {
1388   case 0:
1389     xOverX0 = fgLayers[layerindex].GetThickness(t->GetY(),t->GetZ(),x0);
1390     lengthTimesMeanDensity = xOverX0*x0;
1391     lengthTimesMeanDensity *= dir;
1392     // Bring the track beyond the material
1393     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1394     break;
1395   case 1:
1396     nsteps = (Int_t)(TMath::Abs(xOld-xToGo)/AliITSReconstructor::GetRecoParam()->GetStepSizeTGeo())+1;
1397     if (!t->PropagateToTGeo(xToGo,nsteps)) return 0; // cross the material and apply correction
1398     break;
1399   case 2:
1400     if(fxOverX0Layer[layerindex]<0) BuildMaterialLUT("Layers");
1401     xOverX0 = fxOverX0Layer[layerindex];
1402     lengthTimesMeanDensity = fxTimesRhoLayer[layerindex];
1403     lengthTimesMeanDensity *= dir;
1404     // Bring the track beyond the material
1405     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1406     break;
1407   case 3:
1408     nsteps = (Int_t)(TMath::Abs(xOld-xToGo)/AliITSReconstructor::GetRecoParam()->GetStepSizeTGeo())+1;
1409     if (!t->PropagateToTGeo(xToGo,nsteps,xOverX0,lengthTimesMeanDensity)) return 0;
1410     Double_t angle=TMath::Sqrt((1.+t->GetTgl()*t->GetTgl())/
1411                                ((1.-t->GetSnp())*(1.+t->GetSnp())));
1412     double xOverX0LayerTrks = TMath::Abs(xOverX0)/angle;
1413     double xTimesRhoLayerTrks = TMath::Abs(lengthTimesMeanDensity)/angle;
1414     xOverX0 = xOverX0LayerTrks;
1415     lengthTimesMeanDensity = xTimesRhoLayerTrks;
1416     lengthTimesMeanDensity *= dir;
1417     if (!t->PropagateTo(xToGo,xOverX0,lengthTimesMeanDensity/xOverX0)) return 0;
1418     break;
1419   }
1420
1421
1422   return 1;
1423 }
1424
1425
1426 //------------------------------------------------------------------------
1427 Int_t AliITStrackerHLT::CheckDeadZone(AliHLTITSTrack *track,
1428                                      Int_t ilayer,Int_t idet,
1429                                      Double_t dz,Double_t dy,
1430                                      Bool_t noClusters) const {
1431   //-----------------------------------------------------------------
1432   // This method is used to decide whether to allow a prolongation
1433   // without clusters, because there is a dead zone in the road.
1434   // In this case the return value is > 0:
1435   // return 1: dead zone at z=0,+-7cm in SPD
1436   // return 2: all road is "bad" (dead or noisy) from the OCDB
1437   // return 3: something "bad" (dead or noisy) from the OCDB
1438   //-----------------------------------------------------------------
1439
1440   // check dead zones at z=0,+-7cm in the SPD
1441   if (ilayer<2 && !AliITSReconstructor::GetRecoParam()->GetAddVirtualClustersInDeadZone()) {
1442     Double_t zmindead[3]={fSPDdetzcentre[0] + 0.5*AliITSRecoParam::GetSPDdetzlength(),
1443                           fSPDdetzcentre[1] + 0.5*AliITSRecoParam::GetSPDdetzlength(),
1444                           fSPDdetzcentre[2] + 0.5*AliITSRecoParam::GetSPDdetzlength()};
1445     Double_t zmaxdead[3]={fSPDdetzcentre[1] - 0.5*AliITSRecoParam::GetSPDdetzlength(),
1446                           fSPDdetzcentre[2] - 0.5*AliITSRecoParam::GetSPDdetzlength(),
1447                           fSPDdetzcentre[3] - 0.5*AliITSRecoParam::GetSPDdetzlength()};
1448     for (Int_t i=0; i<3; i++)
1449       if (track->GetZ()-dz<zmaxdead[i] && track->GetZ()+dz>zmindead[i]) {
1450         AliDebug(2,Form("crack SPD %d",ilayer));
1451         return 1;
1452       }
1453   }
1454
1455   // check bad zones from OCDB
1456   if (!AliITSReconstructor::GetRecoParam()->GetUseBadZonesFromOCDB()) return 0;
1457
1458   if (idet<0) return 0;
1459
1460   AliHLTITSDetector &det=fgLayers[ilayer].GetDetector(idet);
1461
1462   Int_t detType=-1;
1463   Float_t detSizeFactorX=0.0001,detSizeFactorZ=0.0001;
1464   if (ilayer==0 || ilayer==1) {        // ----------  SPD
1465     detType = 0;
1466   } else if (ilayer==2 || ilayer==3) { // ----------  SDD
1467     detType = 1;
1468     detSizeFactorX *= 2.;
1469   } else if (ilayer==4 || ilayer==5) { // ----------  SSD
1470     detType = 2;
1471   }
1472   AliITSsegmentation *segm = (AliITSsegmentation*)fkDetTypeRec->GetSegmentationModel(detType);
1473   if (detType==2) segm->SetLayer(ilayer+1);
1474   Float_t detSizeX = detSizeFactorX*segm->Dx();
1475   Float_t detSizeZ = detSizeFactorZ*segm->Dz();
1476
1477   // check if the road overlaps with bad chips
1478   Float_t xloc,zloc;
1479   LocalModuleCoord(ilayer,idet,track,xloc,zloc);
1480   Float_t zlocmin = zloc-dz;
1481   Float_t zlocmax = zloc+dz;
1482   Float_t xlocmin = xloc-dy;
1483   Float_t xlocmax = xloc+dy;
1484   Int_t chipsInRoad[100];
1485
1486   // check if road goes out of detector
1487   Bool_t touchNeighbourDet=kFALSE;
1488   if (TMath::Abs(xlocmin)>0.5*detSizeX) {xlocmin=-0.5*detSizeX; touchNeighbourDet=kTRUE;}
1489   if (TMath::Abs(xlocmax)>0.5*detSizeX) {xlocmax=+0.5*detSizeX; touchNeighbourDet=kTRUE;}
1490   if (TMath::Abs(zlocmin)>0.5*detSizeZ) {zlocmin=-0.5*detSizeZ; touchNeighbourDet=kTRUE;}
1491   if (TMath::Abs(zlocmax)>0.5*detSizeZ) {zlocmax=+0.5*detSizeZ; touchNeighbourDet=kTRUE;}
1492   AliDebug(2,Form("layer %d det %d zmim zmax %f %f xmin xmax %f %f   %f %f",ilayer,idet,zlocmin,zlocmax,xlocmin,xlocmax,detSizeZ,detSizeX));
1493
1494   // check if this detector is bad
1495   if (det.IsBad()) {
1496     AliDebug(2,Form("lay %d  bad detector %d",ilayer,idet));
1497     if(!touchNeighbourDet) {
1498       return 2; // all detectors in road are bad
1499     } else {
1500       return 3; // at least one is bad
1501     }
1502   }
1503
1504   Int_t nChipsInRoad = segm->GetChipsInLocalWindow(chipsInRoad,zlocmin,zlocmax,xlocmin,xlocmax);
1505   AliDebug(2,Form("lay %d nChipsInRoad %d",ilayer,nChipsInRoad));
1506   if (!nChipsInRoad) return 0;
1507
1508   Bool_t anyBad=kFALSE,anyGood=kFALSE;
1509   for (Int_t iCh=0; iCh<nChipsInRoad; iCh++) {
1510     if (chipsInRoad[iCh]<0 || chipsInRoad[iCh]>det.GetNChips()-1) continue;
1511     AliDebug(2,Form("  chip %d bad %d",chipsInRoad[iCh],(Int_t)det.IsChipBad(chipsInRoad[iCh])));
1512     if (det.IsChipBad(chipsInRoad[iCh])) {
1513       anyBad=kTRUE;
1514     } else {
1515       anyGood=kTRUE;
1516     }
1517   }
1518
1519   if (!anyGood) {
1520     if(!touchNeighbourDet) {
1521       AliDebug(2,"all bad in road");
1522       return 2;  // all chips in road are bad
1523     } else {
1524       return 3; // at least a bad chip in road
1525     }
1526   }
1527
1528   if (anyBad) {
1529     AliDebug(2,"at least a bad in road");
1530     return 3; // at least a bad chip in road
1531   }
1532
1533
1534   if (!AliITSReconstructor::GetRecoParam()->GetUseSingleBadChannelsFromOCDB()
1535       || ilayer==4 || ilayer==5     // SSD
1536       || !noClusters) return 0;
1537
1538   // There are no clusters in road: check if there is at least
1539   // a bad SPD pixel or SDD anode
1540
1541   Int_t idetInITS=idet;
1542   for(Int_t l=0;l<ilayer;l++) idetInITS+=AliITSgeomTGeo::GetNLadders(l+1)*AliITSgeomTGeo::GetNDetectors(l+1);
1543
1544   if (fITSChannelStatus->AnyBadInRoad(idetInITS,zlocmin,zlocmax,xlocmin,xlocmax)) {
1545     AliDebug(2,Form("Bad channel in det %d of layer %d\n",idet,ilayer));
1546     return 3;
1547   }
1548   //if (fITSChannelStatus->FractionOfBadInRoad(idet,zlocmin,zlocmax,xlocmin,xlocmax) > AliITSReconstructor::GetRecoParam()->GetMinFractionOfBadInRoad()) return 3;
1549
1550   return 0;
1551 }
1552 //------------------------------------------------------------------------
1553 Bool_t AliITStrackerHLT::LocalModuleCoord(Int_t ilayer,Int_t idet,
1554                                        const AliHLTITSTrack *track,
1555                                        Float_t &xloc,Float_t &zloc) const {
1556   //-----------------------------------------------------------------
1557   // Gives position of track in local module ref. frame
1558   //-----------------------------------------------------------------
1559
1560   xloc=0.;
1561   zloc=0.;
1562
1563   if(idet<0) return kFALSE;
1564
1565   Int_t ndet=AliITSgeomTGeo::GetNDetectors(ilayer+1); // layers from 1 to 6
1566
1567   Int_t lad = Int_t(idet/ndet) + 1;
1568
1569   Int_t det = idet - (lad-1)*ndet + 1;
1570
1571   Double_t xyzGlob[3],xyzLoc[3];
1572
1573   AliHLTITSDetector &detector = fgLayers[ilayer].GetDetector(idet);
1574   // take into account the misalignment: xyz at real detector plane
1575   if(!track->GetXYZAt(detector.GetRmisal(),GetBz(),xyzGlob)) return kFALSE;
1576
1577   if(!AliITSgeomTGeo::GlobalToLocal(ilayer+1,lad,det,xyzGlob,xyzLoc)) return kFALSE;
1578
1579   xloc = (Float_t)xyzLoc[0];
1580   zloc = (Float_t)xyzLoc[2];
1581
1582   return kTRUE;
1583 }
1584
1585 //------------------------------------------------------------------------
1586 Bool_t AliITStrackerHLT::ComputeRoad(AliHLTITSTrack* track,Int_t ilayer,Int_t idet,Double_t &zmin,Double_t &zmax,Double_t &ymin,Double_t &ymax) const {
1587   //--------------------------------------------------------------------
1588   // This function computes the rectangular road for this track
1589   //--------------------------------------------------------------------
1590
1591
1592   AliHLTITSDetector &det = fgLayers[ilayer].GetDetector(idet);
1593   // take into account the misalignment: propagate track to misaligned detector plane
1594   if (!TransportToPhiX( track, det.GetPhi(),det.GetRmisal() ) ) return kFALSE;
1595
1596   Double_t dz=AliITSReconstructor::GetRecoParam()->GetNSigmaRoadZ()*
1597                     TMath::Sqrt(track->GetSigmaZ2() +
1598                     AliITSReconstructor::GetRecoParam()->GetNSigmaZLayerForRoadZ()*
1599                     AliITSReconstructor::GetRecoParam()->GetNSigmaZLayerForRoadZ()*
1600                     AliITSReconstructor::GetRecoParam()->GetSigmaZ2(ilayer));
1601   Double_t dy=AliITSReconstructor::GetRecoParam()->GetNSigmaRoadY()*
1602                     TMath::Sqrt(track->GetSigmaY2() +
1603                     AliITSReconstructor::GetRecoParam()->GetNSigmaYLayerForRoadY()*
1604                     AliITSReconstructor::GetRecoParam()->GetNSigmaYLayerForRoadY()*
1605                     AliITSReconstructor::GetRecoParam()->GetSigmaY2(ilayer));
1606
1607   // track at boundary between detectors, enlarge road
1608   Double_t boundaryWidth=AliITSRecoParam::GetBoundaryWidth();
1609   if ( (track->GetY()-dy < det.GetYmin()+boundaryWidth) ||
1610        (track->GetY()+dy > det.GetYmax()-boundaryWidth) ||
1611        (track->GetZ()-dz < det.GetZmin()+boundaryWidth) ||
1612        (track->GetZ()+dz > det.GetZmax()-boundaryWidth) ) {
1613     Float_t tgl = TMath::Abs(track->GetTgl());
1614     if (tgl > 1.) tgl=1.;
1615     Double_t deltaXNeighbDets=AliITSRecoParam::GetDeltaXNeighbDets();
1616     dz = TMath::Sqrt(dz*dz+deltaXNeighbDets*deltaXNeighbDets*tgl*tgl);
1617     Float_t snp = TMath::Abs(track->GetSnp());
1618     if (snp > AliITSReconstructor::GetRecoParam()->GetMaxSnp()) return kFALSE;
1619     dy = TMath::Sqrt(dy*dy+deltaXNeighbDets*deltaXNeighbDets*snp*snp);
1620   } // boundary
1621
1622   // add to the road a term (up to 2-3 mm) to deal with misalignments
1623   dy = TMath::Sqrt(dy*dy + AliITSReconstructor::GetRecoParam()->GetRoadMisal()*AliITSReconstructor::GetRecoParam()->GetRoadMisal());
1624   dz = TMath::Sqrt(dz*dz + AliITSReconstructor::GetRecoParam()->GetRoadMisal()*AliITSReconstructor::GetRecoParam()->GetRoadMisal());
1625
1626   Double_t r = fgLayers[ilayer].GetR();
1627   zmin = track->GetZ() - dz;
1628   zmax = track->GetZ() + dz;
1629   ymin = track->GetY() + r*det.GetPhi() - dy;
1630   ymax = track->GetY() + r*det.GetPhi() + dy;
1631
1632   // bring track back to idead detector plane
1633   if (!TransportToPhiX( track, det.GetPhi(),det.GetR())) return kFALSE;
1634
1635   return kTRUE;
1636 }