]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDseedV1.cxx
Setters added (Magnus)
[u/mrichter/AliRoot.git] / TRD / AliTRDseedV1.cxx
index 123f8680da0ddfcbb85911f4e462d5f1666b4c02..3061dcd826838582246031483ba42f7c66ff596b 100644 (file)
 ////////////////////////////////////////////////////////////////////////////
 
 #include "TMath.h"
-#include "TLinearFitter.h"
-#include "TClonesArray.h" // tmp
 #include <TTreeStream.h>
 
 #include "AliLog.h"
 #include "AliMathBase.h"
+#include "AliRieman.h"
 #include "AliCDBManager.h"
-#include "AliTracker.h"
 
+#include "AliTRDReconstructor.h"
 #include "AliTRDpadPlane.h"
 #include "AliTRDcluster.h"
 #include "AliTRDseedV1.h"
 #include "Cal/AliTRDCalROC.h"
 #include "Cal/AliTRDCalDet.h"
 
-ClassImp(AliTRDseedV1)
+class AliTracker;
 
-TLinearFitter *AliTRDseedV1::fgFitterY = NULL;
-TLinearFitter *AliTRDseedV1::fgFitterZ = NULL;
+ClassImp(AliTRDseedV1)
 
 //____________________________________________________________________
 AliTRDseedV1::AliTRDseedV1(Int_t det) 
@@ -77,6 +75,7 @@ AliTRDseedV1::AliTRDseedV1(Int_t det)
   ,fDiffL(0.)
   ,fDiffT(0.)
   ,fClusterIdx(0)
+  ,fErrorMsg(0)
   ,fN(0)
   ,fDet(det)
   ,fPt(0.)
@@ -87,7 +86,6 @@ AliTRDseedV1::AliTRDseedV1(Int_t det)
   ,fZ(0.)
   ,fS2Y(0.)
   ,fS2Z(0.)
