]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDtrackV1.cxx
fix a memory leak coming from the TRD stand alone track finder (Markus)
[u/mrichter/AliRoot.git] / TRD / AliTRDtrackV1.cxx
index fca6528c32e7bc05b75aa6dec1727015d44b70e7..8d550fd2fb8113afb03213f8f03942f9b4fc0c1c 100644 (file)
 
 /* $Id$ */
 
+#include "AliESDtrack.h"
+#include "AliTracker.h"
+
 #include "AliTRDtrackV1.h"
 #include "AliTRDcluster.h"
 #include "AliTRDcalibDB.h"
+#include "AliTRDReconstructor.h"
 #include "AliTRDrecoParam.h"
 
-#include "AliESDtrack.h"
-
 ClassImp(AliTRDtrackV1)
 
 ///////////////////////////////////////////////////////////////////////////////
@@ -36,61 +38,103 @@ ClassImp(AliTRDtrackV1)
 ///////////////////////////////////////////////////////////////////////////////
 
 //_______________________________________________________________
-AliTRDtrackV1::AliTRDtrackV1() 
-  :AliTRDtrack()
-  ,fRecoParam(0x0)
+AliTRDtrackV1::AliTRDtrackV1() : AliKalmanTrack()
+  ,fPIDquality(0)
+  ,fDE(0.)
+  ,fBackupTrack(0x0)
 {
   //
   // Default constructor
   //
-       
-       for(int ip=0; ip<6; ip++){
-               fTrackletIndex[ip] = -1;
-               fTracklet[ip].Reset();
-       }
+  //printf("AliTRDtrackV1::AliTRDtrackV1()\n");
+
+  for(int i =0; i<3; i++) fBudget[i] = 0.;
+  
+  Float_t pid = 1./AliPID::kSPECIES;
+  for(int is =0; is<AliPID::kSPECIES; is++) fPID[is] = pid;
+
+  for(int ip=0; ip<kNplane; ip++){
+    fTrackletIndex[ip] = 0xffff;
+    fTracklet[ip]      = 0x0;
+  }
 }
 
 //_______________________________________________________________
-AliTRDtrackV1::AliTRDtrackV1(const AliESDtrack &t) 
-  :AliTRDtrack(t)
-  ,fRecoParam(0x0)
+AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &ref) : AliKalmanTrack(ref)
+  ,fPIDquality(ref.fPIDquality)
+  ,fDE(ref.fDE)
+  ,fBackupTrack(0x0)
 {
   //
-  // Constructor from AliESDtrack
+  // Copy constructor
   //
 
-       //AliInfo(Form("alpha %f", GetAlpha()));
-       t.GetTRDtracklets(&fTrackletIndex[0]);
-       for(int ip=0; ip<6; ip++) fTracklet[ip].Reset();
+  //printf("AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &)\n");
+  SetBit(kOwner, kFALSE);
+  for(int ip=0; ip<kNplane; ip++){ 
+    fTrackletIndex[ip] = ref.fTrackletIndex[ip];
+    fTracklet[ip]      = ref.fTracklet[ip];
+  }
+
+  for (Int_t i = 0; i < 3;i++) fBudget[i]      = ref.fBudget[i];
+
+       for(Int_t is = 0; is<AliPID::kSPECIES; is++) fPID[is] = ref.fPID[is];
+  
+  AliKalmanTrack::SetNumberOfClusters(ref.GetNumberOfClusters());
 }
 
 //_______________________________________________________________
