]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCtracker.cxx
Transition to NewIO
[u/mrichter/AliRoot.git] / TPC / AliTPCtracker.cxx
index bb9ebc648129ba43ffc05564d259f3e29ca5a61c..dee9237690bf5426c5a9b227b000e547d02f8ca0 100644 (file)
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-/*
-$Log$
-Revision 1.11  2001/05/23 08:50:10  hristov
-Weird inline removed
-
-Revision 1.10  2001/05/16 14:57:25  alibrary
-New files for folders and Stack
-
-Revision 1.9  2001/05/11 07:16:56  hristov
-Fix needed on Sun and Alpha
-
-Revision 1.8  2001/05/08 15:00:15  hristov
-Corrections for tracking in arbitrary magnenetic field. Changes towards a concept of global Alice track. Back propagation of reconstructed tracks (Yu.Belikov)
-
-Revision 1.5  2000/12/20 07:51:59  kowal2
-Changes suggested by Alessandra and Paolo to avoid overlapped
-data fields in encapsulated classes.
-
-Revision 1.4  2000/11/02 07:27:16  kowal2
-code corrections
-
-Revision 1.2  2000/06/30 12:07:50  kowal2
-Updated from the TPC-PreRelease branch
-
-Revision 1.1.2.1  2000/06/25 08:53:55  kowal2
-Splitted from AliTPCtracking
-
-*/
+/* $Id$ */
 
 //-------------------------------------------------------
 //          Implementation of the TPC tracker
@@ -51,51 +24,64 @@ Splitted from AliTPCtracking
 #include <TObjArray.h>
 #include <TFile.h>
 #include <TTree.h>
-#include <iostream.h>
+#include <AliRunLoader.h>
+#include <AliLoader.h>
+#include <Riostream.h>
 
 #include "AliTPCtracker.h"
 #include "AliTPCcluster.h"
 #include "AliTPCParam.h"
 #include "AliTPCClustersRow.h"
-
-
-AliTPCtracker::AliTPCtracker(const AliTPCParam *par)
-{;
-//MI change only provisore - need change in the ITS code which depend on it
-}
-
+#include "AliTPCcluster.h"
 
 //_____________________________________________________________________________
-AliTPCtracker::AliTPCtracker(const AliTPCParam *par, Int_t eventn): 
-fkNIS(par->GetNInnerSector()/2), 
-fkNOS(par->GetNOuterSector()/2)
+AliTPCtracker::AliTPCtracker(const AliTPCParam *par, Int_t eventn, const char* evfoldname):
+AliTracker(), fkNIS(par->GetNInnerSector()/2), fkNOS(par->GetNOuterSector()/2),
+fEvFolderName(evfoldname)
 {
   //---------------------------------------------------------------------
   // The main TPC tracker constructor
   //---------------------------------------------------------------------
+  cout<<"fkNIS = "<<fkNIS<<endl;
+  cout<<"fkNOS = "<<fkNOS<<endl;
+
   fInnerSec=new AliTPCSector[fkNIS];         
   fOuterSec=new AliTPCSector[fkNOS];
 
   Int_t i;
+  
   for (i=0; i<fkNIS; i++) fInnerSec[i].Setup(par,0);
   for (i=0; i<fkNOS; i++) fOuterSec[i].Setup(par,1);
 
   fN=0;  fSectors=0;
 
+
+  fEventN = eventn;
+  AliRunLoader* rl = AliRunLoader::GetRunLoader(fEvFolderName);
+  if (rl == 0x0)
+   {
+     Error("AliTPCtracker","Can not get RL from specified folder %s",fEvFolderName.Data());
+     return;
+   }
+  rl->GetEvent(fEventN);
+
+  AliLoader* tpcl = rl->GetLoader("TPCLoader");
+  if (tpcl == 0x0)
+   {
+     Error("AliTPCtracker","Can not get TPC Laoder from Run Loader");
+     return;
+   }
+  
+  if (tpcl->TreeR() == 0x0) tpcl->LoadRecPoints("read");
+  TTree* treeR = tpcl->TreeR();
+  if (treeR == 0x0)
+   {
+     cout<<"Error: Can not get TreeR\n";
+   }
   fClustersArray.Setup(par);
   fClustersArray.SetClusterType("AliTPCcluster");
+  fClustersArray.ConnectTree(treeR);
 
-  char   cname[100];
-  if (eventn==-1) {
-    sprintf(cname,"TreeC_TPC");
-  }
-  else {
-    sprintf(cname,"TreeC_TPC_%d",eventn);
-  }
-
-  fClustersArray.ConnectTree(cname);
-
-  fEventN = eventn;
   fSeeds=0;
 }
 