-  ,fC(0.)
   ,fChi2(0.)
 {
   //
@@ -95,7 +93,7 @@ AliTRDseedV1::AliTRDseedV1(Int_t det)
   //
   memset(fIndexes,0xFF,kNclusters*sizeof(fIndexes[0]));
   memset(fClusters, 0, kNclusters*sizeof(AliTRDcluster*));
-  memset(fPad, 0, 3*sizeof(Float_t));
+  memset(fPad, 0, 4*sizeof(Float_t));
   fYref[0] = 0.; fYref[1] = 0.; 
   fZref[0] = 0.; fZref[1] = 0.; 
   fYfit[0] = 0.; fYfit[1] = 0.; 
@@ -105,6 +103,8 @@ AliTRDseedV1::AliTRDseedV1(Int_t det)
   fLabels[0]=-1; fLabels[1]=-1; // most freq MC labels
   fLabels[2]=0;  // number of different labels for tracklet
   memset(fRefCov, 0, 7*sizeof(Double_t));
+  // stand alone curvature
+  fC[0] = 0.; fC[1] = 0.; 
   // covariance matrix [diagonal]
   // default sy = 200um and sz = 2.3 cm 
   fCov[0] = 4.e-4; fCov[1] = 0.; fCov[2] = 5.3; 
@@ -123,6 +123,7 @@ AliTRDseedV1::AliTRDseedV1(const AliTRDseedV1 &ref)
   ,fDiffL(0.)
   ,fDiffT(0.)
   ,fClusterIdx(0)
+  ,fErrorMsg(0)
   ,fN(0)
   ,fDet(-1)
   ,fPt(0.)
@@ -133,7 +134,6 @@ AliTRDseedV1::AliTRDseedV1(const AliTRDseedV1 &ref)
   ,fZ(0.)
   ,fS2Y(0.)
   ,fS2Z(0.)
-  ,fC(0.)
   ,fChi2(0.)
 {
   //
@@ -200,6 +200,7 @@ void AliTRDseedV1::Copy(TObject &ref) const
   target.fDiffL         = fDiffL;
   target.fDiffT         = fDiffT;
   target.fClusterIdx    = 0;
+  target.fErrorMsg      = fErrorMsg;
   target.fN             = fN;
   target.fDet           = fDet;
   target.fPt            = fPt;
@@ -210,12 +211,11 @@ void AliTRDseedV1::Copy(TObject &ref) const
   target.fZ             = fZ;
   target.fS2Y           = fS2Y;
   target.fS2Z           = fS2Z;
-  target.fC             = fC;
   target.fChi2          = fChi2;
   
   memcpy(target.fIndexes, fIndexes, kNclusters*sizeof(Int_t));
   memcpy(target.fClusters, fClusters, kNclusters*sizeof(AliTRDcluster*));
-  memcpy(target.fPad, fPad, 3*sizeof(Float_t));
+  memcpy(target.fPad, fPad, 4*sizeof(Float_t));
   target.fYref[0] = fYref[0]; target.fYref[1] = fYref[1]; 
   target.fZref[0] = fZref[0]; target.fZref[1] = fZref[1]; 
   target.fYfit[0] = fYfit[0]; target.fYfit[1] = fYfit[1]; 
@@ -224,12 +224,35 @@ void AliTRDseedV1::Copy(TObject &ref) const
   memcpy(target.fProb, fProb, AliPID::kSPECIES*sizeof(Float_t)); 
   memcpy(target.fLabels, fLabels, 3*sizeof(Int_t)); 
   memcpy(target.fRefCov, fRefCov, 7*sizeof(Double_t)); 
+  target.fC[0] = fC[0]; target.fC[1] = fC[1];
   memcpy(target.fCov, fCov, 3*sizeof(Double_t)); 
   
   TObject::Copy(ref);
 }
 
 
+//____________________________________________________________
+void AliTRDseedV1::Init(const AliRieman *rieman)
+{
+// Initialize this tracklet using the riemann fit information
+
+
+  fZref[0] = rieman->GetZat(fX0);
+  fZref[1] = rieman->GetDZat(fX0);
+  fYref[0] = rieman->GetYat(fX0);
+  fYref[1] = rieman->GetDYat(fX0);
+  if(fkReconstructor && fkReconstructor->IsHLT()){
+    fRefCov[0] = 1;
+    fRefCov[2] = 10;
+  }else{
+    fRefCov[0] = rieman->GetErrY(fX0);
+    fRefCov[2] = rieman->GetErrZ(fX0);
+  }
+  fC[0]    = rieman->GetC(); 
+  fChi2    = rieman->GetChi2();
+}
+
+
 //____________________________________________________________
 Bool_t AliTRDseedV1::Init(AliTRDtrackV1 *track)
 {
@@ -255,24 +278,28 @@ Bool_t AliTRDseedV1::Init(AliTRDtrackV1 *track)
 
 
 //_____________________________________________________________________________
-void AliTRDseedV1::Reset()
+void AliTRDseedV1::Reset(Option_t *opt)
 {
-  //
-  // Reset seed
-  //
+//
+// Reset seed. If option opt="c" is given only cluster arrays are cleared.
+//
+  for(Int_t ic=kNclusters; ic--;) fIndexes[ic] = -1;
+  memset(fClusters, 0, kNclusters*sizeof(AliTRDcluster*));
+  fN=0; SetBit(kRowCross, kFALSE);
+  if(strcmp(opt, "c")==0) return;
+
   fExB=0.;fVD=0.;fT0=0.;fS2PRF=0.;
   fDiffL=0.;fDiffT=0.;
   fClusterIdx=0;
-  fN=0;
+  fErrorMsg = 0;
   fDet=-1;
   fPt=0.;
   fdX=0.;fX0=0.; fX=0.; fY=0.; fZ=0.;
   fS2Y=0.; fS2Z=0.;
-  fC=0.; fChi2 = 0.;
+  fC[0]=0.; fC[1]=0.; 
+  fChi2 = 0.;
 
-  for(Int_t ic=kNclusters; ic--;) fIndexes[ic] = -1;
-  memset(fClusters, 0, kNclusters*sizeof(AliTRDcluster*));
-  memset(fPad, 0, 3*sizeof(Float_t));
+  memset(fPad, 0, 4*sizeof(Float_t));
   fYref[0] = 0.; fYref[1] = 0.; 
   fZref[0] = 0.; fZref[1] = 0.; 
   fYfit[0] = 0.; fYfit[1] = 0.; 
@@ -295,7 +322,7 @@ void AliTRDseedV1::Update(const AliTRDtrackV1 *trk)
   Double_t fSnp = trk->GetSnp();
   Double_t fTgl = trk->GetTgl();
   fPt = trk->Pt();
-  Double_t norm =1./TMath::Sqrt(1. - fSnp*fSnp); 
+  Double_t norm =1./TMath::Sqrt((1.-fSnp)*(1.+fSnp)); 
   fYref[1] = fSnp*norm;
   fZref[1] = fTgl*norm;
   SetCovRef(trk->GetCovariance());
@@ -386,44 +413,32 @@ void AliTRDseedV1::CookdEdx(Int_t nslices)
 // 3. cluster size
 //
 
-  Int_t nclusters[kNslices]; 
-  memset(nclusters, 0, kNslices*sizeof(Int_t));
   memset(fdEdx, 0, kNslices*sizeof(Float_t));
-
   const Double_t kDriftLength = (.5 * AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
 
-  AliTRDcluster *c = NULL;
+  AliTRDcluster *c(NULL);
   for(int ic=0; ic<AliTRDtrackerV1::GetNTimeBins(); ic++){
     if(!(c = fClusters[ic]) && !(c = fClusters[ic+kNtb])) continue;
     Float_t dx = TMath::Abs(fX0 - c->GetX());
-    
+
     // Filter clusters for dE/dx calculation
-    
+
     // 1.consider calibration effects for slice determination
     Int_t slice;
-    if(dx<kDriftLength){ // TODO should be replaced by c->IsInChamber() 
+    if(dx<kDriftLength){ // TODO should be replaced by c->IsInChamber()
       slice = Int_t(dx * nslices / kDriftLength);
     } else slice = c->GetX() < fX0 ? nslices-1 : 0;
 
 
     // 2. take sharing into account
     Float_t w = /*c->IsShared() ? .5 :*/ 1.;
-    
+
     // 3. take into account large clusters TODO
     //w *= c->GetNPads() > 3 ? .8 : 1.;
-    
+
     //CHECK !!!
     fdEdx[slice]   += w * GetdQdl(ic); //fdQdl[ic];
-    nclusters[slice]++;
   } // End of loop over clusters
-
-  //if(fkReconstructor->GetPIDMethod() == AliTRDReconstructor::kLQPID){
-  if(nslices == AliTRDpidUtil::kLQslices){
-  // calculate mean charge per slice (only LQ PID)
-    for(int is=0; is<nslices; is++){ 
-      if(nclusters[is]) fdEdx[is] /= nclusters[is];
-    }
-  }
 }
 
 //_____________________________________________________________________________
@@ -451,6 +466,77 @@ void AliTRDseedV1::CookLabels()
   if ((fLabels[2]  > 1) && (out[3] > 1)) fLabels[1] = out[2];
 }
 
+//____________________________________________________________
+Float_t AliTRDseedV1::GetAnodeWireOffset(Float_t zt)
+{
+// Find position inside the amplification cell for reading drift velocity map
+
+  Float_t d = fPad[3] - zt;
+  if(d<0.){
+    AliError(Form("Fail AnodeWireOffset calculation z0[%+7.2f] zt[%+7.2f] d[%+7.2f].", fPad[3], zt, d));
+    return 0.125;
+  } 
+  d -= ((Int_t)(2 * d)) / 2.0;
+  if(d > 0.25) d = 0.5 - d;
+  return d;
+}
+
+
+//____________________________________________________________________
+Float_t AliTRDseedV1::GetCharge(Bool_t useOutliers)
+{
+// Computes total charge attached to tracklet. If "useOutliers" is set clusters 
+// which are not in chamber are also used (default false)
+
+  AliTRDcluster *c(NULL); Float_t qt(0.);
+  for(int ic=0; ic<kNclusters; ic++){
+    if(!(c=fClusters[ic])) continue;
+    if(c->IsInChamber() && !useOutliers) continue;
+    qt += TMath::Abs(c->GetQ());
+  }
+  return qt;
+}
+
+//____________________________________________________________________
+Bool_t AliTRDseedV1::GetEstimatedCrossPoint(Float_t &x, Float_t &z) const
+{
+// Algorithm to estimate cross point in the x-z plane for pad row cross tracklets.
+// Returns true in case of success.
+  if(!IsRowCross()) return kFALSE;
+
+  x=0.; z=0.;
+  AliTRDcluster *c(NULL);
+  // Find radial range for first row
+  Float_t x1[] = {0., 1.e3};
+  for(int ic=0; ic<kNtb; ic++){
+    if(!(c=fClusters[ic])) continue;
+    if(!c->IsInChamber()) continue;
+    if(c->GetX() <= x1[1]) x1[1] = c->GetX();
+    if(c->GetX() >= x1[0]) x1[0] = c->GetX();
+    z=c->GetZ();
+  }
+  if((x1[0] - x1[1])<1.e-5) return kFALSE;
+
+  // Find radial range for second row
+  Bool_t kZ(kFALSE);
+  Float_t x2[] = {0., 1.e3};
+  for(int ic=kNtb; ic<kNclusters; ic++){
+    if(!(c=fClusters[ic])) continue;
+    if(!c->IsInChamber()) continue;
+    if(c->GetX() <= x2[1]) x2[1] = c->GetX();
+    if(c->GetX() >= x2[0]) x2[0] = c->GetX();
+    if(!kZ){
+      z+=c->GetZ();
+      z*=0.5;
+      kZ=kTRUE;
+    }
+  }
+  if((x2[0] - x2[1])<1.e-5) return kFALSE;
+
+  // Find intersection of the 2 radial regions
+  x = 0.5*((x1[0]+x1[1] > x2[0]+x2[1]) ? (x1[1]+x2[0]) : (x1[0]+x2[1]));
+  return kTRUE;
+}
 
 //____________________________________________________________________
 Float_t AliTRDseedV1::GetdQdl(Int_t ic, Float_t *dl) const
@@ -481,7 +567,8 @@ Float_t AliTRDseedV1::GetdQdl(Int_t ic, Float_t *dl) const
   if(fClusters[ic] && fClusters[ic]->IsInChamber()){
     hasClusterInChamber = kTRUE;
     dq += TMath::Abs(fClusters[ic]->GetQ());
-  }else if(fClusters[ic+kNtb] && fClusters[ic+kNtb]->IsInChamber()){
+  }
+  if(fClusters[ic+kNtb] && fClusters[ic+kNtb]->IsInChamber()){
     hasClusterInChamber = kTRUE;
     dq += TMath::Abs(fClusters[ic+kNtb]->GetQ());
   }
@@ -506,7 +593,8 @@ Float_t AliTRDseedV1::GetdQdl(Int_t ic, Float_t *dl) const
   }
   dx *= TMath::Sqrt(1. + fYfit[1]*fYfit[1] + fZref[1]*fZref[1]);
   if(dl) (*dl) = dx;
-  return dq/dx;
+  if(dx>1.e-9) return dq/dx;
+  else return 0.;
 }
 
 //____________________________________________________________
@@ -541,6 +629,20 @@ Float_t AliTRDseedV1::GetMomentum(Float_t *err) const
   return p;
 }
 
+//____________________________________________________________________
+Float_t AliTRDseedV1::GetOccupancyTB() const
+{
+// Returns procentage of TB occupied by clusters
+
+  Int_t n(0);
+  AliTRDcluster *c(NULL);
+  for(int ic=0; ic<AliTRDtrackerV1::GetNTimeBins(); ic++){
+    if(!(c = fClusters[ic]) && !(c = fClusters[ic+kNtb])) continue;
+    n++;
+  }
+
+  return Float_t(n)/AliTRDtrackerV1::GetNTimeBins();
+}
 
 //____________________________________________________________________
 Float_t* AliTRDseedV1::GetProbability(Bool_t force)