-AliTRDtrackV1::AliTRDtrackV1(const AliTRDtrackV1 &ref) 
-  :AliTRDtrack(ref)
-  ,fRecoParam(ref.fRecoParam)
+AliTRDtrackV1::AliTRDtrackV1(const AliESDtrack &t) : AliKalmanTrack()
+  ,fPIDquality(0)
+  ,fDE(0.)
+  ,fBackupTrack(0x0)
 {
   //
-  // Copy constructor
+  // Constructor from AliESDtrack
   //
 
-       for(int ip=0; ip<6; ip++){ 
-               fTrackletIndex[ip] = ref.fTrackletIndex[ip];
-               fTracklet[ip]      = ref.fTracklet[ip];
-       }
+  SetLabel(t.GetLabel());
+  SetChi2(0.0);
+  SetMass(t.GetMass());
+  AliKalmanTrack::SetNumberOfClusters(t.GetTRDncls()); 
+  Int_t ti[kNplane]; t.GetTRDtracklets(&ti[0]);
+  for(int ip=0; ip<kNplane; ip++){ 
+    fTrackletIndex[ip] = ti[ip] < 0 ? 0xffff : ti[ip];
+    fTracklet[ip]      = 0x0;
+  }
+  for(int i =0; i<3; i++) fBudget[i] = 0.;
+  
+  Float_t pid = 1./AliPID::kSPECIES;
+  for(int is =0; is<AliPID::kSPECIES; is++) fPID[is] = pid;
+
+  const AliExternalTrackParam *par = &t;
+  if (t.GetStatus() & AliESDtrack::kTRDbackup) { 
+    par = t.GetOuterParam();
+    if (!par) {
+      AliError("No backup info!"); 
+      par = &t;
+    }
+  }
+  Set(par->GetX() 
+     ,par->GetAlpha()
+     ,par->GetParameter()
+     ,par->GetCovariance());
+
+  if(t.GetStatus() & AliESDtrack::kTIME) {
+    StartTimeIntegral();
+    Double_t times[10]; 
+    t.GetIntegratedTimes(times); 
+    SetIntegratedTimes(times);
+    SetIntegratedLength(t.GetIntegratedLength());
+  }
 }
 
 //_______________________________________________________________
