ITS/UPGRADE/AliITSUTrackerGlo.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, 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 /*
  17   TO CHECK
  18   GetBz usage
  19
  20  */
  21
  22 //-------------------------------------------------------------------------
  23 //  Implementation of the ITS Upgrade tracker mother class.
  24 //-------------------------------------------------------------------------
  25 #include <TTree.h>
  26 #include <Riostream.h>
  27 #include <TMath.h>
  28
  29 #include "AliITSUTrackerGlo.h"
  30 #include "AliESDEvent.h"
  31 #include "AliESDtrack.h"
  32 #include "AliITSURecoDet.h"
  33 #include "AliITSURecoSens.h"
  34 #include "AliITSUReconstructor.h"
  35 #include "AliITSReconstructor.h"
  36 #include "AliITSUSeed.h"
  37 #include "AliITSUAux.h"
  38 #include "AliITSUClusterPix.h"
  39 using namespace AliITSUAux;
  40 using namespace TMath;
  41
  42
  43 //----------------- tmp stuff -----------------
  44
  45 ClassImp(AliITSUTrackerGlo)
  46 //_________________________________________________________________________
  47 AliITSUTrackerGlo::AliITSUTrackerGlo(AliITSUReconstructor* rec)
  48 :  fReconstructor(rec)
  49   ,fITS(0)
  50   ,fCurrESDtrack(0)
  51   ,fCurrMass(kPionMass)
  52   ,fHypStore(100)
  53   ,fCurrHyp(0)
  54   ,fSeedsPool("AliITSUSeed",0)
  55   ,fTrCond()
  56 {
  57   // Default constructor
  58   if (rec) Init(rec);
  59 }
  60
  61 //_________________________________________________________________________
  62 AliITSUTrackerGlo::~AliITSUTrackerGlo()
  63 {
  64  // Default destructor
  65  //
  66   delete fITS;
  67   //
  68 }
  69
  70 //_________________________________________________________________________
  71 void AliITSUTrackerGlo::Init(AliITSUReconstructor* rec)
  72 {
  73   // init with external reconstructor
  74   //
  75   fITS = new AliITSURecoDet(rec->GetGeom(),"ITSURecoInterface");
  76   for (int ilr=fITS->GetNLayersActive();ilr--;) {
  77     fITS->GetLayerActive(ilr)->SetClusters(rec->GetClusters(ilr));
  78   }
  79   //
  80   fSeedsPool.ExpandCreateFast(1000); // RS TOCHECK
  81   //
  82   fTrCond.SetNLayers(fITS->GetNLayersActive());
  83   fTrCond.AddNewCondition(5);
  84   fTrCond.AddGroupPattern( (0x1<<0)|(0x1<<1) );
  85   fTrCond.AddGroupPattern( (0x1<<3)|(0x1<<4) );
  86   fTrCond.AddGroupPattern( (0x1<<5)|(0x1<<6) );
  87   //
  88   fTrCond.AddNewCondition(5);
  89   fTrCond.AddGroupPattern( (0x1<<0)|(0x1<<2) );
  90   fTrCond.AddGroupPattern( (0x1<<3)|(0x1<<4) );
  91   fTrCond.AddGroupPattern( (0x1<<5)|(0x1<<6) );
  92   //
  93   fTrCond.AddNewCondition(5);
  94   fTrCond.AddGroupPattern( (0x1<<1)|(0x1<<2) );
  95   fTrCond.AddGroupPattern( (0x1<<3)|(0x1<<4) );
  96   fTrCond.AddGroupPattern( (0x1<<5)|(0x1<<6) );
  97   //
  98   printf("Tracking Conditions: ");
  99   fTrCond.Print();
 100 }
 101
 102 //_________________________________________________________________________
 103 Int_t AliITSUTrackerGlo::Clusters2Tracks(AliESDEvent *esdEv)
 104 {
 105   //
 106   //
 107   AliITSUReconstructor::GetRecoParam()->Print();
 108   int nTrESD = esdEv->GetNumberOfTracks();
 109   fHypStore.Delete();
 110   if (fHypStore.GetSize()<nTrESD) fHypStore.Expand(nTrESD+100);
 111   //
 112   fITS->ProcessClusters();
 113   // select ESD tracks to propagate
 114   for (int itr=0;itr<nTrESD;itr++) {
 115     AliESDtrack *esdTr = esdEv->GetTrack(itr);
 116     AliInfo(Form("Processing track %d | MCLabel: %d",itr,esdTr->GetTPCLabel()));
 117     FindTrack(esdTr, itr);
 118   }
 119   //
 120   printf("Hypotheses for current event (N seeds in pool: %d, size: %d)\n",fSeedsPool.GetEntriesFast(),fSeedsPool.GetSize());
 121   fHypStore.Print();
 122   FinalizeHypotheses();
 123   //
 124   return 0;
 125 }
 126
 127 //_________________________________________________________________________
 128 Int_t AliITSUTrackerGlo::PropagateBack(AliESDEvent * /*event*/)
 129 {
 130   //
 131   // To be implemented
 132   //
 133
 134  Info("PropagateBack","To be implemented");
 135   return 0;
 136 }
 137
 138 //_________________________________________________________________________
 139 Int_t AliITSUTrackerGlo::RefitInward(AliESDEvent * /*event*/)
 140 {
 141   //
 142   // To be implemented
 143   //
 144
 145   Info("RefitInward","To be implemented");
 146   return 0;
 147 }
 148
 149 //_________________________________________________________________________
 150 Int_t AliITSUTrackerGlo::LoadClusters(TTree * treeRP)
 151 {
 152   // read from tree (if pointer provided) or directly from the ITS reco interface
 153   //
 154   return fReconstructor->LoadClusters(treeRP);
 155 }
 156
 157 //_________________________________________________________________________
 158 void AliITSUTrackerGlo::UnloadClusters()
 159 {
 160   //
 161   // To be implemented
 162   //
 163
 164   Info("UnloadClusters","To be implemented");
 165 }
 166 //_________________________________________________________________________
 167 AliCluster * AliITSUTrackerGlo::GetCluster(Int_t /*index*/) const
 168 {
 169   //
 170   // To be implemented
 171   //
 172
 173   Info("GetCluster","To be implemented");
 174   return 0x0;
 175 }
 176
 177 //_________________________________________________________________________
 178 Bool_t AliITSUTrackerGlo::NeedToProlong(AliESDtrack* esdTr)
 179 {
 180   // do we need to match this track to ITS?
 181   //
 182   static double bz = GetBz();
 183   if (!esdTr->IsOn(AliESDtrack::kTPCin) ||
 184       esdTr->IsOn(AliESDtrack::kTPCout) ||
 185       esdTr->IsOn(AliESDtrack::kITSin)  ||
 186       esdTr->GetKinkIndex(0)>0) return kFALSE;
 187   //
 188   if (esdTr->Pt()<AliITSUReconstructor::GetRecoParam()->GetMinPtForProlongation()) return kFALSE;
 189   //
 190   float dtz[2];
 191   esdTr->GetDZ(GetX(),GetY(),GetZ(),bz,dtz);
 192   // if track is not V0 candidata but has large offset wrt IP, reject it. RS TOCHECK
 193   if ( !(esdTr->GetV0Index(0)>0 && dtz[0]>AliITSUReconstructor::GetRecoParam()->GetMaxDforV0dghtrForProlongation())
 194        && (Abs(dtz[0])>AliITSUReconstructor::GetRecoParam()->GetMaxDForProlongation() ||
 195            Abs(dtz[1])>AliITSUReconstructor::GetRecoParam()->GetMaxDZForProlongation())) return kFALSE;
 196   //
 197   return kTRUE;
 198 }
 199
 200 //_________________________________________________________________________
 201 void AliITSUTrackerGlo::FindTrack(AliESDtrack* esdTr, Int_t esdID)
 202 {
 203   // find prolongaion candidates finding for single seed
 204   //
 205   AliITSUSeed seedUC;  // copy of the seed from the upper layer
 206   AliITSUSeed seedT;   // transient seed between the seedUC and new prolongation hypothesis
 207   //
 208   if (!NeedToProlong(esdTr)) return;  // are we interested in this track?
 209   if (!InitHypothesis(esdTr,esdID)) return;  // initialize prolongations hypotheses tree
 210   //
 211   AliITSURecoSens *hitSens[AliITSURecoSens::kNNeighbors+1];
 212   //
 213   TObjArray clArr; // container for transfer of clusters matching to seed
 214   //
 215   int nLrActive = fITS->GetNLayersActive();
 216   for (int ila=nLrActive;ila--;) {
 217     int ilaUp = ila+1;                         // prolong seeds from layer above
 218     //
 219     // for the outermost layer the seed is created from the ESD track
 220     int nSeedsUp = (ilaUp==nLrActive) ? 1 : fCurrHyp->GetNSeeds(ilaUp);
 221     //
 222     for (int isd=0;isd<nSeedsUp;isd++) {
 223       AliITSUSeed* seedU;
 224       if (ilaUp==nLrActive) {
 225         seedU = 0;
 226         seedUC.InitFromESDTrack(esdTr);
 227       }
 228       else {
 229         seedU = fCurrHyp->GetSeed(ilaUp,isd);  // seed on prev.active layer to prolong
 230         seedUC = *seedU;                       // its copy will be prolonged
 231         seedUC.SetParent(seedU);
 232       }
 233       seedUC.ResetFMatrix();                    // reset the matrix for propagation to next layer
 234       // go till next active layer
 235       AliInfo(Form("working on Lr:%d Seed:%d of %d",ila,isd,nSeedsUp));
 236       if (!TransportToLayer(&seedUC, fITS->GetLrIDActive(ilaUp), fITS->GetLrIDActive(ila)) ) {
 237         //
 238         AliInfo("Transport failed");
 239         // Check if the seed satisfies to track definition
 240         if (NeedToKill(&seedUC,kTransportFailed) && seedU) seedU->Kill();
 241         continue; // RS TODO: decide what to do with tracks stopped on higher layers w/o killing
 242       }
 243       AliITSURecoLayer* lrA = fITS->GetLayerActive(ila);
 244       if (!GetRoadWidth(&seedUC, ila)) { // failed to find road width on the layer
 245         if (NeedToKill(&seedUC,kRWCheckFailed) && seedU) seedU->Kill();
 246         continue;
 247       }
 248       int nsens = lrA->FindSensors(&fTrImpData[kTrPhi0], hitSens);  // find detectors which may be hit by the track
 249       AliInfo(Form("Will check %d sensors on lr:%d ",nsens,ila));
 250       //
 251       double bz = GetBz();
 252       for (int isn=nsens;isn--;) {
 253         seedT = seedUC;
 254         AliITSURecoSens* sens = hitSens[isn];
 255         //
 256         // We need to propagate the seed to sensor on lrA staying the frame of the sensor from prev.layer,
 257         // since the transport matrix should be defined in this frame.
 258         double xs; // X in the TF of current seed, corresponding to intersection with sensor plane
 259         if (!seedT.GetTrackingXAtXAlpha(sens->GetXTF(),sens->GetPhiTF(),bz, xs)) continue;
 260         if (!seedT.PropagateToX(xs,bz)) continue;
 261         //      if (!seedT.Rotate(sens->GetPhiTF())) continue;
 262         if (!seedT.RotateToAlpha(sens->GetPhiTF())) continue;
 263         //
 264         int clID0 = sens->GetFirstClusterId();
 265         for (int icl=sens->GetNClusters();icl--;) {
 266           int res = CheckCluster(&seedT,ila,clID0+icl);
 267           //
 268           if (res==kStopSearchOnSensor) break;     // stop looking on this sensor
 269           if (res==kClusterNotMatching) continue;  // cluster does not match
 270           // cluster is matching and it was added to the hypotheses tree
 271         }
 272       }
 273       // cluster search is done. Do we need ta have a version of this seed skipping current layer
 274       seedT.SetLr(ila);
 275       if (!NeedToKill(&seedT,kMissingCluster)) {
 276         AliITSUSeed* seedSkp = NewSeedFromPool(&seedT);
 277         double penalty = -AliITSUReconstructor::GetRecoParam()->GetMissPenalty(ila);
 278         // to do: make penalty to account for probability to miss the cluster for good reason
 279         seedSkp->SetChi2Cl(penalty);
 280         AddProlongationHypothesis(seedSkp,ila);
 281       }
 282     }
 283     ((TObjArray*)fCurrHyp->GetLayerSeeds(ila))->Sort();
 284     printf(">>> All hypotheses on lr %d: \n",ila);
 285     for (int ih=0;ih<fCurrHyp->GetNSeeds(ila);ih++) {
 286       printf(" #%3d ",ih); fCurrHyp->GetSeed(ila,ih)->Print();
 287       //
 288       /*
 289       if (ila!=0) continue;
 290       double vecL[5] = {0};
 291       double matL[15] = {0};
 292       AliITSUSeed* sp = fCurrHyp->GetSeed(ila,ih);
 293       while(sp->GetParent()) {
 294         sp->Smooth(vecL,matL);
 295         if (sp->GetLayerID()>=fITS->GetNLayersActive()-1) break;
 296         sp = (AliITSUSeed*)sp->GetParent();
 297       }
 298       */
 299     }
 300   }
 301   //
 302   SaveCurrentTrackHypotheses();
 303   //
 304 }
 305
 306 //_________________________________________________________________________
 307 Bool_t AliITSUTrackerGlo::InitHypothesis(AliESDtrack *esdTr, Int_t esdID)
 308 {
 309   // init prolongaion candidates finding for single seed
 310   fCurrHyp = GetTrackHyp(esdID);
 311   if (fCurrHyp) return kTRUE;
 312   //
 313   fCurrMass = esdTr->GetMass();
 314   fCurrESDtrack = esdTr;
 315   if (fCurrMass<kPionMass*0.9) fCurrMass = kPionMass; // don't trust to mu, e identification from TPCin
 316   //
 317   fCurrHyp = new AliITSUTrackHyp(fITS->GetNLayersActive());
 318   fCurrHyp->SetESDTrack(esdTr);
 319   fCurrHyp->SetUniqueID(esdID);
 320   SetTrackHyp(fCurrHyp,esdID);
 321   //
 322   return kTRUE;
 323   // TO DO
 324 }
 325
 326 //_________________________________________________________________________
 327 Bool_t AliITSUTrackerGlo::TransportToLayer(AliITSUSeed* seed, Int_t lFrom, Int_t lTo)
 328 {
 329   // transport seed from layerFrom to the entrance of layerTo
 330   //
 331   const double kToler = 1e-6; // tolerance for layer on-surface check
 332   //
 333   int dir = lTo > lFrom ? 1:-1;
 334   AliITSURecoLayer* lrFr = fITS->GetLayer(lFrom); // this can be 0 when extrapolation from TPC to ITS is requested
 335   Bool_t checkFirst = kTRUE;
 336   while(lFrom!=lTo) {
 337     double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
 338     if (lrFr) {
 339       Bool_t doLayer = kTRUE;
 340       double xToGo = dir>0 ? lrFr->GetRMax() : lrFr->GetRMin();
 341       if (checkFirst) { // do we need to track till the surface of the current layer ?
 342         checkFirst = kFALSE;
 343         if      (dir>0) { if (curR2-xToGo*xToGo>kToler) doLayer = kFALSE; } // on the surface or outside of the layer
 344         else if (dir<0) { if (xToGo*xToGo-curR2>kToler) doLayer = kFALSE; } // on the surface or outside of the layer
 345       }
 346       if (doLayer) {
 347         if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
 348         // go via layer to its boundary, applying material correction.
 349         if (!PropagateSeed(seed,xToGo,fCurrMass, lrFr->GetMaxStep())) return kFALSE;
 350       }
 351     }
 352     AliITSURecoLayer* lrTo =  fITS->GetLayer( (lFrom+=dir) );
 353     if (!lrTo) AliFatal(Form("Layer %d does not exist",lFrom));
 354     //
 355     // go the entrance of the layer, assuming no materials in between
 356     double xToGo = dir>0 ? lrTo->GetRMin() : lrTo->GetRMax();
 357     if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
 358     if (!PropagateSeed(seed,xToGo,fCurrMass,100, kFALSE )) return kFALSE;
 359     lrFr = lrTo;
 360   }
 361   return kTRUE;
 362   //
 363 }
 364
 365 //_________________________________________________________________________
 366 Bool_t AliITSUTrackerGlo::TransportToLayer(AliExternalTrackParam* seed, Int_t lFrom, Int_t lTo)
 367 {
 368   // transport track from layerFrom to the entrance of layerTo
 369   //
 370   const double kToler = 1e-6; // tolerance for layer on-surface check
 371   //
 372   int dir = lTo > lFrom ? 1:-1;
 373   AliITSURecoLayer* lrFr = fITS->GetLayer(lFrom); // this can be 0 when extrapolation from TPC to ITS is requested
 374   Bool_t checkFirst = kTRUE;
 375   while(lFrom!=lTo) {
 376     double curR2 = seed->GetX()*seed->GetX() + seed->GetY()*seed->GetY(); // current radius
 377     if (lrFr) {
 378       Bool_t doLayer = kTRUE;
 379       double xToGo = dir>0 ? lrFr->GetRMax() : lrFr->GetRMin();
 380       if (checkFirst) { // do we need to track till the surface of the current layer ?
 381         checkFirst = kFALSE;
 382         if      (dir>0) { if (curR2-xToGo*xToGo>kToler) doLayer = kFALSE; } // on the surface or outside of the layer
 383         else if (dir<0) { if (xToGo*xToGo-curR2>kToler) doLayer = kFALSE; } // on the surface or outside of the layer
 384       }
 385       if (doLayer) {
 386         if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
 387         // go via layer to its boundary, applying material correction.
 388         if (!PropagateSeed(seed,xToGo,fCurrMass, lrFr->GetMaxStep())) return kFALSE;
 389       }
 390     }
 391     AliITSURecoLayer* lrTo =  fITS->GetLayer( (lFrom+=dir) );
 392     if (!lrTo) AliFatal(Form("Layer %d does not exist",lFrom));
 393     //
 394     // go the entrance of the layer, assuming no materials in between
 395     double xToGo = dir>0 ? lrTo->GetRMin() : lrTo->GetRMax();
 396     if (!seed->GetXatLabR(xToGo,xToGo,GetBz(),dir)) return kFALSE;
 397     if (!PropagateSeed(seed,xToGo,fCurrMass,100, kFALSE )) return kFALSE;
 398     lrFr = lrTo;
 399   }
 400   return kTRUE;
 401   //
 402 }
 403
 404 //_________________________________________________________________________
 405 Bool_t AliITSUTrackerGlo::GetRoadWidth(AliITSUSeed* seed, int ilrA)
 406 {
 407   // calculate road width in terms of phi and z for the track which MUST be on the external radius of the layer
 408   // as well as some aux info
 409   double bz = GetBz();
 410   AliITSURecoLayer* lrA = fITS->GetLayerActive(ilrA);
 411   seed->GetXYZ(&fTrImpData[kTrXIn]);    // lab position at the entrance from above
 412   static AliExternalTrackParam sc;   // seed copy for manipalitions
 413   sc = *seed;
 414   //
 415   fTrImpData[kTrPhiIn] = ATan2(fTrImpData[kTrYIn],fTrImpData[kTrXIn]);
 416   if (!sc.Rotate(fTrImpData[kTrPhiIn])) return kFALSE; // go to the frame of the entry point into the layer
 417   double dr  = lrA->GetDR();                              // approximate X dist at the inner radius
 418   if (!sc.GetXYZAt(sc.GetX()-dr, bz, fTrImpData + kTrXOut)) {
 419     // special case: track does not reach inner radius, might be tangential
 420     double r = sc.GetD(0,0,bz);
 421     double x;
 422     if (!sc.GetXatLabR(r,x,bz,-1)) {
 423       sc.Print();
 424       AliFatal(Form("This should not happen: r=%f",r));
 425     }
 426     dr = Abs(sc.GetX() - x);
 427     if (!sc.GetXYZAt(x, bz, fTrImpData + kTrXOut)) {
 428       sc.Print();
 429       AliFatal(Form("This should not happen: x=%f",x));
 430     }
 431   }
 432   //
 433   fTrImpData[kTrPhiOut] = ATan2(fTrImpData[kTrYOut],fTrImpData[kTrXOut]);
 434   double sgy = sc.GetSigmaY2() + dr*dr*sc.GetSigmaSnp2() + AliITSUReconstructor::GetRecoParam()->GetSigmaY2(ilrA);
 435   double sgz = sc.GetSigmaZ2() + dr*dr*sc.GetSigmaTgl2() + AliITSUReconstructor::GetRecoParam()->GetSigmaZ2(ilrA);
 436   sgy = Sqrt(sgy)*AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadY();
 437   sgz = Sqrt(sgz)*AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadZ();
 438   fTrImpData[kTrPhi0] = 0.5*(fTrImpData[kTrPhiOut]+fTrImpData[kTrPhiIn]);
 439   fTrImpData[kTrZ0]   = 0.5*(fTrImpData[kTrZOut]+fTrImpData[kTrZIn]);
 440   fTrImpData[kTrDPhi] = 0.5*Abs(fTrImpData[kTrPhiOut]-fTrImpData[kTrPhiIn]) + sgy/lrA->GetR();
 441   fTrImpData[kTrDZ]   = 0.5*Abs(fTrImpData[kTrZOut]-fTrImpData[kTrZIn])   + sgz;
 442   //
 443   return kTRUE;
 444 }
 445
 446 //_________________________________________________________________________
 447 AliITSUSeed* AliITSUTrackerGlo::NewSeedFromPool(const AliITSUSeed* src)
 448 {
 449   // create new seed, optionally copying from the source
 450   return src ?
 451     new(fSeedsPool[fSeedsPool.GetEntriesFast()]) AliITSUSeed(*src) :
 452     new(fSeedsPool[fSeedsPool.GetEntriesFast()]) AliITSUSeed();
 453 }
 454
 455 //_________________________________________________________________________
 456 Int_t AliITSUTrackerGlo::CheckCluster(AliITSUSeed* track, Int_t lr, Int_t clID)
 457 {
 458   // Check if the cluster (in tracking frame!) is matching to track.
 459   // The track must be already propagated to sensor tracking frame.
 460   // Returns:  kStopSearchOnSensor if the search on given sensor should be stopped,
 461   //           kClusterMatching    if the cluster is matching
 462   //           kClusterMatching    otherwise
 463   //
 464   const double kTolerX = 5e-4;
 465   AliCluster *cl = fITS->GetLayerActive(lr)->GetCluster(clID);
 466   //
 467   Bool_t goodCl = kFALSE;
 468   int currLabel = Abs(fCurrESDtrack->GetTPCLabel());
 469   //
 470   if (cl->GetLabel(0)>=0) {for (int i=0;i<3;i++) if (cl->GetLabel(i)>=0 && cl->GetLabel(i)==currLabel) {goodCl = kTRUE; break;}}
 471   else goodCl = kTRUE;
 472   //
 473   if (TMath::Abs(cl->GetX())>kTolerX) { // if due to the misalingment X is large, propagate track only
 474     if (!track->PropagateParamOnlyTo(track->GetX()+cl->GetX(),GetBz())) {
 475       if (goodCl) {printf("Loose good cl: Failed propagation. |"); cl->Print();}
 476       return kStopSearchOnSensor; // propagation failed, seedT is intact
 477     }
 478   }
 479   double dy = cl->GetY()-track->GetY();
 480   double dz = cl->GetZ()-track->GetZ();
 481   //
 482   double dy2 = dy*dy;
 483   double tol2 = (track->GetSigmaY2() + AliITSUReconstructor::GetRecoParam()->GetSigmaY2(lr))*
 484     AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadY()*AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadY(); // RS TOOPTIMIZE
 485   if (dy2>tol2) {                          // the clusters are sorted in Z(col) then in Y(row).
 486     if (goodCl) {printf("Loose good cl: dy2=%e > tol2=%e |",dy2,tol2); cl->Print();}
 487     if (dy>0) return kStopSearchOnSensor;  // No chance that other cluster of this sensor will match (all Y's will be even larger)
 488     else      return kClusterNotMatching;   // Other clusters may match
 489   }
 490   double dz2 = dz*dz;
 491   tol2 = (track->GetSigmaZ2() + AliITSUReconstructor::GetRecoParam()->GetSigmaZ2(lr))*
 492     AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadZ()*AliITSUReconstructor::GetRecoParam()->GetNSigmaRoadZ(); // RS TOOPTIMIZE
 493   if (dz2>tol2) {
 494     if (goodCl) {printf("Loose good cl: dz2=%e > tol2=%e |",dz2,tol2); cl->Print();}
 495     return kClusterNotMatching; // Other clusters may match
 496   }
 497   //
 498   // check chi2
 499   Double_t p[2]={cl->GetY(), cl->GetZ()};
 500   Double_t cov[3]={cl->GetSigmaY2(), cl->GetSigmaYZ(), cl->GetSigmaZ2()};
 501   double chi2 = track->GetPredictedChi2(p,cov);
 502   if (chi2>AliITSUReconstructor::GetRecoParam()->GetMaxTr2ClChi2(lr)) {
 503     if (goodCl) {
 504       printf("Loose good cl: Chi2=%e > Chi2Max=%e |dy: %+.3e dz: %+.3e\n",
 505              chi2,AliITSUReconstructor::GetRecoParam()->GetMaxTr2ClChi2(lr),dy,dz);
 506       track->Print("etp");
 507       cl->Print("");
 508     }
 509     return kClusterNotMatching;
 510   }
 511   //
 512   track = NewSeedFromPool(track);  // input track will be reused, use its clone for updates
 513   if (!track->Update()) {
 514     if (goodCl) {printf("Loose good cl: Failed update |"); cl->Print();}
 515     return kClusterNotMatching;
 516   }
 517   track->SetChi2Cl(chi2);
 518   track->SetLrClusterID(lr,clID);
 519   cl->IncreaseClusterUsage();
 520   //
 521   track->SetFake(!goodCl);
 522   //
 523   AliInfo(Form("AddCl(%d) Cl%d lr:%d: dY:%+8.4f dZ:%+8.4f (MC: %5d %5d %5d) |Chi2=%f(%c)",
 524                goodCl,clID,lr,dy,dz2,cl->GetLabel(0),cl->GetLabel(1),cl->GetLabel(2), chi2, track->IsFake() ? '-':'+'));
 525   //
 526   AddProlongationHypothesis(track,lr);
 527   //
 528   return kClusterMatching;
 529 }
 530
 531 //_________________________________________________________________________
 532 Bool_t AliITSUTrackerGlo::NeedToKill(AliITSUSeed *seed, Int_t flag)
 533 {
 534   // check if the seed should not be discarded
 535   const UShort_t kMask = 0xffff;
 536   if (flag==kMissingCluster) {
 537     int lastChecked = seed->GetLayerID();
 538     UShort_t patt = seed->GetHitsPattern();
 539     if (lastChecked) patt |= ~(kMask<<lastChecked); // not all layers were checked, complete unchecked once by potential hits
 540     Bool_t seedOK = fTrCond.CheckPattern(patt);
 541     return !seedOK;
 542   }
 543   return kTRUE;
 544 }
 545
 546 //______________________________________________________________________________
 547 Bool_t AliITSUTrackerGlo::PropagateSeed(AliITSUSeed *seed, Double_t xToGo, Double_t mass, Double_t maxStep, Bool_t matCorr)
 548 {
 549   // propagate seed to given x applying material correction if requested
 550   const Double_t kEpsilon = 1e-5;
 551   Double_t xpos     = seed->GetX();
 552   Int_t dir         = (xpos<xToGo) ? 1:-1;
 553   Double_t xyz0[3],xyz1[3],param[7];
 554   //
 555   if (matCorr) seed->GetXYZ(xyz1);   //starting global position
 556   while ( (xToGo-xpos)*dir > kEpsilon){
 557     Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
 558     Double_t x    = xpos+step;
 559     Double_t bz=GetBz();   // getting the local Bz
 560     if (!seed->PropagateToX(x,bz))  return kFALSE;
 561     if (matCorr) {
 562       xyz0[0]=xyz1[0]; // global pos at the beginning of step
 563       xyz0[1]=xyz1[1];
 564       xyz0[2]=xyz1[2];
 565       seed->GetXYZ(xyz1);    //  // global pos at the end of step
 566       MeanMaterialBudget(xyz0,xyz1,param);
 567       Double_t xrho=param[0]*param[4], xx0=param[1];
 568       if (dir>0) xrho = -xrho; // outward should be negative
 569       if (!seed->ApplyMaterialCorrection(xx0,xrho,mass,kFALSE)) return kFALSE;
 570     }
 571     xpos = seed->GetX();
 572   }
 573   return kTRUE;
 574 }
 575
 576 //______________________________________________________________________________
 577 Bool_t AliITSUTrackerGlo::PropagateSeed(AliExternalTrackParam *seed, Double_t xToGo, Double_t mass, Double_t maxStep, Bool_t matCorr)
 578 {
 579   // propagate seed to given x applying material correction if requested
 580   const Double_t kEpsilon = 1e-5;
 581   Double_t xpos     = seed->GetX();
 582   Int_t dir         = (xpos<xToGo) ? 1:-1;
 583   Double_t xyz0[3],xyz1[3],param[7];
 584   //
 585   Bool_t updTime = dir>0 && seed->IsStartedTimeIntegral();
 586   if (matCorr || updTime) seed->GetXYZ(xyz1);   //starting global position
 587   while ( (xToGo-xpos)*dir > kEpsilon){
 588     Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
 589     Double_t x    = xpos+step;
 590     Double_t bz=GetBz();   // getting the local Bz
 591     if (!seed->PropagateTo(x,bz))  return kFALSE;
 592     double ds = 0;
 593     if (matCorr || updTime) {
 594       xyz0[0]=xyz1[0]; // global pos at the beginning of step
 595       xyz0[1]=xyz1[1];
 596       xyz0[2]=xyz1[2];
 597       seed->GetXYZ(xyz1);    //  // global pos at the end of step
 598       //
 599       if (matCorr) {
 600         MeanMaterialBudget(xyz0,xyz1,param);
 601         Double_t xrho=param[0]*param[4], xx0=param[1];
 602         if (dir>0) xrho = -xrho; // outward should be negative
 603         if (!seed->CorrectForMeanMaterial(xx0,xrho,mass)) return kFALSE;
 604         ds = param[4];
 605       }
 606       else { // matCorr is not requested but time integral is
 607         double d0 = xyz1[0]-xyz0[0];
 608         double d1 = xyz1[1]-xyz0[1];
 609         double d2 = xyz1[2]-xyz0[2];
 610         ds = TMath::Sqrt(d0*d0+d1*d1+d2*d2);
 611       }
 612     }
 613     if (updTime) seed->AddTimeStep(ds);
 614     //
 615     xpos = seed->GetX();
 616   }
 617   return kTRUE;
 618 }
 619
 620 //______________________________________________________________________________
 621 void AliITSUTrackerGlo::SaveCurrentTrackHypotheses()
 622 {
 623   // RS: shall we clean up killed seeds?
 624   fCurrHyp = 0;
 625   // TODO
 626
 627 }
 628
 629 //______________________________________________________________________________
 630 void AliITSUTrackerGlo::FinalizeHypotheses()
 631 {
 632   // select winner for each hypothesis, remove cl. sharing conflicts
 633   AliInfo("TODO");
 634   //
 635   int nh = fHypStore.GetEntriesFast();
 636   for (int ih=0;ih<nh;ih++) {
 637     AliITSUTrackHyp* hyp = (AliITSUTrackHyp*) fHypStore.UncheckedAt(ih);
 638     if (!hyp) continue;
 639     hyp->UpdateESD();
 640   }
 641
 642 }