@@ -586,18 +688,18 @@ Bool_t AliTRDseedV1::CookPID()
     AliError("No access to AliTRDCalPID object");
     return kFALSE;
   }
-  //AliInfo(Form("Method[%d] : %s", fkReconstructor->GetRecoParam() ->GetPIDMethod(), pd->IsA()->GetName()));
 
   // calculate tracklet length TO DO
-  Float_t length = (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick());
-  /// TMath::Sqrt((1.0 - fSnp[iPlane]*fSnp[iPlane]) / (1.0 + fTgl[iPlane]*fTgl[iPlane]));
+  Float_t length = (AliTRDgeometry::AmThick() + AliTRDgeometry::DrThick())/ TMath::Sqrt((1.0 - GetSnp()*GetSnp()) / (1.0 + GetTgl()*GetTgl()));
   
   //calculate dE/dx
-  CookdEdx(fkReconstructor->GetNdEdxSlices());
-  
+  CookdEdx(AliTRDCalPID::kNSlicesNN);
+  AliDebug(4, Form("p=%6.4f[GeV/c] dEdx{%7.2f %7.2f %7.2f %7.2f %7.2f %7.2f %7.2f %7.2f} l=%4.2f[cm]", GetMomentum(), fdEdx[0], fdEdx[1], fdEdx[2], fdEdx[3], fdEdx[4], fdEdx[5], fdEdx[6], fdEdx[7], length));
+
   // Sets the a priori probabilities
+  Bool_t kPIDNN(fkReconstructor->GetPIDMethod()==AliTRDpidUtil::kNN);
   for(int ispec=0; ispec<AliPID::kSPECIES; ispec++)
-    fProb[ispec] = pd->GetProbability(ispec, GetMomentum(), &fdEdx[0], length, GetPlane());
+    fProb[ispec] = pd->GetProbability(ispec, GetMomentum(), &fdEdx[0], length, kPIDNN?GetPlane():fkReconstructor->GetRecoParam()->GetPIDLQslices());
   
   return kTRUE;
 }
@@ -677,18 +779,23 @@ void AliTRDseedV1::GetCovAt(Double_t x, Double_t *cov) const
     sy2 += sys[0];
     sz2 += sys[1];
   }
-  // rotate covariance matrix
-  Double_t t2 = GetTilt()*GetTilt();
-  Double_t correction = 1./(1. + t2);
-  cov[0] = (sy2+t2*sz2)*correction;
-  cov[1] = GetTilt()*(sz2 - sy2)*correction;
-  cov[2] = (t2*sy2+sz2)*correction;
-
-  //printf("C(%6.1f %+6.3f %6.1f)  [%s]\n", 1.e4*TMath::Sqrt(cov[0]), cov[1], 1.e4*TMath::Sqrt(cov[2]), IsRowCross()?" RC ":"-");
+
+  // rotate covariance matrix if no RC
+  if(!IsRowCross()){
+    Double_t t2 = GetTilt()*GetTilt();
+    Double_t correction = 1./(1. + t2);
+    cov[0] = (sy2+t2*sz2)*correction;
+    cov[1] = GetTilt()*(sz2 - sy2)*correction;
+    cov[2] = (t2*sy2+sz2)*correction;
+  } else {
+    cov[0] = sy2; cov[1] = 0.; cov[2] = sy2;
+  }
+
+  AliDebug(4, Form("C(%6.1f %+6.3f %6.1f)  RC[%c]", 1.e4*TMath::Sqrt(cov[0]), cov[1], 1.e4*TMath::Sqrt(cov[2]), IsRowCross()?'y':'n'));
 }
 
 //____________________________________________________________
-Double_t AliTRDseedV1::GetCovSqrt(const Double_t * const c, Double_t *d)
+Int_t AliTRDseedV1::GetCovSqrt(const Double_t * const c, Double_t *d)
 {
 // Helper function to calculate the square root of the covariance matrix. 
 // The input matrix is stored in the vector c and the result in the vector d. 
@@ -707,6 +814,7 @@ Double_t AliTRDseedV1::GetCovSqrt(const Double_t * const c, Double_t *d)
 // Author A.Bercuci <A.Bercuci@gsi.de>
 // Date   Mar 19 2009
 
+  const Double_t kZero(1.e-20);
   Double_t l[2], // eigenvalues
            v[3]; // eigenvectors
   // the secular equation and its solution :
@@ -715,26 +823,33 @@ Double_t AliTRDseedV1::GetCovSqrt(const Double_t * const c, Double_t *d)
   // L12 = [Tr(c) +- sqrt(Tr(c)^2-4*DET(c))]/2
   Double_t tr = c[0]+c[2],           // trace
           det = c[0]*c[2]-c[1]*c[1]; // determinant
-  if(TMath::Abs(det)<1.e-20) return -1.;
+  if(TMath::Abs(det)<kZero) return 1;
   Double_t dd = TMath::Sqrt(tr*tr - 4*det);
-  l[0] = .5*(tr + dd);
-  l[1] = .5*(tr - dd);
-  if(l[0]<0. || l[1]<0.) return -1.;
-
+  l[0] = .5*(tr + dd*(c[0]>c[2]?-1.:1.));
+  l[1] = .5*(tr + dd*(c[0]>c[2]?1.:-1.));
+  if(l[0]<kZero || l[1]<kZero) return 2;
   // the sym V matrix
   // | v00   v10|
   // | v10   v11|
-  Double_t tmp = (l[0]-c[0])/c[1];
-  v[0] = TMath::Sqrt(1./(tmp*tmp+1));
-  v[1] = tmp*v[0];
-  v[2] = v[1]*c[1]/(l[1]-c[2]);
+  Double_t den = (l[0]-c[0])*(l[0]-c[0])+c[1]*c[1];
+  if(den<kZero){ // almost diagonal
+    v[0] = TMath::Sign(0., c[1]);
+    v[1] = TMath::Sign(1., (l[0]-c[0]));
+    v[2] = TMath::Sign(0., c[1]*(l[0]-c[0])*(l[1]-c[2]));
+  } else {
+    Double_t tmp = 1./TMath::Sqrt(den);
+    v[0] = c[1]* tmp;
+    v[1] = (l[0]-c[0])*tmp;
+    if(TMath::Abs(l[1]-c[2])<kZero) v[2] = TMath::Sign(v[0]*(l[0]-c[0])/kZero, (l[1]-c[2]));
+    else v[2] = v[0]*(l[0]-c[0])/(l[1]-c[2]);
+  }
   // the VD^{1/2}V is: 
   l[0] = TMath::Sqrt(l[0]); l[1] = TMath::Sqrt(l[1]);
   d[0] = v[0]*v[0]*l[0]+v[1]*v[1]*l[1];
   d[1] = v[0]*v[1]*l[0]+v[1]*v[2]*l[1];
   d[2] = v[1]*v[1]*l[0]+v[2]*v[2]*l[1];
 
-  return 1.;
+  return 0;
 }
 
 //____________________________________________________________
@@ -760,26 +875,14 @@ Double_t AliTRDseedV1::GetCovInv(const Double_t * const c, Double_t *d)
 //____________________________________________________________________
 UShort_t AliTRDseedV1::GetVolumeId() const
 {
-  Int_t ic=0;
-  while(ic<kNclusters && !fClusters[ic]) ic++;
-  return fClusters[ic] ? fClusters[ic]->GetVolumeId() : 0;
-}
+// Returns geometry volume id by delegation 
 
-//____________________________________________________________________
-TLinearFitter* AliTRDseedV1::GetFitterY()
-{
-  if(!fgFitterY) fgFitterY = new TLinearFitter(1, "pol1");
-  fgFitterY->ClearPoints();
-  return fgFitterY;
+  for(Int_t ic(0);ic<kNclusters; ic++){
+    if(fClusters[ic]) return fClusters[ic]->GetVolumeId();
+  }
+  return 0;
 }
 