-// AliTRDtrackV1::~AliTRDtrackV1()
-// {
-//     
-// }
-       
-//_______________________________________________________________
-AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 *trklts, const Double_t p[5]
-                           , const Double_t cov[15]
-                           , Double_t x, Double_t alpha)
-  :AliTRDtrack()
-  ,fRecoParam(0x0)
+AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 *trklts, const Double_t p[5], const Double_t cov[15]
+             , Double_t x, Double_t alpha) : AliKalmanTrack()
+  ,fPIDquality(0)
+  ,fDE(0.)
+  ,fBackupTrack(0x0)
 {
   //
   // The stand alone tracking constructor
@@ -100,7 +144,8 @@ AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 *trklts, const Double_t p[5]
   // 2. dQdl calculation
   //
 
-  Double_t cnv   = 1.0 / (GetBz() * kB2C);
+  Double_t cnv = GetBz() < 1.e-5 ? 1.e5 : 1.0 / (GetBz() * kB2C);
+  //  Double_t cnv   = 1.0 / (GetBz() * kB2C);
 
   Double_t pp[5] = { p[0]    
                    , p[1]
@@ -119,33 +164,99 @@ AliTRDtrackV1::AliTRDtrackV1(AliTRDseedV1 *trklts, const Double_t p[5]
                     , c20,         c21,         c22
                     , cov[ 6],     cov[ 7],     c32,     cov[ 9]
                     , cov[10]*cnv, cov[11]*cnv, c42*cnv, cov[13]*cnv, cov[14]*cnv*cnv };
-  
+
+  Double_t mostProbablePt=AliExternalTrackParam::GetMostProbablePt();
+  Double_t p0=TMath::Sign(1/mostProbablePt,pp[4]);
+  Double_t w0=cc[14]/(cc[14] + p0*p0), w1=p0*p0/(cc[14] + p0*p0);
+  pp[4] = w0*p0 + w1*pp[4];
+  cc[10]*=w1; cc[11]*=w1; cc[12]*=w1; cc[13]*=w1; cc[14]*=w1;
+
        Set(x,alpha,pp,cc);
-  Double_t s = GetSnp();
-  Double_t t = GetTgl();
-
-       Int_t ncl = 0, nclPlane; AliTRDcluster *c = 0x0;
-       for(int iplane=0; iplane<6; iplane++){
-               fTrackletIndex[iplane] = -1;
-               fTracklet[iplane] = trklts[iplane];
-               nclPlane = 0;
-               for(int ic = 0; ic<AliTRDseed::knTimebins; ic++){
-                       if(!fTracklet[iplane].IsUsable(ic)) continue;
-                       if(!(c = fTracklet[iplane].GetClusters(ic))) continue;
-                       
-                       fIndex[ncl]      = fTracklet[iplane].GetIndexes(ic);
-                       Double_t q = TMath::Abs(c->GetQ());
-                       fClusters[ncl]   = c;
-                       // temporary !!!
-                       // This is not correctly. Has to be updated in the AliTRDtrackerV1::FollowBackProlonagation()
-                       fdQdl[ncl]       = q * (s*s < 1.) ? TMath::Sqrt((1-s*s)/(1+t*t)) : 1.;
-                       ncl++;
-                       nclPlane++;
-               }
-               //printf("%d N clusters plane %d [%d %d].\n", iplane, nclPlane, fTracklet[iplane].GetN2(), trklts[iplane].GetN());
+  Int_t ncls = 0;
+       for(int iplane=0; iplane<kNplane; iplane++){
+    fTrackletIndex[iplane] = 0xffff;
+               if(!trklts[iplane].IsOK()) fTracklet[iplane] = 0x0;
+    else{ 
+      fTracklet[iplane] = new AliTRDseedV1(trklts[iplane]);
+      ncls += fTracklet[iplane]->GetN();
+    }
        }
-       //printf("N clusters in AliTRDtrackV1 %d.\n", ncl);
-       SetNumberOfClusters(ncl);
+  AliKalmanTrack::SetNumberOfClusters(ncls);           
+  for(int i =0; i<3; i++) fBudget[i] = 0.;
+  
+  Float_t pid = 1./AliPID::kSPECIES;
+  for(int is =0; is<AliPID::kSPECIES; is++) fPID[is] = pid;
+
+}
+
+//_______________________________________________________________
+AliTRDtrackV1::~AliTRDtrackV1()
+{
+  //AliInfo("");
+  //printf("I-AliTRDtrackV1::~AliTRDtrackV1() : Owner[%s]\n", TestBit(kOwner)?"YES":"NO");
+  if(fBackupTrack) {
+    delete fBackupTrack;
+  }
+  fBackupTrack = 0x0;
+
+  if(TestBit(kOwner)){
+    for(Int_t ip=0; ip<kNplane; ip++){
+      if(fTracklet[ip]) delete fTracklet[ip];
+      fTracklet[ip] = 0x0;
+    }
+  }
+}
+       
+//_______________________________________________________________
+Bool_t AliTRDtrackV1::CookLabel(Float_t wrong)
+{
+  // set MC label for this track
+  
+  Int_t s[kMAXCLUSTERSPERTRACK][2];
+  for (Int_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
+    s[i][0] = -1;
+    s[i][1] =  0;
+  }
+  
+  Bool_t labelAdded;
+  Int_t label;
+  AliTRDcluster *c    = 0x0;
+  for (Int_t ip = 0; ip < kNplane; ip++) {
+    if(fTrackletIndex[ip] == 0xffff) continue;
+    for (Int_t ic = 0; ic < AliTRDseed::knTimebins; ic++) {
+      if(!(c = fTracklet[ip]->GetClusters(ic))) continue;
+      for (Int_t k = 0; k < 3; k++) { 
+        label      = c->GetLabel(k);
+        labelAdded = kFALSE; 
+        Int_t j = 0;
+        if (label >= 0) {
+          while ((!labelAdded) && (j < kMAXCLUSTERSPERTRACK)) {
+            if ((s[j][0] == label) || 
+                (s[j][1] ==     0)) {
+              s[j][0] = label; 
+              s[j][1]++; 
+              labelAdded = kTRUE;
+            }
+            j++;
+          }
+        }
+      }
+    }
+  }
+  
+  Int_t max = 0;
+  label = -123456789;
+  for (Int_t i = 0; i < kMAXCLUSTERSPERTRACK; i++) {
+    if (s[i][1] <= max) continue;
+    max   = s[i][1]; 
+    label = s[i][0];
+  }
+  
+  if ((1. - Float_t(max)/GetNumberOfClusters()) > wrong) label = -label;
+  
+  SetLabel(label); 
+  
+  return kTRUE;
 }
 
 //_______________________________________________________________
@@ -154,66 +265,85 @@ Bool_t AliTRDtrackV1::CookPID()
   //
   // Cook the PID information
   //
-
-// CookdEdx();  // truncated mean ... do we still need it ?
-
-// CookdEdxTimBin(seed->GetID());
-       
-  // Sets the a priori probabilities
+  
+  // Reset the a priori probabilities
+  Double_t pid = 1. / AliPID::kSPECIES;
   for(int ispec=0; ispec<AliPID::kSPECIES; ispec++) {
-    fPID[ispec] = 1.0 / AliPID::kSPECIES;      
+    fPID[ispec] = pid; 
+  }
+  fPIDquality = 0;
+  
+  // steer PID calculation @ tracklet level
+  Double_t *prob = 0x0;
+  for(int ip=0; ip<kNplane; ip++){
+    if(fTrackletIndex[ip] == 0xffff) continue;
+    if(!fTracklet[ip]->IsOK()) continue;
+    if(!(prob = fTracklet[ip]->GetProbability())) return kFALSE;
+    
+    Int_t nspec = 0; // quality check of tracklet dEdx
+    for(int ispec=0; ispec<AliPID::kSPECIES; ispec++){
+      if(prob[ispec] < 0.) continue;
+      fPID[ispec] *= prob[ispec];
+      nspec++;
+    }
+    if(!nspec) continue;
+    
+    fPIDquality++;
   }
-       
-       // steer PID calculation @ tracklet level
-       Double_t *prob = 0x0;
-       fPIDquality = 0;
-       for(int itrklt=0; itrklt<AliESDtrack::kNPlane; itrklt++){
-    //for (Int_t iSlice = 0; iSlice < AliESDtrack::kNSlice; iSlice++) fdEdxPlane[itrklt][iSlice] = -1.;
-
-               if(fTrackletIndex[itrklt]<0) continue;
-               if(!fTracklet[itrklt].IsOK()) continue;
-               if(!(prob = fTracklet[itrklt].GetProbability())) return kFALSE;
-               
-               Int_t nspec = 0; // quality check of tracklet dEdx
-               for(int ispec=0; ispec<AliPID::kSPECIES; ispec++){
-                       if(prob[ispec] < 0.) continue;
-                       fPID[ispec] *= prob[ispec];
-                       nspec++;
-               }
-               if(!nspec) continue;
-               
-               fPIDquality++;
-       }
   
   // no tracklet found for PID calculations
   if(!fPIDquality) return kTRUE;
-
-       // slot for PID calculation @ track level
-       
-
+  
+  // slot for PID calculation @ track level
+  
+  
   // normalize probabilities
   Double_t probTotal = 0.0;
   for (Int_t is = 0; is < AliPID::kSPECIES; is++) probTotal += fPID[is];
-
-
+  
+  
   if (probTotal <= 0.0) {
     AliWarning("The total probability over all species <= 0. This may be caused by some error in the reference data.");
     return kFALSE;
   }
-
+  
   for (Int_t iSpecies = 0; iSpecies < AliPID::kSPECIES; iSpecies++) fPID[iSpecies] /= probTotal;
-
+  
   return kTRUE;
 }
 
-//_______________________________________________________________
-Float_t AliTRDtrackV1::GetMomentum(Int_t plane) const
+//_____________________________________________________________________________
+Double_t AliTRDtrackV1::GetBz() const 
 {
   //
-  // Get the momentum at a given plane
+  // Returns Bz component of the magnetic field (kG)
   //
 
-       return plane >=0 && plane < 6 && fTrackletIndex[plane] >= 0 ? fTracklet[plane].GetMomentum() : -1.;
+  if (AliTracker::UniformField()) return AliTracker::GetBz();
+
+  Double_t r[3]; 
+  GetXYZ(r);
+  return AliTracker::GetBz(r);
+}
+
+//_______________________________________________________________
+Int_t  AliTRDtrackV1::GetClusterIndex(Int_t id) const
+{
+  Int_t n = 0;
+  for(Int_t ip=0; ip<kNplane; ip++){
+    if(!fTracklet[ip]) continue;
+    if(n+fTracklet[ip]->GetN() < id){ 
+      n+=fTracklet[ip]->GetN();
+      continue;
+    }
+    for(Int_t ic=34; ic>=0; ic--){
+      if(!fTracklet[ip]->GetClusters(ic)) continue;
+      n++;
+      if(n<id) continue;
+      return fTracklet[ip]->GetIndexes(ic);
+    }
+  }
+  return -1;
 }
 
 //_______________________________________________________________
@@ -222,108 +352,308 @@ Double_t AliTRDtrackV1::GetPredictedChi2(const AliTRDseedV1 *trklt) const
   //
   // Get the predicted chi2
   //
-
+  
   Double_t x      = trklt->GetX0();
   Double_t p[2]   = { trklt->GetYat(x)
                     , trklt->GetZat(x) };
   Double_t cov[3];
-       trklt->GetCovAt(x, cov);
-
+  trklt->GetCovAt(x, cov);
+  
   return AliExternalTrackParam::GetPredictedChi2(p, cov);
+}
+
+       
+//_____________________________________________________________________________
+void AliTRDtrackV1::MakeBackupTrack()
+{
+  //
+  // Creates a backup track
+  //
+
+  if(fBackupTrack) {
+    fBackupTrack->~AliTRDtrackV1();
+    new(fBackupTrack) AliTRDtrackV1((AliTRDtrackV1&)(*this));
+    return;
+  }
+  fBackupTrack = new AliTRDtrackV1((AliTRDtrackV1&)(*this));
+}
+
+//_____________________________________________________________________________
+Int_t AliTRDtrackV1::GetProlongation(Double_t xk, Double_t &y, Double_t &z)
+{
+  //
+  // Find a prolongation at given x
+  // Return 0 if it does not exist
+  //  
+
+  Double_t bz = GetBz();
+  if (!AliExternalTrackParam::GetYAt(xk,bz,y)) return 0;
+  if (!AliExternalTrackParam::GetZAt(xk,bz,z)) return 0;
+
+  return 1;  
 
 }
 