@@ -106,7 +92,10 @@ AliTPCtracker::~AliTPCtracker() {
   //------------------------------------------------------------------
   delete[] fInnerSec;
   delete[] fOuterSec;
-  delete fSeeds;
+  if (fSeeds) {
+    fSeeds->Delete(); 
+    delete fSeeds;
+  }
 }
 
 //_____________________________________________________________________________
@@ -207,7 +196,14 @@ void AliTPCtracker::LoadOuterSectors() {
   // This function fills outer TPC sectors with clusters.
   //-----------------------------------------------------------------
   UInt_t index;
-  Int_t j=Int_t(fClustersArray.GetTree()->GetEntries());
+  TTree* tree = fClustersArray.GetTree();
+  if (tree == 0x0)
+   {
+    Error("LoadOuterSectors","Can not get tree from fClustersArray");
+    return;
+   }
+  Int_t j=Int_t(tree->GetEntries());
+  cout<<"fClustersArray.GetTree()->GetEntries() = "<<j<<endl;
   for (Int_t i=0; i<j; i++) {
       AliSegmentID *s=fClustersArray.LoadEntry(i);
       Int_t sec,row;
@@ -304,7 +300,8 @@ Int_t AliTPCtracker::FollowProlongation(AliTPCseed& t, Int_t rf) {
   if (alpha < 0.            ) alpha += 2.*TMath::Pi();  
   Int_t s=Int_t(alpha/fSectors->GetAlpha())%fN;
 
-  for (Int_t nr=fSectors->GetRowNumber(xt)-1; nr>=rf; nr--) {
+  Int_t nrows=fSectors->GetRowNumber(xt)-1;
+  for (Int_t nr=nrows; nr>=rf; nr--) {
     Double_t x=fSectors->GetX(nr), ymax=fSectors->GetMaxY(nr);
     if (!t.PropagateTo(x)) return 0;
 
@@ -326,20 +323,21 @@ Int_t AliTPCtracker::FollowProlongation(AliTPCseed& t, Int_t rf) {
 
     if (krow) {
       for (Int_t i=krow.Find(y-road); i<krow; i++) {
-       AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
-       if (c->GetY() > y+road) break;
-       if (c->IsUsed()) continue;
-       if ((c->GetZ()-z)*(c->GetZ()-z) > 16.*(t.GetSigmaZ2()+sz2)) continue;
-       Double_t chi2=t.GetPredictedChi2(c);
-       if (chi2 > maxchi2) continue;
-       maxchi2=chi2;
-       cl=c;
+       AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
+       if (c->GetY() > y+road) break;
+       if (c->IsUsed()) continue;
+       if ((c->GetZ()-z)*(c->GetZ()-z) > 16.*(t.GetSigmaZ2()+sz2)) continue;
+       Double_t chi2=t.GetPredictedChi2(c);
+       if (chi2 > maxchi2) continue;
+       maxchi2=chi2;
+       cl=c;
         index=krow.GetIndex(i);       
       }
     }
     if (cl) {
       Float_t l=fSectors->GetPadPitchWidth();
-      t.SetSampledEdx(cl->GetQ()/l,t.GetNumberOfClusters());
+      Float_t corr=1.; if (nr>63) corr=0.67; // new (third) pad response !
+      t.SetSampledEdx(cl->GetQ()/l*corr,t.GetNumberOfClusters());
       if (!t.Update(cl,maxchi2,index)) {
          if (!tryAgain--) return 0;
       } else tryAgain=kSKIP;
@@ -431,14 +429,14 @@ Int_t AliTPCtracker::FollowBackProlongation
        if (accepted>27)
        if (krow) {
           for (Int_t i=krow.Find(y-road); i<krow; i++) {
-           AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
-           if (c->GetY() > y+road) break;
-           if (c->IsUsed()) continue;
-        if ((c->GetZ()-z)*(c->GetZ()-z)>16.*(seed.GetSigmaZ2()+sz2)) continue;
-           Double_t chi2=seed.GetPredictedChi2(c);
-           if (chi2 > maxchi2) continue;
-           maxchi2=chi2;
-           cl=c;
+           AliTPCcluster *c=(AliTPCcluster*)(krow[i]);
+           if (c->GetY() > y+road) break;
+           if (c->IsUsed()) continue;
+        if ((c->GetZ()-z)*(c->GetZ()-z)>16.*(seed.GetSigmaZ2()+sz2)) continue;
+           Double_t chi2=seed.GetPredictedChi2(c);
+           if (chi2 > maxchi2) continue;
+           maxchi2=chi2;
+           cl=c;
             index=krow.GetIndex(i);
           }
        }
@@ -446,7 +444,8 @@ Int_t AliTPCtracker::FollowBackProlongation
 
     if (cl) {
       Float_t l=fSectors->GetPadPitchWidth();
-      seed.SetSampledEdx(cl->GetQ()/l,seed.GetNumberOfClusters());
+      Float_t corr=1.; if (i>63) corr=0.67; // new (third) pad response !
+      seed.SetSampledEdx(cl->GetQ()/l*corr,seed.GetNumberOfClusters());
       seed.Update(cl,maxchi2,index);
     }
 
@@ -462,6 +461,7 @@ void AliTPCtracker::MakeSeeds(Int_t i1, Int_t i2) {
   //-----------------------------------------------------------------
   // This function creates track seeds.
   //-----------------------------------------------------------------
+  cout<<" Making Seeds i1="<<i1<<"  i2="<<i2<<"\n";
   if (fSeeds==0) fSeeds=new TObjArray(15000);
 
   Double_t x[5], c[15];
@@ -480,61 +480,62 @@ void AliTPCtracker::MakeSeeds(Int_t i1, Int_t i2) {
     for (Int_t is=0; is < kr1; is++) {
       Double_t y1=kr1[is]->GetY(), z1=kr1[is]->GetZ();
       for (Int_t js=0; js < nl+nm+nu; js++) {
-       const AliTPCcluster *kcl;
+       const AliTPCcluster *kcl;
         Double_t x2,   y2,   z2;
-        Double_t x3=0.,y3=0.;
+        Double_t x3=GetX(), y3=GetY(), z3=GetZ();
 
-       if (js<nl) {
-         const AliTPCRow& kr2=fOuterSec[(ns-1+fkNOS)%fkNOS][i2];
-         kcl=kr2[js];
+       if (js<nl) {
+         const AliTPCRow& kr2=fOuterSec[(ns-1+fkNOS)%fkNOS][i2];
+         kcl=kr2[js];
           y2=kcl->GetY(); z2=kcl->GetZ();
           x2= xx2*cs+y2*sn;
           y2=-xx2*sn+y2*cs;
-       } else 
-         if (js<nl+nm) {
-           const AliTPCRow& kr2=fOuterSec[ns][i2];
-           kcl=kr2[js-nl];
+       } else 
+         if (js<nl+nm) {
+           const AliTPCRow& kr2=fOuterSec[ns][i2];
+           kcl=kr2[js-nl];
             x2=xx2; y2=kcl->GetY(); z2=kcl->GetZ();
-         } else {
-           const AliTPCRow& kr2=fOuterSec[(ns+1)%fkNOS][i2];
-           kcl=kr2[js-nl-nm];
+         } else {
+           const AliTPCRow& kr2=fOuterSec[(ns+1)%fkNOS][i2];
+           kcl=kr2[js-nl-nm];
             y2=kcl->GetY(); z2=kcl->GetZ();
             x2=xx2*cs-y2*sn;
             y2=xx2*sn+y2*cs;
-         }
+         }
 
-        Double_t zz=z1 - z1/x1*(x1-x2); 
+        Double_t zz=z1 - (z1-z3)/(x1-x3)*(x1-x2); 
         if (TMath::Abs(zz-z2)>5.) continue;
 
         Double_t d=(x2-x1)*(0.-y2)-(0.-x2)*(y2-y1);
         if (d==0.) {cerr<<"MakeSeeds warning: Straight seed !\n"; continue;}
 
-       x[0]=y1;
-       x[1]=z1;
-       x[4]=f1(x1,y1,x2,y2,x3,y3);
-       if (TMath::Abs(x[4]) >= 0.0066) continue;
-       x[2]=f2(x1,y1,x2,y2,x3,y3);
-       //if (TMath::Abs(x[4]*x1-x[2]) >= 0.99999) continue;
-       x[3]=f3(x1,y1,x2,y2,z1,z2);
-       if (TMath::Abs(x[3]) > 1.2) continue;
-       Double_t a=asin(x[2]);
-       Double_t zv=z1 - x[3]/x[4]*(a+asin(x[4]*x1-x[2]));
-       if (TMath::Abs(zv)>10.) continue; 
+       x[0]=y1;
+       x[1]=z1;
+       x[4]=f1(x1,y1,x2,y2,x3,y3);
+       if (TMath::Abs(x[4]) >= 0.0066) continue;
+       x[2]=f2(x1,y1,x2,y2,x3,y3);
+       //if (TMath::Abs(x[4]*x1-x[2]) >= 0.99999) continue;
+       x[3]=f3(x1,y1,x2,y2,z1,z2);
+       if (TMath::Abs(x[3]) > 1.2) continue;
+       Double_t a=asin(x[2]);
+       Double_t zv=z1 - x[3]/x[4]*(a+asin(x[4]*x1-x[2]));
+       if (TMath::Abs(zv-z3)>10.) continue; 
 
         Double_t sy1=kr1[is]->GetSigmaY2(), sz1=kr1[is]->GetSigmaZ2();
         Double_t sy2=kcl->GetSigmaY2(),     sz2=kcl->GetSigmaZ2();
-       Double_t sy3=100*0.025, sy=0.1, sz=0.1;
-
-       Double_t f40=(f1(x1,y1+sy,x2,y2,x3,y3)-x[4])/sy;
-       Double_t f42=(f1(x1,y1,x2,y2+sy,x3,y3)-x[4])/sy;
-       Double_t f43=(f1(x1,y1,x2,y2,x3,y3+sy)-x[4])/sy;
-       Double_t f20=(f2(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy;
-       Double_t f22=(f2(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy;
-       Double_t f23=(f2(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
-       Double_t f30=(f3(x1,y1+sy,x2,y2,z1,z2)-x[3])/sy;
-       Double_t f31=(f3(x1,y1,x2,y2,z1+sz,z2)-x[3])/sz;
-       Double_t f32=(f3(x1,y1,x2,y2+sy,z1,z2)-x[3])/sy;
-       Double_t f34=(f3(x1,y1,x2,y2,z1,z2+sz)-x[3])/sz;
+       //Double_t sy3=400*3./12., sy=0.1, sz=0.1;
+        Double_t sy3=25000*x[4]*x[4]+0.1, sy=0.1, sz=0.1;
+
+       Double_t f40=(f1(x1,y1+sy,x2,y2,x3,y3)-x[4])/sy;
+       Double_t f42=(f1(x1,y1,x2,y2+sy,x3,y3)-x[4])/sy;
+       Double_t f43=(f1(x1,y1,x2,y2,x3,y3+sy)-x[4])/sy;
+       Double_t f20=(f2(x1,y1+sy,x2,y2,x3,y3)-x[2])/sy;
+       Double_t f22=(f2(x1,y1,x2,y2+sy,x3,y3)-x[2])/sy;
+       Double_t f23=(f2(x1,y1,x2,y2,x3,y3+sy)-x[2])/sy;
+       Double_t f30=(f3(x1,y1+sy,x2,y2,z1,z2)-x[3])/sy;
+       Double_t f31=(f3(x1,y1,x2,y2,z1+sz,z2)-x[3])/sz;
+       Double_t f32=(f3(x1,y1,x2,y2+sy,z1,z2)-x[3])/sy;
+       Double_t f34=(f3(x1,y1,x2,y2,z1,z2+sz)-x[3])/sz;
 
         c[0]=sy1;
         c[1]=0.;       c[2]=sz1;
@@ -546,13 +547,17 @@ void AliTPCtracker::MakeSeeds(Int_t i1, Int_t i2) {
         c[14]=f40*sy1*f40+f42*sy2*f42+f43*sy3*f43;
 
         UInt_t index=kr1.GetIndex(is);
-       AliTPCseed *track=new AliTPCseed(index, x, c, x1, ns*alpha+shift);
+       AliTPCseed *track=new AliTPCseed(index, x, c, x1, ns*alpha+shift);
         Float_t l=fOuterSec->GetPadPitchWidth();
         track->SetSampledEdx(kr1[is]->GetQ()/l,0);
 
         Int_t rc=FollowProlongation(*track, i2);
         if (rc==0 || track->GetNumberOfClusters()<(i1-i2)/2) delete track;
-        else fSeeds->AddLast(track); 
+        else 
+         {
+           fSeeds->AddLast(track); 
+           cout<<"Adding seed   "<<fSeeds->GetEntries()<<"\r";
+         }
       }
     }
   }
@@ -590,33 +595,64 @@ Int_t AliTPCtracker::ReadSeeds(const TFile *inp) {
   }
   
   delete seed;
+
+  delete seedTree; //Thanks to Mariana Bondila
+
   savedir->cd();
 
   return 0;
 }
 
 //_____________________________________________________________________________
-Int_t AliTPCtracker::Clusters2Tracks(const TFile *inp, TFile *out) {
+Int_t AliTPCtracker::Clusters2Tracks() 
+{
   //-----------------------------------------------------------------
   // This is a track finder.
   //-----------------------------------------------------------------
-  TDirectory *savedir=gDirectory; 
+  Int_t retval = 0;
 
-  if (inp) {
-     TFile *in=(TFile*)inp;
-     if (!in->IsOpen()) {
-        cerr<<"AliTPCtracker::Clusters2Tracks(): input file is not open !\n";
-        return 1;
+  AliRunLoader* rl = AliRunLoader::GetRunLoader(fEvFolderName);
+  if (rl == 0x0)
+   {
+     Error("Clusters2Tracks","Can not get RL from specified folder %s",fEvFolderName.Data());
+     return 1;
+   }
+  
+  retval = rl->GetEvent(fEventN);
+  if (retval)
+   {
+     Error("Clusters2Tracks","Error while getting event %d",fEventN);
+     return 1;
+   }
+  
+  AliLoader* tpcl = rl->GetLoader("TPCLoader");
+  if (tpcl == 0x0)
+   {
+     Error("Clusters2Tracks","Can not get TPC Laoder from Run Loader");
+     return 1;
+   }
+  if ( tpcl->TreeR() == 0x0)
+   { 
+    retval = tpcl->LoadRecPoints("READ");
+    if (retval)
+     {
+       Error("Clusters2Tracks","Error while loading Reconstructed Points");
+       return 1;
      }
-  }
+   }
+   
+  if (tpcl->TreeT() == 0x0) tpcl->MakeTree("T");
 
-  if (!out->IsOpen()) {
-     cerr<<"AliTPCtracker::Clusters2Tracks(): output file is not open !\n";
-     return 2;
-  }
-
-  out->cd();
-  TTree tracktree("TPCf","Tree with TPC tracks");
+  TTree &tracktree = *(tpcl->TreeT());
+  
+  TBranch* br= tracktree.GetBranch("tracks");
+  if (br)
+   {
+     Error("Clusters2Tracks","Branch \"tracks\" already exists in TreeT for TPC");
+     return 1;
+   }
+  
   AliTPCtrack *iotrack=0;
   tracktree.Branch("tracks","AliTPCtrack",&iotrack,32000,0);
 
@@ -634,6 +670,7 @@ Int_t AliTPCtracker::Clusters2Tracks(const TFile *inp, TFile *out) {
 
   //tracking in outer sectors
   Int_t nseed=fSeeds->GetEntriesFast();
+  cout<<"nseed="<<nseed<<endl;
   Int_t i;
   for (i=0; i<nseed; i++) {
     AliTPCseed *pt=(AliTPCseed*)fSeeds->UncheckedAt(i), &t=*pt;
@@ -643,7 +680,7 @@ Int_t AliTPCtracker::Clusters2Tracks(const TFile *inp, TFile *out) {
     }
     delete fSeeds->RemoveAt(i);
   }  
-  UnloadOuterSectors();
+  //UnloadOuterSectors();
 
   //tracking in inner sectors
   LoadInnerSectors();
@@ -664,6 +701,7 @@ Int_t AliTPCtracker::Clusters2Tracks(const TFile *inp, TFile *out) {
        if (FollowProlongation(t)) {
           if (t.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
              t.CookdEdx();
+             CookLabel(pt,0.1); //For comparison only
              iotrack=pt;
              tracktree.Fill();
              UseClusters(&t,nc);
@@ -674,58 +712,81 @@ Int_t AliTPCtracker::Clusters2Tracks(const TFile *inp, TFile *out) {
     delete fSeeds->RemoveAt(i); 
   }  
   UnloadInnerSectors();
+  UnloadOuterSectors();
 
-  char   tname[100];
-  if (fEventN==-1) {
-    sprintf(tname,"TreeT_TPC");
-  }
-  else {
-    sprintf(tname,"TreeT_TPC_%d",fEventN);
-  }
-
-
-  tracktree.Write(tname);
-
+  tpcl->WriteTracks("OVERWRITE");
+  
   cerr<<"Number of found tracks : "<<found<<endl;
 
-  savedir->cd();
-
   return 0;
 }
 
 //_____________________________________________________________________________
-Int_t AliTPCtracker::PropagateBack(const TFile *inp, TFile *out) {
+Int_t AliTPCtracker::PropagateBack() 
+ {
   //-----------------------------------------------------------------
   // This function propagates tracks back through the TPC.
   //-----------------------------------------------------------------
-  fSeeds=new TObjArray(15000);
+  
+  cout<<"This method is not converted to NewIO yet\n";
+  return 1;
 
-  TFile *in=(TFile*)inp;
-  TDirectory *savedir=gDirectory; 
+  fSeeds=new TObjArray(15000);
+  Int_t retval = 0;
 
-  if (!in->IsOpen()) {
-     cerr<<"AliTPCtracker::PropagateBack(): ";
-     cerr<<"file with back propagated ITS tracks is not open !\n";
+  AliRunLoader* rl = AliRunLoader::GetRunLoader(fEvFolderName);
+  if (rl == 0x0)
+   {
+     Error("Clusters2Tracks","Can not get RL from specified folder %s",fEvFolderName.Data());
      return 1;
-  }
+   }
+  
+  retval = rl->GetEvent(fEventN);
+  if (retval)
+   {
+     Error("Clusters2Tracks","Error while getting event %d",fEventN);
+     return 1;
+   }
+  
+  AliLoader* tpcl = rl->GetLoader("TPCLoader");
+  if (tpcl == 0x0)
+   {
+     Error("Clusters2Tracks","Can not get TPC Laoder from Run Loader");
+     return 1;
+   }
 
-  if (!out->IsOpen()) {
-     cerr<<"AliTPCtracker::PropagateBack(): ";
-     cerr<<"file for back propagated TPC tracks is not open !\n";
-     return 2;
-  }
+  AliLoader* itsl = rl->GetLoader("ITSLoader");
+  if (tpcl == 0x0)
+   {
+     Error("Clusters2Tracks","Can not get ITS Laoder from Run Loader");
+     return 1;
+   }
+  
+  if (itsl->TreeT() == 0x0) itsl->LoadTracks();
 
-  in->cd();
-  TTree *bckTree=(TTree*)in->Get("ITSb");
+  TTree *bckTree=itsl->TreeT();
   if (!bckTree) {
      cerr<<"AliTPCtracker::PropagateBack() ";
      cerr<<"can't get a tree with back propagated ITS tracks !\n";
      return 3;
   }
+  
   AliTPCtrack *bckTrack=new AliTPCtrack; 
   bckTree->SetBranchAddress("tracks",&bckTrack);
 
-  TTree *tpcTree=(TTree*)in->Get("TPCf");
+
+  TFile* in = 0x0;
+  TFile* out = 0x0;
+  cout<<"And NOW there will be a segmentation violation!!!!\n";
+  bckTree=(TTree*)in->Get("TreeT_ITSb_0");
+  if (!bckTree) {
+     cerr<<"AliTPCtracker::PropagateBack() ";
+     cerr<<"can't get a tree with back propagated ITS tracks !\n";
+     return 3;
+  }
+
+
+  TTree *tpcTree=(TTree*)in->Get("TreeT_TPC_0");
   if (!tpcTree) {
      cerr<<"AliTPCtracker::PropagateBack() ";
      cerr<<"can't get a tree with TPC tracks !\n";
@@ -759,7 +820,7 @@ Int_t AliTPCtracker::PropagateBack(const TFile *inp, TFile *out) {
   }
 
   out->cd();
-  TTree backTree("TPCb","Tree with back propagated TPC tracks");
+  TTree backTree("TreeT_TPCb_0","Tree with back propagated TPC tracks");
   AliTPCtrack *otrack=0;
   backTree.Branch("tracks","AliTPCtrack",&otrack,32000,0);
 
@@ -774,7 +835,7 @@ Int_t AliTPCtracker::PropagateBack(const TFile *inp, TFile *out) {
      s.SetLabel(nc-1); //set number of the cluster to start with
 
      if (FollowBackProlongation(s,t)) {
-       UseClusters(&s);
+       UseClusters(&s);
         continue;
      }
      delete fSeeds->RemoveAt(i);
@@ -800,23 +861,23 @@ Int_t AliTPCtracker::PropagateBack(const TFile *inp, TFile *out) {
 
      if (s.Rotate(alpha)) {
         if (FollowBackProlongation(s,t)) {
-          if (s.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
-             s.CookdEdx();
+          if (s.GetNumberOfClusters() >= Int_t(0.4*nrows)) {
+             s.CookdEdx();
               s.SetLabel(t.GetLabel());
-             UseClusters(&s,nc);
+              UseClusters(&s,nc);
               otrack=ps;
               backTree.Fill();
-             cerr<<found++<<'\r';
+              cerr<<found++<<'\r';
               continue;
-          }
-       }
+          }
+       }
      }
      delete fSeeds->RemoveAt(i);
   }  
   UnloadOuterSectors();
 
   backTree.Write();
-  savedir->cd();
+
   cerr<<"Number of seeds: "<<nseed<<endl;
   cerr<<"Number of back propagated ITS tracks: "<<bckN<<endl;
   cerr<<"Number of back propagated TPC tracks: "<<found<<endl;
@@ -824,6 +885,9 @@ Int_t AliTPCtracker::PropagateBack(const TFile *inp, TFile *out) {
   delete bckTrack;
   delete tpcTrack;
 
+  delete bckTree; //Thanks to Mariana Bondila
+  delete tpcTree; //Thanks to Mariana Bondila
+
   return 0;
 }
 
@@ -906,20 +970,25 @@ void AliTPCtracker::AliTPCSector::Setup(const AliTPCParam *par, Int_t f) {
      fAlpha=par->GetInnerAngle();
      fAlphaShift=par->GetInnerAngleShift();
      fPadPitchWidth=par->GetInnerPadPitchWidth();
-     fPadPitchLength=par->GetInnerPadPitchLength();
+     f1PadPitchLength=par->GetInnerPadPitchLength();
+     f2PadPitchLength=f1PadPitchLength;
 
      fN=par->GetNRowLow();
+//     cout<<"par->GetNRowLow() = "<<par->GetNRowLow()<<"  fN = "<<fN<<endl;
      fRow=new AliTPCRow[fN];
      for (Int_t i=0; i<fN; i++) fRow[i].SetX(par->GetPadRowRadiiLow(i));
   } else {
      fAlpha=par->GetOuterAngle();
      fAlphaShift=par->GetOuterAngleShift();
      fPadPitchWidth=par->GetOuterPadPitchWidth();
-     fPadPitchLength=par->GetOuterPadPitchLength();
+     f1PadPitchLength=par->GetOuter1PadPitchLength();
+     f2PadPitchLength=par->GetOuter2PadPitchLength();
 
      fN=par->GetNRowUp();
      fRow=new AliTPCRow[fN];
-     for (Int_t i=0; i<fN; i++) fRow[i].SetX(par->GetPadRowRadiiUp(i));
+     for (Int_t i=0; i<fN; i++){ 
+     fRow[i].SetX(par->GetPadRowRadiiUp(i));
+     }
   } 
 }
 
@@ -944,6 +1013,7 @@ Int_t AliTPCtracker::AliTPCRow::Find(Double_t y) const {
   //-----------------------------------------------------------------------
   // Return the index of the nearest cluster 
   //-----------------------------------------------------------------------
+  if(fN<=0) return 0;
   if (y <= fClusters[0]->GetY()) return 0;
   if (y > fClusters[fN-1]->GetY()) return fN;
   Int_t b=0, e=fN-1, m=(b+e)/2;
@@ -961,23 +1031,34 @@ void AliTPCtracker::AliTPCseed::CookdEdx(Double_t low, Double_t up) {
   //-----------------------------------------------------------------
   Int_t i;
   Int_t nc=GetNumberOfClusters();
-
-  Int_t swap;//stupid sorting
-  do {
-    swap=0;
-    for (i=0; i<nc-1; i++) {
-      if (fdEdxSample[i]<=fdEdxSample[i+1]) continue;
-      Float_t tmp=fdEdxSample[i];
-      fdEdxSample[i]=fdEdxSample[i+1]; fdEdxSample[i+1]=tmp;
-      swap++;
-    }
-  } while (swap);
+  Int_t * index = new Int_t[nc];
+  TMath::Sort(nc, fdEdxSample,index,kFALSE);
 
   Int_t nl=Int_t(low*nc), nu=Int_t(up*nc);
   Float_t dedx=0;
-  for (i=nl; i<=nu; i++) dedx += fdEdxSample[i];
+  for (i=nl; i<=nu; i++) dedx += fdEdxSample[index[i]];
   dedx /= (nu-nl+1);
   SetdEdx(dedx);
+
+  delete [] index;
+
+  //Very rough PID
+  Double_t p=TMath::Sqrt((1.+ GetTgl()*GetTgl())/(Get1Pt()*Get1Pt()));
+
+  Double_t log1=TMath::Log(p+0.45), log2=TMath::Log(p+0.12);
+  if (p<0.6) {
+    if (dedx < 34 + 30/(p+0.45)/(p+0.45) + 24*log1) {SetMass(0.13957); return;}
+    if (dedx < 34 + 30/(p+0.12)/(p+0.12) + 24*log2) {SetMass(0.49368); return;}
+    SetMass(0.93827); return;
+  }
+
+  if (p<1.2) {
+    if (dedx < 34 + 30/(p+0.12)/(p+0.12) + 24*log2) {SetMass(0.13957); return;}
+    SetMass(0.93827); return;
+  }
+
+  SetMass(0.13957); return;
+
 }