-//____________________________________________________________________
-TLinearFitter* AliTRDseedV1::GetFitterZ()
-{
-  if(!fgFitterZ) fgFitterZ = new TLinearFitter(1, "pol1");
-  fgFitterZ->ClearPoints();
-  return fgFitterZ;
-}
 
 //____________________________________________________________________
 void AliTRDseedV1::Calibrate()
@@ -825,12 +928,15 @@ void AliTRDseedV1::Calibrate()
     }
   }
 
-  fT0    = t0Det->GetValue(fDet) + t0ROC->GetValue(col,row);
+  fT0    = (t0Det->GetValue(fDet) + t0ROC->GetValue(col,row)) / AliTRDCommonParam::Instance()->GetSamplingFrequency();
   fVD    = vdDet->GetValue(fDet) * vdROC->GetValue(col, row);
   fS2PRF = calib->GetPRFWidth(fDet, col, row); fS2PRF *= fS2PRF;
   fExB   = AliTRDCommonParam::Instance()->GetOmegaTau(fVD);
   AliTRDCommonParam::Instance()->GetDiffCoeff(fDiffL,
   fDiffT, fVD);
+  AliDebug(4, Form("Calibration params for Det[%3d] Col[%3d] Row[%2d]\n  t0[%f]  vd[%f]  s2PRF[%f]  ExB[%f]  Dl[%f]  Dt[%f]", fDet, col, row, fT0, fVD, fS2PRF, fExB, fDiffL, fDiffT));
+
+
   SetBit(kCalib, kTRUE);
 }
 
@@ -848,13 +954,13 @@ void AliTRDseedV1::SetOwner()
 }
 
 //____________________________________________________________
-void AliTRDseedV1::SetPadPlane(AliTRDpadPlane *p)
+void AliTRDseedV1::SetPadPlane(AliTRDpadPlane * const p)
 {
 // Shortcut method to initialize pad geometry.
-  if(!p) return;
-  SetTilt(TMath::Tan(TMath::DegToRad()*p->GetTiltingAngle()));
-  SetPadLength(p->GetLengthIPad());
-  SetPadWidth(p->GetWidthIPad());
+  fPad[0] = p->GetLengthIPad();
+  fPad[1] = p->GetWidthIPad();
+  fPad[2] = TMath::Tan(TMath::DegToRad()*p->GetTiltingAngle());
+  fPad[3] = p->GetRow0() + p->GetAnodeWireOffset();
 }
 
 
@@ -888,42 +994,45 @@ Bool_t    AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *const chamber, Bool_t
 // Author : Alexandru Bercuci <A.Bercuci@gsi.de>
 // Debug  : level >3
 
-  Bool_t kPRINT = kFALSE;
-  if(!fkReconstructor->GetRecoParam() ){
-    AliError("Seed can not be used without a valid RecoParam.");
+  const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam(); //the dynamic cast in GetRecoParam is slow, so caching the pointer to it
+
+  if(!recoParam){
+    AliError("Tracklets can not be used without a valid RecoParam.");
     return kFALSE;
   }
   // Initialize reco params for this tracklet
   // 1. first time bin in the drift region
   Int_t t0 = 14;
-  Int_t kClmin = Int_t(fkReconstructor->GetRecoParam() ->GetFindableClusters()*AliTRDtrackerV1::GetNTimeBins());
+  Int_t kClmin = Int_t(recoParam->GetFindableClusters()*AliTRDtrackerV1::GetNTimeBins());
 
-  Double_t sysCov[5]; fkReconstructor->GetRecoParam()->GetSysCovMatrix(sysCov);        
+  Double_t sysCov[5]; recoParam->GetSysCovMatrix(sysCov);      
   Double_t s2yTrk= fRefCov[0], 
            s2yCl = 0., 
            s2zCl = GetPadLength()*GetPadLength()/12., 
            syRef = TMath::Sqrt(s2yTrk),
            t2    = GetTilt()*GetTilt();
   //define roads
-  Double_t kroady = 1., //fkReconstructor->GetRecoParam() ->GetRoad1y();
-           kroadz = GetPadLength() * fkReconstructor->GetRecoParam()->GetRoadzMultiplicator() + 1.;
+  Double_t kroady = 1., //recoParam->GetRoad1y();
+           kroadz = GetPadLength() * recoParam->GetRoadzMultiplicator() + 1.;
   // define probing cluster (the perfect cluster) and default calibration
   Short_t sig[] = {0, 0, 10, 30, 10, 0,0};
   AliTRDcluster cp(fDet, 6, 75, 0, sig, 0);
-  if(fkReconstructor->IsHLT())cp.SetRPhiMethod(AliTRDcluster::kCOG);
-  Calibrate();
+  if(fkReconstructor->IsHLT()) cp.SetRPhiMethod(AliTRDcluster::kCOG);
+  if(!IsCalibrated()) Calibrate();
 
-  if(kPRINT) printf("AttachClusters() sy[%f] road[%f]\n", syRef, kroady);
+  AliDebug(4, "");
+  AliDebug(4, Form("syKalman[%f] rY[%f] rZ[%f]", syRef, kroady, kroadz));
 
   // working variables
   const Int_t kNrows = 16;
   const Int_t kNcls  = 3*kNclusters; // buffer size
   AliTRDcluster *clst[kNrows][kNcls];
+  Bool_t blst[kNrows][kNcls];
   Double_t cond[4], dx, dy, yt, zt, yres[kNrows][kNcls];
   Int_t idxs[kNrows][kNcls], ncl[kNrows], ncls = 0;
   memset(ncl, 0, kNrows*sizeof(Int_t));
   memset(yres, 0, kNrows*kNcls*sizeof(Double_t));
-  memset(clst, 0, kNrows*kNcls*sizeof(AliTRDcluster*));
+  memset(blst, 0, kNrows*kNcls*sizeof(Bool_t));   //this is 8 times faster to memset than "memset(clst, 0, kNrows*kNcls*sizeof(AliTRDcluster*))"
 
   // Do cluster projection
   AliTRDcluster *c = NULL;
@@ -943,7 +1052,9 @@ Bool_t     AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *const chamber, Bool_t
     // get estimated road
     kroady = 3.*TMath::Sqrt(12.*(s2yTrk + s2yCl));
 
-    if(kPRINT) printf("  %2d dx[%f] yt[%f] zt[%f] sT[um]=%6.2f sy[um]=%6.2f syTilt[um]=%6.2f yRoad[mm]=%f\n", it, dx, yt, zt, 1.e4*TMath::Sqrt(s2yTrk), 1.e4*TMath::Sqrt(cp.GetSigmaY2()), 1.e4*TMath::Sqrt(s2yCl), 1.e1*kroady);
+    AliDebug(5, Form("  %2d x[%f] yt[%f] zt[%f]", it, dx, yt, zt));
+
+    AliDebug(5, Form("  syTrk[um]=%6.2f syCl[um]=%6.2f syClTlt[um]=%6.2f Ry[mm]=%f", 1.e4*TMath::Sqrt(s2yTrk), 1.e4*TMath::Sqrt(cp.GetSigmaY2()), 1.e4*TMath::Sqrt(s2yCl), 1.e1*kroady));
 
     // select clusters
     cond[0] = yt; cond[2] = kroady;
@@ -960,122 +1071,185 @@ Bool_t AliTRDseedV1::AttachClusters(AliTRDtrackingChamber *const chamber, Bool_t
         continue;
       }*/
       Int_t r = c->GetPadRow();
-      if(kPRINT) printf("\t\t%d dy[%f] yc[%f] r[%d]\n", ic, TMath::Abs(dy), c->GetY(), r);
+      AliDebug(5, Form("   -> dy[%f] yc[%f] r[%d]", TMath::Abs(dy), c->GetY(), r));
       clst[r][ncl[r]] = c;
+      blst[r][ncl[r]] = kTRUE;
       idxs[r][ncl[r]] = idx[ic];
       yres[r][ncl[r]] = dy;
       ncl[r]++; ncls++;
 
       if(ncl[r] >= kNcls) {
-        AliWarning(Form("Cluster candidates reached buffer limit %d. Some may be lost.", kNcls));
+        AliWarning(Form("Cluster candidates row[%d] reached buffer limit[%d]. Some may be lost.", r, kNcls));
         kBUFFER = kTRUE;
         break;
       }
     }
     if(kBUFFER) break;
   }