-//_______________________________________________________________
-Bool_t AliTRDtrackV1::IsOwner() const
+//_____________________________________________________________________________
+Bool_t AliTRDtrackV1::PropagateTo(Double_t xk, Double_t xx0, Double_t xrho)
 {
   //
-  // Check whether track owns the tracklets
+  // Propagates this track to a reference plane defined by "xk" [cm] 
+  // correcting for the mean crossed material.
   //
+  // "xx0"  - thickness/rad.length [units of the radiation length] 
+  // "xrho" - thickness*density    [g/cm^2] 
+  // 
 
-       for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
-               if(fTrackletIndex[ip] < 0) continue;
-               if(!fTracklet[ip].IsOwner()) return kFALSE;
-       }
-       return kTRUE;
+  if (xk == GetX()) {
+    return kTRUE;
+  }
+
+  Double_t oldX = GetX();
+  Double_t oldY = GetY();
+  Double_t oldZ = GetZ();
+
+  Double_t bz   = GetBz();
+
+  if (!AliExternalTrackParam::PropagateTo(xk,bz)) {
+    return kFALSE;
+  }
+
+  Double_t x = GetX();
+  Double_t y = GetY();
+  Double_t z = GetZ();
+
+  if (oldX < xk) {
+    xrho = -xrho;
+    if (IsStartedTimeIntegral()) {
+      Double_t l2  = TMath::Sqrt((x-oldX)*(x-oldX) 
+                               + (y-oldY)*(y-oldY) 
+                               + (z-oldZ)*(z-oldZ));
+      Double_t crv = AliExternalTrackParam::GetC(bz);
+      if (TMath::Abs(l2*crv) > 0.0001) {
+        // Make correction for curvature if neccesary
+        l2 = 0.5 * TMath::Sqrt((x-oldX)*(x-oldX) 
+                             + (y-oldY)*(y-oldY));
+        l2 = 2.0 * TMath::ASin(l2 * crv) / crv;
+        l2 = TMath::Sqrt(l2*l2 + (z-oldZ)*(z-oldZ));
+      }
+      AddTimeStep(l2);
+    }
+  }
+
+  if (!AliExternalTrackParam::CorrectForMeanMaterial(xx0,xrho,GetMass())) { 
+    return kFALSE;
+  }
+
+  {
+
+    // Energy losses
+    Double_t p2    = (1.0 + GetTgl()*GetTgl()) / (GetSigned1Pt()*GetSigned1Pt());
+    Double_t beta2 = p2 / (p2 + GetMass()*GetMass());
+    if ((beta2 < 1.0e-10) || 
+        ((5940.0 * beta2/(1.0 - beta2 + 1.0e-10) - beta2) < 0.0)) {
+      return kFALSE;
+    }
+
+    Double_t dE    = 0.153e-3 / beta2 
+                   * (TMath::Log(5940.0 * beta2/(1.0 - beta2 + 1.0e-10)) - beta2)
+                   * xrho;
+    fBudget[0] += xrho;
+
+    /*
+    // Suspicious part - think about it ?
+    Double_t kinE =  TMath::Sqrt(p2);
+    if (dE > 0.8*kinE) dE = 0.8 * kinE;  //      
+    if (dE < 0)        dE = 0.0;         // Not valid region for Bethe bloch 
+    */
+    fDE += dE;
+
+    /*
+    // Suspicious ! I.B.
+    Double_t sigmade = 0.07 * TMath::Sqrt(TMath::Abs(dE));   // Energy loss fluctuation 
+    Double_t sigmac2 = sigmade*sigmade*fC*fC*(p2+GetMass()*GetMass())/(p2*p2);
+    fCcc += sigmac2;
+    fCee += fX*fX * sigmac2;  
+    */
+
+  }
+
+  return kTRUE;
+}
+
+//_____________________________________________________________________________
+Int_t   AliTRDtrackV1::PropagateToR(Double_t r,Double_t step)
+{
+  //
+  // Propagate track to the radial position
+  // Rotation always connected to the last track position
+  //
+
+  Double_t xyz0[3];
+  Double_t xyz1[3];
+  Double_t y;
+  Double_t z; 
+
+  Double_t radius = TMath::Sqrt(GetX()*GetX() + GetY()*GetY());
+  // Direction +-
+  Double_t dir    = (radius > r) ? -1.0 : 1.0;   
+
+  for (Double_t x = radius+dir*step; dir*x < dir*r; x += dir*step) {
+
+    GetXYZ(xyz0);      
+    Double_t alpha = TMath::ATan2(xyz0[1],xyz0[0]);
+    Rotate(alpha,kTRUE);
+    GetXYZ(xyz0);      
+    GetProlongation(x,y,z);
+    xyz1[0] = x * TMath::Cos(alpha) + y * TMath::Sin(alpha); 
+    xyz1[1] = x * TMath::Sin(alpha) - y * TMath::Cos(alpha);
+    xyz1[2] = z;
+    Double_t param[7];
+    AliTracker::MeanMaterialBudget(xyz0,xyz1,param);
+    if (param[1] <= 0) {
+      param[1] = 100000000;
+    }
+    PropagateTo(x,param[1],param[0]*param[4]);
+
+  } 
+
+  GetXYZ(xyz0);        
+  Double_t alpha = TMath::ATan2(xyz0[1],xyz0[0]);
+  Rotate(alpha,kTRUE);
+  GetXYZ(xyz0);        
+  GetProlongation(r,y,z);
+  xyz1[0] = r * TMath::Cos(alpha) + y * TMath::Sin(alpha); 
+  xyz1[1] = r * TMath::Sin(alpha) - y * TMath::Cos(alpha);
+  xyz1[2] = z;
+  Double_t param[7];
+  AliTracker::MeanMaterialBudget(xyz0,xyz1,param);
+
+  if (param[1] <= 0) {
+    param[1] = 100000000;
+  }
+  PropagateTo(r,param[1],param[0]*param[4]);
+
+  return 0;
+
+}
+
+//_____________________________________________________________________________
+Bool_t AliTRDtrackV1::Rotate(Double_t alpha, Bool_t absolute)
+{
+  //
+  // Rotates track parameters in R*phi plane
+  // if absolute rotation alpha is in global system
+  // otherwise alpha rotation is relative to the current rotation angle
+  //  
+
+  if (absolute) alpha -= GetAlpha();
+  //else fNRotate++;
+
+  return AliExternalTrackParam::Rotate(GetAlpha()+alpha);
 }