-  if(kPRINT) printf("Found %d clusters\n", ncls);
-  if(ncls<kClmin) return kFALSE;
+  AliDebug(4, Form("Found %d clusters. Processing ...", ncls));
+  if(ncls<kClmin){ 
+    AliDebug(1, Form("CLUSTERS FOUND %d LESS THAN THRESHOLD %d.", ncls, kClmin));
+    SetErrorMsg(kAttachClFound);
+    return kFALSE;
+  }
+
   // analyze each row individualy
-  Double_t mean, syDis;
-  Int_t nrow[] = {0, 0, 0}, nr = 0, lr=-1;
-  for(Int_t ir=kNrows; ir--;){
+  Bool_t kRowSelection(kFALSE);
+  Double_t mean[]={1.e3, 1.e3, 1.3}, syDis[]={1.e3, 1.e3, 1.3};
+  Int_t nrow[] = {0, 0, 0}, rowId[] = {-1, -1, -1}, nr = 0, lr=-1;
+  TVectorD vdy[3];
+  for(Int_t ir=0; ir<kNrows; ir++){
     if(!(ncl[ir])) continue;
-    if(lr>0 && lr-ir != 1){
-      if(kPRINT) printf("W - gap in rows attached !!\n"); 
+    if(lr>0 && ir-lr != 1){ 
+      AliDebug(2, "Rows attached not continuous. Turn on selection."); 
+      kRowSelection=kTRUE;
     }
-    if(kPRINT) printf("\tir[%d] lr[%d] n[%d]\n", ir, lr, ncl[ir]);
-    // Evaluate truncated mean on the y direction
-    if(ncl[ir] > 3) AliMathBase::EvaluateUni(ncl[ir], yres[ir], mean, syDis, Int_t(ncl[ir]*.8));
-    else {
-      mean = 0.; syDis = 0.;
-      continue;
-    } 
 
-    if(fkReconstructor->GetStreamLevel(AliTRDReconstructor::kTracker) > 3){
-      TTreeSRedirector &cstreamer = *fkReconstructor->GetDebugStream(AliTRDReconstructor::kTracker);
-      TVectorD vdy(ncl[ir], yres[ir]);
-      UChar_t stat(0);
-      if(IsKink()) SETBIT(stat, 0);
-      if(IsStandAlone()) SETBIT(stat, 1);
-      cstreamer << "AttachClusters"
-          << "stat="   << stat
-          << "det="    << fDet
-          << "pt="     << fPt
-          << "s2y="    << s2yTrk
-          << "dy="     << &vdy
-          << "m="      << mean
-          << "s="      << syDis
-          << "\n";
-    }
+    AliDebug(5, Form("  r[%d] n[%d]", ir, ncl[ir]));
+    // Evaluate truncated mean on the y direction
+    if(ncl[ir] < 4) continue;
+    AliMathBase::EvaluateUni(ncl[ir], yres[ir], mean[nr], syDis[nr], Int_t(ncl[ir]*.8));
 
     // TODO check mean and sigma agains cluster resolution !!
-    if(kPRINT) printf("\tr[%2d] m[%f %5.3fsigma] s[%f]\n", ir, mean, TMath::Abs(mean/syDis), syDis);
-    // select clusters on a 3 sigmaDistr level
+    AliDebug(4, Form("  m_%d[%+5.3f (%5.3fs)] s[%f]", nr, mean[nr], TMath::Abs(mean[nr]/syDis[nr]), syDis[nr]));
+    // remove outliers based on a 3 sigmaDistr level
     Bool_t kFOUND = kFALSE;
     for(Int_t ic = ncl[ir]; ic--;){
-      if(yres[ir][ic] - mean > 3. * syDis){ 
-        clst[ir][ic] = NULL; continue;
+      if(yres[ir][ic] - mean[nr] > 3. * syDis[nr]){ 
+        blst[ir][ic] = kFALSE; continue;
       }
-      nrow[nr]++; kFOUND = kTRUE;
+      nrow[nr]++; rowId[nr]=ir; kFOUND = kTRUE;
+    }
+    if(kFOUND){ 
+      vdy[nr].Use(nrow[nr], yres[ir]);
+      nr++; 
     }
-    // exit loop
-    if(kFOUND) nr++; 
     lr = ir; if(nr>=3) break;
   }
-  if(kPRINT) printf("lr[%d] nr[%d] nrow[0]=%d nrow[1]=%d nrow[2]=%d\n", lr, nr, nrow[0], nrow[1], nrow[2]);
-
-  // classify cluster rows
-  Int_t row = -1;
-  switch(nr){
-  case 1:
-    row = lr;
-    break;
-  case 2:
-    SetBit(kRowCross, kTRUE); // mark pad row crossing
-    if(nrow[0] > nrow[1]){ row = lr+1; lr = -1;}
-    else{ 
-      row = lr; lr = 1;
-      nrow[2] = nrow[1];
-      nrow[1] = nrow[0];
-      nrow[0] = nrow[2];
+  if(recoParam->GetStreamLevel(AliTRDrecoParam::kTracker) > 3 && fkReconstructor->IsDebugStreaming()){
+    TTreeSRedirector &cstreamer = *fkReconstructor->GetDebugStream(AliTRDrecoParam::kTracker);
+    UChar_t stat(0);
+    if(IsKink()) SETBIT(stat, 1);
+    if(IsStandAlone()) SETBIT(stat, 2);
+    cstreamer << "AttachClusters"
+        << "stat="   << stat
+        << "det="    << fDet
+        << "pt="     << fPt
+        << "s2y="    << s2yTrk
+        << "r0="     << rowId[0]
+        << "dy0="    << &vdy[0]
+        << "m0="     << mean[0]
+        << "s0="     << syDis[0]
+        << "r1="     << rowId[1]
+        << "dy1="    << &vdy[1]
+        << "m1="     << mean[1]
+        << "s1="     << syDis[1]
+        << "r2="     << rowId[2]
+        << "dy2="    << &vdy[2]
+        << "m2="     << mean[2]
+        << "s2="     << syDis[2]
+        << "\n";
+  }
+
+
+  // analyze gap in rows attached 
+  if(kRowSelection){
+    SetErrorMsg(kAttachRowGap);
+    Int_t rowRemove(-1); 
+    if(nr==2){ // select based on minimum distance to track projection
+      if(TMath::Abs(mean[0])<TMath::Abs(mean[1])){ 
+        if(nrow[1]>nrow[0]) AliDebug(2, Form("Conflicting mean[%f < %f] but ncl[%d < %d].", TMath::Abs(mean[0]), TMath::Abs(mean[1]), nrow[0], nrow[1]));
+      }else{
+        if(nrow[1]<nrow[0]) AliDebug(2, Form("Conflicting mean[%f > %f] but ncl[%d > %d].", TMath::Abs(mean[0]), TMath::Abs(mean[1]), nrow[0], nrow[1]));
+        Swap(nrow[0],nrow[1]); Swap(rowId[0],rowId[1]);
+        Swap(mean[0],mean[1]); Swap(syDis[0],syDis[1]);
+      }
+      rowRemove=1; nr=1; 
+    } else if(nr==3){ // select based on 2 consecutive rows
+      if(rowId[1]==rowId[0]+1 && rowId[1]!=rowId[2]-1){ 
+        nr=2;rowRemove=2;
+      } else if(rowId[1]!=rowId[0]+1 && rowId[1]==rowId[2]-1){ 
+        Swap(nrow[0],nrow[2]); Swap(rowId[0],rowId[2]);
+        Swap(mean[0],mean[2]); Swap(syDis[0],syDis[2]);
+        nr=2; rowRemove=2;
+      }
     }
-    break;
-  case 3:
+    if(rowRemove>0){nrow[rowRemove]=0; rowId[rowRemove]=-1;}
+  }
+  AliDebug(4, Form("  Ncl[%d[%d] + %d[%d] + %d[%d]]", nrow[0], rowId[0],  nrow[1], rowId[1], nrow[2], rowId[2]));
+
+  if(nr==3){
     SetBit(kRowCross, kTRUE); // mark pad row crossing
-    break;
+    SetErrorMsg(kAttachRow);
+    const Float_t am[]={TMath::Abs(mean[0]), TMath::Abs(mean[1]), TMath::Abs(mean[2])};
+    AliDebug(4, Form("complex row configuration\n"
+      "  r[%d] n[%d] m[%6.3f] s[%6.3f]\n"
+      "  r[%d] n[%d] m[%6.3f] s[%6.3f]\n"
+      "  r[%d] n[%d] m[%6.3f] s[%6.3f]\n"
+      , rowId[0], nrow[0], am[0], syDis[0]
+      , rowId[1], nrow[1], am[1], syDis[1]
+      , rowId[2], nrow[2], am[2], syDis[2]));
+    Int_t id[]={0,1,2}; TMath::Sort(3, am, id, kFALSE);
+    // backup
+    Int_t rnn[3]; memcpy(rnn, nrow, 3*sizeof(Int_t));
+    Int_t rid[3]; memcpy(rid, rowId, 3*sizeof(Int_t));
+    Double_t rm[3]; memcpy(rm, mean, 3*sizeof(Double_t));
+    Double_t rs[3]; memcpy(rs, syDis, 3*sizeof(Double_t));
+    nrow[0]=rnn[id[0]]; rowId[0]=rid[id[0]]; mean[0]=rm[id[0]]; syDis[0]=rs[id[0]];
+    nrow[1]=rnn[id[1]]; rowId[1]=rid[id[1]]; mean[1]=rm[id[1]]; syDis[1]=rs[id[1]];
+    nrow[2]=0;          rowId[2]=-1; mean[2] = 1.e3; syDis[2] = 1.e3;
+    AliDebug(4, Form("solved configuration\n"
+      "  r[%d] n[%d] m[%+6.3f] s[%6.3f]\n"
+      "  r[%d] n[%d] m[%+6.3f] s[%6.3f]\n"
+      "  r[%d] n[%d] m[%+6.3f] s[%6.3f]\n"
+      , rowId[0], nrow[0], mean[0], syDis[0]
+      , rowId[1], nrow[1], mean[1], syDis[1]
+      , rowId[2], nrow[2], mean[2], syDis[2]));
+    nr=2;
+  } else if(nr==2) {
+    SetBit(kRowCross, kTRUE); // mark pad row crossing
+    if(nrow[1] > nrow[0]){ // swap row order
+      Swap(nrow[0],nrow[1]); Swap(rowId[0],rowId[1]);
+      Swap(mean[0],mean[1]); Swap(syDis[0],syDis[1]);
+    }
   }
-  if(kPRINT) printf("\trow[%d] n[%d]\n\n", row, nrow[0]);
-  if(row<0) return kFALSE;
 
   // Select and store clusters 
   // We should consider here :
   //  1. How far is the chamber boundary
   //  2. How big is the mean
-  Int_t n = 0;
+  Int_t n(0); Float_t dyc[kNclusters]; memset(dyc,0,kNclusters*sizeof(Float_t));
   for (Int_t ir = 0; ir < nr; ir++) {
-    Int_t jr = row + ir*lr; 
-    if(kPRINT) printf("\tattach %d clusters for row %d\n", ncl[jr], jr);
+    Int_t jr(rowId[ir]);
+    AliDebug(4, Form("  Attaching Ncl[%d]=%d ...", jr, ncl[jr]));
     for (Int_t ic = 0; ic < ncl[jr]; ic++) {
-      if(!(c = clst[jr][ic])) continue;
-      Int_t it = c->GetPadTime();
+      if(!blst[jr][ic])continue;
+      c = clst[jr][ic];
+      Int_t it(c->GetPadTime());
+      Int_t idx(it+kNtb*ir);
+      if(fClusters[idx]){
+        AliDebug(4, Form("Many cluster candidates on row[%2d] tb[%2d].", jr, it));
+        // TODO should save also the information on where the multiplicity happened and its size
+        SetErrorMsg(kAttachMultipleCl);
+        // TODO should also compare with mean and sigma for this row
+        if(yres[jr][ic] > dyc[idx]) continue;
+      }
+
       // TODO proper indexing of clusters !!
-      fIndexes[it+kNtb*ir]  = chamber->GetTB(it)->GetGlobalIndex(idxs[jr][ic]);
-      fClusters[it+kNtb*ir] = c;
-  
-      //printf("\tid[%2d] it[%d] idx[%d]\n", ic, it, fIndexes[it]);
-  
+      fIndexes[idx]  = chamber->GetTB(it)->GetGlobalIndex(idxs[jr][ic]);
+      fClusters[idx] = c;
+      dyc[idx]        = yres[jr][ic];
       n++;
     }
-  }  
+  }
+  SetN(n);
 
   // number of minimum numbers of clusters expected for the tracklet
-  if (n < kClmin){
-    //AliWarning(Form("Not enough clusters to fit the tracklet %d [%d].", n, kClmin));
+  if (GetN() < kClmin){
+    AliDebug(1, Form("NOT ENOUGH CLUSTERS %d ATTACHED TO THE TRACKLET [min %d] FROM FOUND %d.", GetN(), kClmin, n));
+    SetErrorMsg(kAttachClAttach);
     return kFALSE;
   }
-  SetN(n);
 
   // Load calibration parameters for this tracklet  
   Calibrate();
@@ -1109,10 +1283,8 @@ void AliTRDseedV1::Bootstrap(const AliTRDReconstructor *rec)
 //
   fkReconstructor = rec;
   AliTRDgeometry g;
-  AliTRDpadPlane *pp = g.GetPadPlane(fDet);
-  fPad[0] = pp->GetLengthIPad();
-  fPad[1] = pp->GetWidthIPad();
-  fPad[3] = TMath::Tan(TMath::DegToRad()*pp->GetTiltingAngle());
+  SetPadPlane(g.GetPadPlane(fDet));
+
   //fSnp = fYref[1]/TMath::Sqrt(1+fYref[1]*fYref[1]);
   //fTgl = fZref[1];
   Int_t n = 0, nshare = 0, nused = 0;
@@ -1131,16 +1303,17 @@ void AliTRDseedV1::Bootstrap(const AliTRDReconstructor *rec)
 
 
 //____________________________________________________________________
-Bool_t AliTRDseedV1::Fit(Bool_t tilt, Bool_t zcorr)
+Bool_t AliTRDseedV1::Fit(UChar_t opt)
 {
 //
 // Linear fit of the clusters attached to the tracklet
 //
 // Parameters :
-//   - tilt : switch for tilt pad correction of cluster y position based on 
-//            the z, dzdx info from outside [default false].
-//   - zcorr : switch for using z information to correct for anisochronity 
-//            and a finner error parameterization estimation [default false]  
+//   - opt : switch for tilt pad correction of cluster y position. Options are
+//           0 no correction [default]
+//           1 full tilt correction [dz/dx and z0]
+//           2 pseudo tilt correction [dz/dx from pad-chamber geometry]
+//
 // Output :
 //  True if successful
 //
@@ -1225,33 +1398,63 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Bool_t zcorr)
 // Author 
 // A.Bercuci <A.Bercuci@gsi.de>
 
+  if(!fkReconstructor){
+    AliError("The tracklet needs the reconstruction setup. Please initialize by SetReconstructor().");
+    return kFALSE;
+  }
   if(!IsCalibrated()) Calibrate();
+  if(opt>2){
+    AliWarning(Form("Option [%d] outside range [0, 2]. Using default",opt));
+    opt=0;
+  }
 
   const Int_t kClmin = 8;
-
+  const Float_t kScalePulls = 10.; // factor to scale y pulls - NOT UNDERSTOOD
   // get track direction
   Double_t y0   = fYref[0];
   Double_t dydx = fYref[1]; 
   Double_t z0   = fZref[0];
   Double_t dzdx = fZref[1];
-  Double_t yt, zt;
 
-  //AliTRDtrackerV1::AliTRDLeastSquare fitterZ;
-  TLinearFitter& fitterY=*GetFitterY();
-  TLinearFitter& fitterZ=*GetFitterZ();
+  AliTRDtrackerV1::AliTRDLeastSquare fitterY;
+  AliTRDtrackerV1::AliTRDLeastSquare fitterZ;
 
   // book cluster information
   Double_t qc[kNclusters], xc[kNclusters], yc[kNclusters], zc[kNclusters], sy[kNclusters];
 
-  Int_t n = 0;
-  AliTRDcluster *c=NULL, **jc = &fClusters[0];
-  for (Int_t ic=0; ic<kNtb; ic++, ++jc) {
-    xc[ic]  = -1.;
-    yc[ic]  = 999.;
-    zc[ic]  = 999.;
-    sy[ic]  = 0.;
+  Bool_t tilt(opt==1)       // full tilt correction
+        ,pseudo(opt==2)     // pseudo tilt correction
+        ,rc(IsRowCross())   // row cross candidate 
+        ,kDZDX(IsPrimary());// switch dzdx calculation for barrel primary tracks
+  Int_t n(0);   // clusters used in fit 
+  AliTRDcluster *c(NULL), *cc(NULL), **jc = &fClusters[0];
+  const AliTRDrecoParam* const recoParam = fkReconstructor->GetRecoParam(); //the dynamic cast in GetRecoParam is slow, so caching the pointer to it
+
+  const Char_t *tcName[]={"NONE", "FULL", "HALF"};
+  AliDebug(2, Form("Options : TC[%s] dzdx[%c]", tcName[opt], kDZDX?'Y':'N'));
+
+  for (Int_t ic=0; ic<kNclusters; ic++, ++jc) {
+    xc[ic]  = -1.; yc[ic]  = 999.; zc[ic]  = 999.; sy[ic]  = 0.;
     if(!(c = (*jc))) continue;
     if(!c->IsInChamber()) continue;
+    // compute pseudo tilt correction
+    if(kDZDX){ 
+      fZfit[0] = c->GetZ();
+      if(rc){
+        for(Int_t kc=AliTRDseedV1::kNtb; kc<AliTRDseedV1::kNclusters; kc++){
+          if(!(cc=fClusters[kc])) continue;
+          if(!cc->IsInChamber()) continue;
+          fZfit[0] += cc->GetZ(); fZfit[0] *= 0.5;
+          break;
+        }
+      }
+      fZfit[1] = fZfit[0]/fX0;
+      if(rc){
+        fZfit[0] += fZfit[1]*0.5*AliTRDgeometry::CdrHght();
+        fZfit[1] = fZfit[0]/fX0;
+      }
+      kDZDX=kFALSE;
+    }
 
     Float_t w = 1.;
     if(c->GetNPads()>4) w = .5;
@@ -1261,48 +1464,64 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Bool_t zcorr)
     qc[n]   = TMath::Abs(c->GetQ());
     // pad row of leading 
 
-    // Radial cluster position
-    //Int_t jc = TMath::Max(fN-3, 0);
-    //xc[fN]   = c->GetXloc(fT0, fVD, &qc[jc], &xc[jc]/*, z0 - c->GetX()*dzdx*/);
     xc[n]   = fX0 - c->GetX();
 
-    // extrapolated track to cluster position
-    yt = y0 - xc[n]*dydx; 
-    zt = z0 - xc[n]*dzdx; 
-
     // Recalculate cluster error based on tracking information
-    c->SetSigmaY2(fS2PRF, fDiffT, fExB, xc[n], zcorr?zt:-1., dydx);
+    c->SetSigmaY2(fS2PRF, fDiffT, fExB, xc[n], -1./*zcorr?zt:-1.*/, dydx);
+    c->SetSigmaZ2(fPad[0]*fPad[0]/12.); // for HLT
     sy[n]  = TMath::Sqrt(c->GetSigmaY2());
 
-    yc[n]   = fkReconstructor->UseGAUS() ? 
+    yc[n]  = recoParam->UseGAUS() ? 
       c->GetYloc(y0, sy[n], GetPadWidth()): c->GetY();
     zc[n]   = c->GetZ();
-    //optional tilt correction
-    if(tilt) yc[n] -= (GetTilt()*(zc[n] - zt)); 
 
-    fitterY.AddPoint(&xc[n], yc[n], TMath::Sqrt(sy[n]));
-    fitterZ.AddPoint(&xc[n], qc[n], 1.);
+    //optional r-phi correction
+    //printf("   n[%2d] yc[%7.5f] ", n, yc[n]);
+    Float_t correction(0.);
+    if(tilt) correction = fPad[2]*(xc[n]*dzdx + zc[n] - z0);
+    else if(pseudo) correction = fPad[2]*(xc[n]*fZfit[1] + zc[n]-fZfit[0]);
+    yc[n]-=correction;
+    //printf("corr(%s%s)[%7.5f] yc1[%7.5f]\n", (tilt?"TC":""), (zcorr?"PC":""), correction, yc[n]);
+
+    AliDebug(5, Form("  tb[%2d] dx[%6.3f] y[%6.2f+-%6.3f]", c->GetLocalTimeBin(), xc[n], yc[n], sy[n]));
+    fitterY.AddPoint(&xc[n], yc[n], sy[n]);
+    if(rc) fitterZ.AddPoint(&xc[n], qc[n]*(ic<kNtb?1.:-1.), 1.);
     n++;
   }
-  // to few clusters
-  if (n < kClmin) return kFALSE; 
 
+  // to few clusters
+  if (n < kClmin){ 
+    AliDebug(1, Form("Not enough clusters to fit. Clusters: Attached[%d] Fit[%d].", GetN(), n));
+    SetErrorMsg(kFitCl);
+    return kFALSE; 
+  }
   // fit XY
-  fitterY.Eval();
-  fYfit[0] = fitterY.GetParameter(0);
-  fYfit[1] = -fitterY.GetParameter(1);
+  if(!fitterY.Eval()){
+    AliDebug(1, "Fit Y failed.");
+    SetErrorMsg(kFitFailedY);
+    return kFALSE;
+  }
+  fYfit[0] = fitterY.GetFunctionParameter(0);
+  fYfit[1] = -fitterY.GetFunctionParameter(1);
   // store covariance
-  Double_t *p = fitterY.GetCovarianceMatrix();
-  fCov[0] = p[0]; // variance of y0
-  fCov[1] = p[1]; // covariance of y0, dydx
-  fCov[2] = p[3]; // variance of dydx
+  Double_t p[3];
+  fitterY.GetCovarianceMatrix(p);
+  fCov[0] = kScalePulls*p[1]; // variance of y0
+  fCov[1] = kScalePulls*p[2]; // covariance of y0, dydx
+  fCov[2] = kScalePulls*p[0]; // variance of dydx
   // the ref radial position is set at the minimum of 
   // the y variance of the tracklet
   fX   = -fCov[1]/fCov[2];
+  fS2Y = fCov[0] +2.*fX*fCov[1] + fX*fX*fCov[2];
+
+  Float_t xs=fX+.5*AliTRDgeometry::CamHght();
+  if(xs < 0. || xs > AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght()){
+    AliDebug(1, Form("Ref radial position ouside chamber x[%5.2f].", fX));
+    SetErrorMsg(kFitFailedY);
+    return kFALSE;
+  }
 
-  // fit XZ
-  if(IsRowCross()){
-/*    // THE LEADING CLUSTER METHOD
+/*    // THE LEADING CLUSTER METHOD for z fit
     Float_t xMin = fX0;
     Int_t ic=n=kNclusters-1; jc = &fClusters[ic];
     AliTRDcluster *c0 =0x0, **kc = &fClusters[kNtb-1];
@@ -1318,37 +1537,26 @@ Bool_t AliTRDseedV1::Fit(Bool_t tilt, Bool_t zcorr)
     fS2Z     = fdX*fZref[1];
     fS2Z    *= fS2Z; fS2Z    *= 0.2887; //  1/sqrt(12)*/
 
-    // THE FIT X-Q PLANE METHOD 
-    Int_t ic=n=kNclusters-1; jc = &fClusters[ic];
-    for(; ic>kNtb; ic--, --jc){
-      if(!(c = (*jc))) continue;
-      if(!c->IsInChamber()) continue;
-      qc[n]   = TMath::Abs(c->GetQ());
-      xc[n]   = fX0 - c->GetX();
-      zc[n]   = c->GetZ();
-      fitterZ.AddPoint(&xc[n], -qc[n], 1.);
-      n--;
-    }
-    // fit XZ
-    fitterZ.Eval();
-    if(fitterZ.GetParameter(1)!=0.){ 
-      fX = -fitterZ.GetParameter(0)/fitterZ.GetParameter(1);
-      fX=(fX<0.)?0.:fX;
-      Float_t dl = .5*AliTRDgeometry::CamHght()+AliTRDgeometry::CdrHght();
-      fX=(fX> dl)?dl:fX;
-      fX-=.055; // TODO to be understood
+  // fit QZ
+  if(opt!=1 && IsRowCross()){
+    if(!fitterZ.Eval()) SetErrorMsg(kFitFailedZ);
+    if(!HasError(kFitFailedZ) && TMath::Abs(fitterZ.GetFunctionParameter(1))>1.e-10){ 
+      // TODO - one has to recalculate xy fit based on
+      // better knowledge of z position
+//       Double_t x = -fitterZ.GetFunctionParameter(0)/fitterZ.GetFunctionParameter(1);
+//       Double_t z0 = .5*(zc[0]+zc[n-1]);
+//       fZfit[0] = z0 + fZfit[1]*x; 
+//       fZfit[1] = fZfit[0]/fX0; 
+//       redo fit on xy plane
     }
-
-    fZfit[0] = .5*(zc[0]+zc[kNclusters-1]); fZfit[1] = 0.;
     // temporary external error parameterization
     fS2Z     = 0.05+0.4*TMath::Abs(fZref[1]); fS2Z *= fS2Z;
     // TODO correct formula
     //fS2Z     = sigma_x*TMath::Abs(fZref[1]);
   } else {
-    fZfit[0] = zc[0]; fZfit[1] = 0.;
+    //fZfit[0] = zc[0] + dzdx*0.5*AliTRDgeometry::CdrHght();
     fS2Z     = GetPadLength()*GetPadLength()/12.;
   }
-  fS2Y = fCov[0] +2.*fX*fCov[1] + fX*fX*fCov[2];
   return kTRUE;
 }
 
@@ -1615,13 +1823,17 @@ void AliTRDseedV1::Print(Option_t *o) const
   AliInfo(Form("Det[%3d] X0[%7.2f] Pad{L[%5.2f] W[%5.2f] Tilt[%+6.2f]}", fDet, fX0, GetPadLength(), GetPadWidth(), GetTilt()));
   AliInfo(Form("N[%2d] Nused[%2d] Nshared[%2d] [%d]", GetN(), GetNUsed(), GetNShared(), fN));
   AliInfo(Form("FLAGS : RC[%c] Kink[%c] SA[%c]", IsRowCross()?'y':'n', IsKink()?'y':'n', IsStandAlone()?'y':'n'));
+  AliInfo(Form("CALIB PARAMS :  T0[%5.2f]  Vd[%5.2f]  s2PRF[%5.2f]  ExB[%5.2f]  Dl[%5.2f]  Dt[%5.2f]", fT0, fVD, fS2PRF, fExB, fDiffL, fDiffT));
 
   Double_t cov[3], x=GetX();
   GetCovAt(x, cov);
   AliInfo("    |  x[cm]  |      y[cm]       |      z[cm]      |  dydx |  dzdx |");
   AliInfo(Form("Fit | %7.2f | %7.2f+-%7.2f | %7.2f+-%7.2f| %5.2f | ----- |", x, GetY(), TMath::Sqrt(cov[0]), GetZ(), TMath::Sqrt(cov[2]), fYfit[1]));
   AliInfo(Form("Ref | %7.2f | %7.2f+-%7.2f | %7.2f+-%7.2f| %5.2f | %5.2f |", x, fYref[0]-fX*fYref[1], TMath::Sqrt(fRefCov[0]), fZref[0]-fX*fYref[1], TMath::Sqrt(fRefCov[2]), fYref[1], fZref[1]))
-
+  AliInfo(Form("P / Pt [GeV/c] = %f / %f", GetMomentum(), fPt));
+  if(IsStandAlone()) AliInfo(Form("C Rieman / Vertex [1/cm] = %f / %f", fC[0], fC[1]));
+  AliInfo(Form("dEdx [a.u.]    = %f / %f / %f / %f / %f/ %f / %f / %f", fdEdx[0], fdEdx[1], fdEdx[2], fdEdx[3], fdEdx[4], fdEdx[5], fdEdx[6], fdEdx[7]));
+  AliInfo(Form("PID            = %5.3f / %5.3f / %5.3f / %5.3f / %5.3f", fProb[0], fProb[1], fProb[2], fProb[3], fProb[4]));
 
   if(strcmp(o, "a")!=0) return;
 
@@ -1648,9 +1860,9 @@ Bool_t AliTRDseedV1::IsEqual(const TObject *o) const
     if ( fZref[i] != inTracklet->fZref[i] ) return kFALSE;
   }
   
-  if ( fS2Y != inTracklet->fS2Y ) return kFALSE;
-  if ( GetTilt() != inTracklet->GetTilt() ) return kFALSE;
-  if ( GetPadLength() != inTracklet->GetPadLength() ) return kFALSE;
+  if ( TMath::Abs(fS2Y - inTracklet->fS2Y)>1.e-10 ) return kFALSE;
+  if ( TMath::Abs(GetTilt() - inTracklet->GetTilt())>1.e-10 ) return kFALSE;
+  if ( TMath::Abs(GetPadLength() - inTracklet->GetPadLength())>1.e-10 ) return kFALSE;
   
   for (Int_t i = 0; i < kNclusters; i++){
 //     if ( fX[i] != inTracklet->GetX(i) ) return kFALSE;
@@ -1673,14 +1885,14 @@ Bool_t AliTRDseedV1::IsEqual(const TObject *o) const
   //if ( fFreq != inTracklet->GetFreq() ) return kFALSE;
   //if ( fNChange != inTracklet->GetNChange() ) return kFALSE;
    
-  if ( fC != inTracklet->fC ) return kFALSE;
+  if ( TMath::Abs(fC[0] - inTracklet->fC[0])>1.e-10 ) return kFALSE;
   //if ( fCC != inTracklet->GetCC() ) return kFALSE;
-  if ( fChi2 != inTracklet->fChi2 ) return kFALSE;
+  if ( TMath::Abs(fChi2 - inTracklet->fChi2)>1.e-10 ) return kFALSE;
   //  if ( fChi2Z != inTracklet->GetChi2Z() ) return kFALSE;
 
   if ( fDet != inTracklet->fDet ) return kFALSE;
-  if ( fPt != inTracklet->fPt ) return kFALSE;
-  if ( fdX != inTracklet->fdX ) return kFALSE;
+  if ( TMath::Abs(fPt - inTracklet->fPt)>1.e-10 ) return kFALSE;
+  if ( TMath::Abs(fdX - inTracklet->fdX)>1.e-10 ) return kFALSE;
   
   for (Int_t iCluster = 0; iCluster < kNclusters; iCluster++){
     AliTRDcluster *curCluster = fClusters[iCluster];