+
+//___________________________________________________________
+void AliTRDtrackV1::SetNumberOfClusters() 
+{
+// Calculate the number of clusters attached to this track
+       
+  Int_t ncls = 0;
+  for(int ip=0; ip<kNplane; ip++){
+    if(fTracklet[ip] && fTrackletIndex[ip] != 0xffff) ncls += fTracklet[ip]->GetN();
+  }
+  AliKalmanTrack::SetNumberOfClusters(ncls);   
+}
+
        
 //_______________________________________________________________
-void AliTRDtrackV1::SetOwner(Bool_t own)
+void AliTRDtrackV1::SetOwner()
 {
   //
   // Toggle ownership of tracklets
   //
 
-       for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
-               if(fTrackletIndex[ip] < 0) continue;
-               //AliInfo(Form("p[%d] index[%d]", ip, fTrackletIndex[ip]));
-               fTracklet[ip].SetOwner(own);
-       }
+  if(TestBit(kOwner)) return;
+  for (Int_t ip = 0; ip < kNplane; ip++) {
+    if(fTrackletIndex[ip] == 0xffff) continue;
+    fTracklet[ip] = new AliTRDseedV1(*fTracklet[ip]);
+    fTracklet[ip]->SetOwner();
+  }
+  SetBit(kOwner);
 }
 
 //_______________________________________________________________
-void AliTRDtrackV1::SetTracklet(AliTRDseedV1 *trklt, Int_t plane, Int_t index)
+void AliTRDtrackV1::SetTracklet(AliTRDseedV1 *trklt, Int_t index)
 {
   //
   // Set the tracklets
   //
+  Int_t plane = trklt->GetPlane();
+  if(fTrackletIndex[plane]==0xffff && fTracklet[plane]){ 
+    delete fTracklet[plane];
+  }
+  fTracklet[plane]      = trklt;
+  fTrackletIndex[plane] = index;
+}
 
-       if(plane < 0 || plane >=6) return;
-       fTracklet[plane]      = (*trklt);
-       fTrackletIndex[plane] = index;
+//_______________________________________________________________
+void AliTRDtrackV1::UnsetTracklet(Int_t plane)
+{
+  if(plane<0 && plane >= kNplane) return;
+  fTrackletIndex[plane] = 0xffff;
+  fTracklet[plane] = 0x0;
 }
 
+
 //_______________________________________________________________
 Bool_t  AliTRDtrackV1::Update(AliTRDseedV1 *trklt, Double_t chisq)
 {
   //
-  // Update track parameters
+  // Update track and tracklet parameters 
   //
-
+  
   Double_t x      = trklt->GetX0();
   Double_t p[2]   = { trklt->GetYat(x)
                     , trklt->GetZat(x) };
   Double_t cov[3];
-       trklt->GetCovAt(x, cov);
-       
-//     Print();
-//     AliInfo(Form("cov[%f %f %f]", cov[0], cov[1], cov[2]));
-
+  trklt->GetCovAt(x, cov);
+  
   if(!AliExternalTrackParam::Update(p, cov)) return kFALSE;
-       //Print();
-
+  
+  AliTRDcluster *c = 0x0;
+  Int_t ic = 0; while(!(c = trklt->GetClusters(ic))) ic++;
+  AliTracker::FillResiduals(this, p, cov, c->GetVolumeId());
+  
+  
   // Register info to track
-//   Int_t n      = GetNumberOfClusters();
-//   fIndex[n]    = index;
-//   fClusters[n] = c;
-  SetNumberOfClusters(GetNumberOfClusters()+trklt->GetN());
+  SetNumberOfClusters();
   SetChi2(GetChi2() + chisq);
-       
-       // update tracklet
-       trklt->SetMomentum(GetP());
-  Double_t s = GetSnp(), t = GetTgl();
-       trklt->SetdQdl(TMath::Sqrt((1.0 - s*s) / (1.0 + t*t)));
-       return kTRUE;
+  
+  // update tracklet
+  trklt->SetMomentum(GetP());
+  trklt->SetSnp(GetSnp());
+  trklt->SetTgl(GetTgl());
+  return kTRUE;
 }
 
 //_______________________________________________________________
 void AliTRDtrackV1::UpdateESDtrack(AliESDtrack *track)
 {
   //
-  // Update the ESD track
+  // Update the TRD PID information in the ESD track
   //
-       
-       // copy dEdx to ESD
-       Float_t *dedx = 0x0;
-       for (Int_t ip = 0; ip < AliESDtrack::kNPlane; ip++) {
-               if(fTrackletIndex[ip] < 0) continue;
-               fTracklet[ip].CookdEdx(AliESDtrack::kNSlice);
-               dedx = fTracklet[ip].GetdEdx();
-               for (Int_t js = 0; js < AliESDtrack::kNSlice; js++) track->SetTRDsignals(dedx[js], ip, js);
-               //track->SetTRDTimBin(fTimBinPlane[i], i);
-       }
 
-       // copy PID to ESD
-       track->SetTRDpid(fPID);
-       track->SetTRDpidQuality(fPIDquality);
+  track->SetNumberOfTRDslices(kNslice);
+
+  for (Int_t ip = 0; ip < kNplane; ip++) {
+    if(fTrackletIndex[ip] == 0xffff) continue;
+    fTracklet[ip]->CookdEdx(kNslice);
+    Float_t *dedx = fTracklet[ip]->GetdEdx();
+    for (Int_t js = 0; js < kNslice; js++) track->SetTRDslice(dedx[js], ip, js);
+  }
+
+  // copy PID to ESD
+  if(!fPIDquality) return;
+  track->SetTRDpid(fPID);
+  track->SetTRDpidQuality(fPIDquality);
 }