A skeleton of the "cooked covariance" standalone tracker

diff --git a/ITS/UPGRADE/AliITSUTrackCooked.cxx b/ITS/UPGRADE/AliITSUTrackCooked.cxx
new file mode 100644 (file)
index 0000000..d807940
--- /dev/null
+++ b/ITS/UPGRADE/AliITSUTrackCooked.cxx
@@ -0,0 +1,117 @@
+//-------------------------------------------------------------------------
+//               Implementation of the ITSU track class
+//             based on the "cooked covariance" approach
+//-------------------------------------------------------------------------
+#include "AliITSUTrackCooked.h"
+#include "AliCluster.h"
+#include "AliESDtrack.h"
+
+ClassImp(AliITSUTrackCooked)
+
+AliITSUTrackCooked::AliITSUTrackCooked(): 
+AliKalmanTrack()
+{
+  //--------------------------------------------------------------------
+  // This default constructor needs to be provided
+  //--------------------------------------------------------------------
+    for (Int_t i=0; i<2*AliITSUTrackerCooked::kNLayers; i++) {
+        fIndex[i]=0;
+    }
+}
+
+AliITSUTrackCooked::AliITSUTrackCooked(const AliESDtrack &t):
+AliKalmanTrack()
+{
+    //--------------------------------------------------------------------
+    // Constructior from an ESD track
+    //--------------------------------------------------------------------
+    Set(t.GetX(), t.GetAlpha(), t.GetParameter(), t.GetCovariance());
+    SetLabel(t.GetITSLabel());
+    SetChi2(t.GetITSchi2());
+    SetNumberOfClusters(t.GetITSclusters(fIndex));
+}
+
+AliITSUTrackCooked::~AliITSUTrackCooked()
+{
+  //--------------------------------------------------------------------
+  // Virtual destructor
+  //--------------------------------------------------------------------
+}
+
+void AliITSUTrackCooked::SetClusterIndex(Int_t l, Int_t i)
+{
+    //--------------------------------------------------------------------
+    // Set the cluster index
+    //--------------------------------------------------------------------
+    Int_t idx = (l<<28) + i;
+    Int_t n=GetNumberOfClusters();
+    fIndex[n]=idx;
+    SetNumberOfClusters(n+1);
+}
+
+Double_t AliITSUTrackCooked::GetPredictedChi2(const AliCluster *c) const {
+  //-----------------------------------------------------------------
+  // This function calculates a predicted chi2 increment.
+  //-----------------------------------------------------------------
+  Double_t p[2]={c->GetY(), c->GetZ()};
+  Double_t cov[3]={c->GetSigmaY2(), 0., c->GetSigmaZ2()};
+  return AliExternalTrackParam::GetPredictedChi2(p,cov);
+}
+
+Bool_t AliITSUTrackCooked::PropagateTo(Double_t xk, Double_t t,Double_t x0rho) {
+  //------------------------------------------------------------------
+  // This function propagates a track
+  // t is the material thicknes in units X/X0
+  // x0rho is the material X0*density
+  //------------------------------------------------------------------
+  Double_t xOverX0,xTimesRho; 
+  xOverX0 = t; xTimesRho = t*x0rho;
+  if (!CorrectForMeanMaterial(xOverX0,xTimesRho,GetMass(),kTRUE)) return kFALSE;
+
+  Double_t bz=GetBz();
+  if (!AliExternalTrackParam::PropagateTo(xk,bz)) return kFALSE;
+  //Double_t b[3]; GetBxByBz(b);
+  //if (!AliExternalTrackParam::PropagateToBxByBz(xk,b)) return kFALSE;
+
+  return kTRUE;
+}
+
+Bool_t AliITSUTrackCooked::
+Propagate(const AliCluster *c, Double_t t, Double_t x0rho) {
+  //------------------------------------------------------------------
+  // This function propagates a track to the plane the cluster belong to
+  // t is the material thicknes in units X/X0
+  // x0rho is the material X0*density
+  //------------------------------------------------------------------
+  Double_t xOverX0,xTimesRho; 
+  xOverX0 = t; xTimesRho = t*x0rho;
+  if (!CorrectForMeanMaterial(xOverX0,xTimesRho,GetMass(),kTRUE)) return kFALSE;
+
+  Float_t xk, alpha;
+  if (!c->GetXAlphaRefPlane(xk,alpha)) return kFALSE;
+
+  Double_t bz=GetBz();
+  if (!AliExternalTrackParam::Propagate(alpha, xk, bz)) return kFALSE;
+  //Double_t b[3]; GetBxByBz(b);
+  //if (!AliExternalTrackParam::PropagateBxByBz(alpha, xk, b)) return kFALSE;
+
+  return kTRUE;
+}
+
+Bool_t AliITSUTrackCooked::Update(const AliCluster *c, Double_t chi2, Int_t idx)
+{
+  //--------------------------------------------------------------------
+  // Update track params
+  //--------------------------------------------------------------------
+  Double_t p[2]={c->GetY(), c->GetZ()};
+  Double_t cov[3]={c->GetSigmaY2(), c->GetSigmaYZ(), c->GetSigmaZ2()};
+
+  if (!AliExternalTrackParam::Update(p,cov)) return kFALSE;
+
+  Int_t n=GetNumberOfClusters();
+  fIndex[n]=idx;
+  SetNumberOfClusters(n+1);
+  SetChi2(GetChi2()+chi2);
+
+  return kTRUE;
+}

diff --git a/ITS/UPGRADE/AliITSUTrackCooked.h b/ITS/UPGRADE/AliITSUTrackCooked.h
new file mode 100644 (file)
index 0000000..e067e4a
--- /dev/null
+++ b/ITS/UPGRADE/AliITSUTrackCooked.h
@@ -0,0 +1,40 @@
+#ifndef ALIITSUTRACKCooked_H
+#define ALIITSUTRACKCooked_H
+
+//-------------------------------------------------------------------------
+//                   The stand-alone ITSU track
+//           Created with a "cooked covariance" approach
+//-------------------------------------------------------------------------
+
+#include "AliKalmanTrack.h"
+#include "AliITSUTrackerCooked.h"
+
+class AliCluster;
+class AliESDtrack;
+
+//-------------------------------------------------------------------------
+class AliITSUTrackCooked : public AliKalmanTrack {
+public:
+    AliITSUTrackCooked();
+    AliITSUTrackCooked(const AliESDtrack &t);
+    virtual ~AliITSUTrackCooked();
+    
+//These functions must be provided
+    Double_t GetPredictedChi2(const AliCluster *c) const;
+    Bool_t PropagateTo(Double_t xr, Double_t thick, Double_t rho);
+    Bool_t Update(const AliCluster *c, Double_t chi2, Int_t idx);
+
+//Other functions    
+    Int_t GetClusterIndex(Int_t i) const { return fIndex[i]; }
+    Double_t GetPIDsignal() const { return 0; }
+
+    void SetClusterIndex(Int_t layer, Int_t index);
+    Bool_t Propagate(const AliCluster *c, Double_t thick, Double_t x0rho);
+    
+private:
+    Int_t fIndex[2*AliITSUTrackerCooked::kNLayers]; // indices of associated clusters
+  
+    ClassDef(AliITSUTrackCooked,1)   //ITSU stand-alone track
+};
+
+#endif

diff --git a/ITS/UPGRADE/AliITSUTrackerCooked.cxx b/ITS/UPGRADE/AliITSUTrackerCooked.cxx
new file mode 100644 (file)
index 0000000..8109801
--- /dev/null
+++ b/ITS/UPGRADE/AliITSUTrackerCooked.cxx
@@ -0,0 +1,429 @@
+//-------------------------------------------------------------------------
+//               Implementation of the ITS tracker class
+//    The pattern recongintion based on the "cooked covariance" approach
+//-------------------------------------------------------------------------
+
+#include <Riostream.h>
+#include <TTree.h>
+#include <TClonesArray.h>
+
+#include "AliLog.h"
+#include "AliESDEvent.h"
+#include "AliITSUClusterPix.h"
+#include "AliITSUGeomTGeo.h"
+#include "AliITSUTrackerCooked.h"
+#include "AliITSUTrackCooked.h" 
+
+ClassImp(AliITSUTrackerCooked)
+
+//************************************************
+// Constants hardcoded for the moment:
+//************************************************
+// radial positions of layers: default contructor
+const 
+Double_t klRadius[7]={2.34, 3.15, 3.93, 19.61, 24.55, 34.39, 39.34}; //tdr6
+// seed "windows" in z and phi: MakeSeeds
+const Double_t kzWin=0.33, kpWin=3.14/4;
+// Maximal accepted impact parameters for the seeds 
+const Double_t kmaxDCAxy=7.;
+const Double_t kmaxDCAz= 2.;
+
+//************************************************
+// TODO:
+//************************************************
+// Seeding:
+// Precalculate cylidnrical (r,phi) for the clusters;
+// use exact r's for the clusters
+
+
+AliITSUTrackerCooked::AliITSUlayer
+              AliITSUTrackerCooked::fgLayers[AliITSUTrackerCooked::kNLayers];
+
+AliITSUTrackerCooked::AliITSUlayer::AliITSUlayer():
+  fR(0),
+  fN(0),
+  fNsel(0) 
+{
+  //--------------------------------------------------------------------
+  // This default constructor needs to be provided
+  //--------------------------------------------------------------------
+  for (Int_t i=0; i<kMaxClusterPerLayer; i++) {
+    fClusters[i]=0;
+    fIndex[i]=-1;
+  }
+}
+
+void AliITSUTrackerCooked::AliITSUlayer::InsertClusters(TClonesArray *clusters)
+{
+  //--------------------------------------------------------------------
+  // Load clusters to this layer
+  //--------------------------------------------------------------------
+  Int_t ncl=clusters->GetEntriesFast();
+ 
+  while (ncl--) {
+     AliITSUClusterPix *c=(AliITSUClusterPix*)clusters->UncheckedAt(ncl);
+     c->GoToFrameGlo();
+     //if (!c->Misalign()) AliWarning("Can't misalign this cluster !");
+     InsertCluster(new AliITSUClusterPix(*c));
+  }
+
+}
+
+void AliITSUTrackerCooked::AliITSUlayer::DeleteClusters()
+{
+  //--------------------------------------------------------------------
+  // Load clusters to this layer
+  //--------------------------------------------------------------------
+  for (Int_t i=0; i<fN; i++) delete fClusters[i];
+  fN=0;
+}
+
+Int_t 
+AliITSUTrackerCooked::AliITSUlayer::InsertCluster(AliCluster *c) {
+  //--------------------------------------------------------------------
+  // This function inserts a cluster to this layer in increasing
+  // order of the cluster's fZ
+  //--------------------------------------------------------------------
+  if (fN>=kMaxClusterPerLayer) {
+     ::Error("InsertCluster","Too many clusters !\n");
+     return 1;
+  }
+  if (fN==0) fClusters[0]=c;
+  else {
+     Int_t i=FindClusterIndex(c->GetZ());
+     Int_t k=fN-i;
+     memmove(fClusters+i+1 ,fClusters+i,k*sizeof(AliCluster*));
+     fClusters[i]=c;
+  }
+  fN++;
+  return 0;
+}
+
+Int_t 
+AliITSUTrackerCooked::AliITSUlayer::FindClusterIndex(Double_t z) const {
+  //--------------------------------------------------------------------
+  // This function returns the index of the first 
+  // with its fZ >= "z". 
+  //--------------------------------------------------------------------
+  if (fN==0) return 0;
+
+  Int_t b=0;
+  if (z <= fClusters[b]->GetZ()) return b;
+
+  Int_t e=b+fN-1;
+  if (z > fClusters[e]->GetZ()) return e+1;
+
+  Int_t m=(b+e)/2;
+  for (; b<e; m=(b+e)/2) {
+    if (z > fClusters[m]->GetZ()) b=m+1;
+    else e=m; 
+  }
+  return m;
+}
+
+
+AliITSUTrackerCooked::AliITSUTrackerCooked(): 
+  AliTracker(),fSeeds(0) 
+{
+  //--------------------------------------------------------------------
+  // This default constructor needs to be provided
+  //--------------------------------------------------------------------
+
+  AliITSUGeomTGeo *gm  = new AliITSUGeomTGeo(kTRUE,kTRUE);
+  AliITSUClusterPix::SetGeom(gm);
+
+  for (Int_t i=0; i<kNLayers; i++) fgLayers[i].SetR(klRadius[i]);
+
+  // Some default primary vertex
+  Double_t xyz[]={0.,0.,0.};
+  Double_t ers[]={2.,2.,2.};
+
+  SetVertex(xyz,ers);
+
+}
+
+AliITSUTrackerCooked::~AliITSUTrackerCooked() 
+{
+  //--------------------------------------------------------------------
+  // Virtual destructor
+  //--------------------------------------------------------------------
+
+  UnloadClusters();
+
+  if (fSeeds) {fSeeds->Delete(); delete fSeeds;} 
+  fSeeds=0;
+}
+
+static Double_t 
+f1(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3)
+{
+    //-----------------------------------------------------------------
+    // Initial approximation of the track curvature
+    //-----------------------------------------------------------------
+    Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
+    Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
+                    (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
+    Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
+                    (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
+    
+    Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
+    
+    return -xr*yr/sqrt(xr*xr+yr*yr);
+}
+
+static Double_t 
+f2(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3)
+{
+    //-----------------------------------------------------------------
+    // Initial approximation of the track curvature times center of curvature
+    //-----------------------------------------------------------------
+    Double_t d=(x2-x1)*(y3-y2)-(x3-x2)*(y2-y1);
+    Double_t a=0.5*((y3-y2)*(y2*y2-y1*y1+x2*x2-x1*x1)-
+                    (y2-y1)*(y3*y3-y2*y2+x3*x3-x2*x2));
+    Double_t b=0.5*((x2-x1)*(y3*y3-y2*y2+x3*x3-x2*x2)-
+                    (x3-x2)*(y2*y2-y1*y1+x2*x2-x1*x1));
+    
+    Double_t xr=TMath::Abs(d/(d*x1-a)), yr=d/(d*y1-b);
+    
+    return -a/(d*y1-b)*xr/sqrt(xr*xr+yr*yr);
+}
+
+static Double_t 
+f3(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t z1, Double_t z2)
+{
+    //-----------------------------------------------------------------
+    // Initial approximation of the tangent of the track dip angle
+    //-----------------------------------------------------------------
+    return (z1 - z2)/sqrt((x1-x2)*(x1-x2)+(y1-y2)*(y1-y2));
+}
+
+Bool_t AliITSUTrackerCooked::
+AddCookedSeed(const Float_t r1[3], Int_t l1, Int_t i1, 
+              const Float_t r2[3], Int_t l2, Int_t i2,
+              const AliCluster *c3,Int_t l3, Int_t i3) 
+{
+    //--------------------------------------------------------------------
+    // This is the main cooking function.
+    // Creates seed parameters out of provided clusters.
+    //--------------------------------------------------------------------
+    Float_t x,a;
+    if (!c3->GetXAlphaRefPlane(x,a)) return kFALSE;
+
+    Double_t par[5];
+    par[0]=c3->GetY();
+    par[1]=c3->GetZ();
+    Float_t r3[3]; c3->GetGlobalXYZ(r3);
+    Double_t crv=f1(r1[0],r1[1],r2[0],r2[1],r3[0],r3[1]); //curvature
+    Double_t cx0=f2(r1[0],r1[1],r2[0],r2[1],r3[0],r3[1]); //curvature*x0
+    Double_t tgl12=f3(r1[0],r1[1],r2[0],r2[1],r1[2],r2[2]);
+    Double_t tgl23=f3(r2[0],r2[1],r3[0],r3[1],r2[2],r3[2]);
+    par[2]=TMath::Sin(TMath::ASin(r3[0]*crv - cx0) - a); //FIXME
+    if (TMath::Abs(par[2]) >= kAlmost1) return kFALSE;
+    par[3]=0.5*(tgl12 + tgl23);
+    Double_t bz=GetBz();
+    par[4]=(TMath::Abs(bz) < kAlmost0Field) ? par[4]=kAlmost0 : crv/(bz*kB2C);
+
+    Double_t cov[15];
+    for (Int_t i=0; i<15; i++) cov[i]=0.;
+    cov[0] =0.0005*0.0005*10;
+    cov[2] =0.0005*0.0005*10;
+    cov[5] =0.007*0.007*10;   //FIXME all these lines
+    cov[9] =0.007*0.007*10;
+    cov[14]=0.1*0.1*10;
+
+    AliITSUTrackCooked *seed=new AliITSUTrackCooked();
+    seed->Set(Double_t(x), Double_t(a), par, cov);
+
+    Float_t dz[2]; 
+    seed->GetDZ(GetX(),GetY(),GetZ(),GetBz(),dz);
+    if (TMath::Abs(dz[0]) > kmaxDCAxy) {delete seed; return kFALSE;} 
+    if (TMath::Abs(dz[1]) > kmaxDCAz ) {delete seed; return kFALSE;} 
+
+    seed->SetClusterIndex(l1,i1);
+    seed->SetClusterIndex(l2,i2);
+    seed->SetClusterIndex(l3,i3);
+
+    fSeeds->AddLast(seed);
+
+    return kTRUE;
+}
+
+Int_t AliITSUTrackerCooked::MakeSeeds(Int_t l1, Int_t l2) {
+  //--------------------------------------------------------------------
+  // This is the main pattern recongition function.
+  // Creates seeds out of two clusters and another point.
+  //--------------------------------------------------------------------
+   Int_t l3=5;
+
+   if (fSeeds) {fSeeds->Delete(); delete fSeeds;}
+   fSeeds=new TObjArray(77777);
+
+   Double_t zv=GetZ();
+
+   AliITSUlayer &layer1=fgLayers[l1];
+   AliITSUlayer &layer2=fgLayers[l2];
+   AliITSUlayer &layer3=fgLayers[l3];
+   Double_t r1=layer1.GetR();
+   Double_t r2=layer2.GetR();
+   Double_t r3=layer3.GetR();
+
+   Int_t nClusters1=layer1.GetNumberOfClusters();
+   Int_t nClusters2=layer2.GetNumberOfClusters();
+   Int_t nClusters3=layer3.GetNumberOfClusters();
+   for (Int_t n1=0; n1<nClusters1; n1++) {
+     AliCluster *c1=layer1.GetCluster(n1);
+     //
+     //Int_t lab=c1->GetLabel(0);
+     //
+     Double_t z1=c1->GetZ();
+     Float_t xyz1[3]; c1->GetGlobalXYZ(xyz1);
+     Double_t phi1=TMath::ATan2(xyz1[1],xyz1[0]);
+     Double_t zr2=zv + r2/r1*(z1-zv);
+     Int_t start2=layer2.FindClusterIndex(zr2-kzWin);
+     for (Int_t n2=start2; n2<nClusters2; n2++) {
+         AliCluster *c2=layer2.GetCluster(n2);
+         //
+         //if (c2->GetLabel(0)!=lab) continue;
+        //
+         Double_t z2=c2->GetZ();
+         
+         if (z2 > (zr2+kzWin)) break;  //check in Z
+         Float_t xyz2[3]; c2->GetGlobalXYZ(xyz2);
+         Double_t phi2=TMath::ATan2(xyz2[1],xyz2[0]);
+         if (TMath::Abs(phi2-phi1) > kpWin) continue;  //check in Phi
+ 
+         Double_t zr3=z1 + (r3-r1)/(r2-r1)*(z2-z1);
+         Int_t start3=layer3.FindClusterIndex(zr3-kzWin/2);
+         for (Int_t n3=start3; n3<nClusters3; n3++) {
+             AliCluster *c3=layer3.GetCluster(n3);
+             //
+             //if (c3->GetLabel(0)!=lab) continue;
+             //
+             Double_t z3=c3->GetZ();
+         
+             if (z3 > (zr3+kzWin/2)) break;  //check in Z
+             Float_t xyz3[3]; c3->GetGlobalXYZ(xyz3);
+             Double_t phi3=TMath::ATan2(xyz3[1],xyz3[0]);
+             if (TMath::Abs(phi2-phi3) > kpWin/2) continue;  //check in Phi
+
+             AliITSUClusterPix cc(*((AliITSUClusterPix*)c2));
+             cc.GoToFrameTrk();
+             AddCookedSeed(xyz1,l1,n1, xyz3,l3,n3, &cc,l2,n2);
+
+        }
+     }
+   }
+
+   return fSeeds->GetEntriesFast();
+}
+
+Int_t AliITSUTrackerCooked::Clusters2Tracks(AliESDEvent *event) {
+  //--------------------------------------------------------------------
+  // This is the main tracking function
+  // The clusters must already be loaded
+  //--------------------------------------------------------------------
+  Int_t l1=6, l2=4;
+
+  // Possibly, create the track "seeds" (combinatorial)
+  Int_t nSeeds=MakeSeeds(l1,l2);
+
+  // Possibly, icrement the seeds with additional clusters (Kalman)
+
+  // Possibly, (re)fit the found tracks 
+
+
+  Int_t nClusters1=fgLayers[l1].GetNumberOfClusters();
+  Int_t nClusters2=fgLayers[l2].GetNumberOfClusters();
+  cout<<nClusters1<<' '<<nClusters2<<endl;
+
+    Int_t ngood=0;
+    for (Int_t s=0; s<nSeeds; s++) {
+        AliITSUTrackCooked *track=(AliITSUTrackCooked*)fSeeds->At(s);
+        CookLabel(track,0.);
+        if (track->GetLabel()>0) {ngood++;
+            //cout<<track->GetLabel()<<endl;
+        }
+        AliESDtrack iotrack;
+        iotrack.UpdateTrackParams(track,AliESDtrack::kITSin);
+        iotrack.SetLabel(track->GetLabel());
+        event->AddTrack(&iotrack);
+    }
+    cout<<"Good tracks "<<ngood<<endl;
+    cout<<"Good tracks/seeds "<<Float_t(ngood)/nSeeds<<endl;
+
+  if (fSeeds) {fSeeds->Delete(); delete fSeeds;}
+  fSeeds=0;
+    
+  return 0;
+}
+
+Int_t AliITSUTrackerCooked::PropagateBack(AliESDEvent *event) {
+  //--------------------------------------------------------------------
+  // Here, we implement the Kalman smoother ?
+  // The clusters must already be loaded
+  //--------------------------------------------------------------------
+    Int_t n=event->GetNumberOfTracks();
+    for (Int_t i=0; i<n; i++) {
+        AliESDtrack *esdTrack=event->GetTrack(i);
+        AliITSUTrackCooked track(*esdTrack);
+        esdTrack->UpdateTrackParams(&track,AliESDtrack::kITSout);
+    }
+    
+  return 0;
+}
+
+Int_t AliITSUTrackerCooked::RefitInward(AliESDEvent *event) {
+  //--------------------------------------------------------------------
+  // Some final refit, after the outliers get removed by the smoother ?  
+  // The clusters must be loaded
+  //--------------------------------------------------------------------
+    Int_t n=event->GetNumberOfTracks();
+    for (Int_t i=0; i<n; i++) {
+        AliESDtrack *esdTrack=event->GetTrack(i);
+        AliITSUTrackCooked track(*esdTrack);
+        esdTrack->UpdateTrackParams(&track,AliESDtrack::kITSrefit);
+    }
+    
+  return 0;
+}
+
+Int_t AliITSUTrackerCooked::LoadClusters(TTree *cTree) {
+  //--------------------------------------------------------------------
+  // This function reads the ITSU clusters from the tree,
+  // sort them, distribute over the internal tracker arrays, etc
+  //--------------------------------------------------------------------
+  if (!cTree) {
+     AliFatal("No cluster tree !");
+     return 1;
+  }
+
+  //This TClonesArray is not the owner of the clusters
+  TClonesArray dummy("AliITSUClusterPix",kMaxClusterPerLayer), *clusters=&dummy;
+
+  for (Int_t i=0; i<kNLayers; i++) {
+      TBranch *br = cTree->GetBranch(Form("ITSRecPoints%d",i));
+      if (!br) AliFatal(Form("No cluster branch for layer %d",i));
+      br->SetAddress(&clusters);
+      br->GetEvent(0);
+      fgLayers[i].InsertClusters(clusters);
+      clusters->Delete();
+  }
+
+  return 0;
+}
+
+void AliITSUTrackerCooked::UnloadClusters() {
+  //--------------------------------------------------------------------
+  // This function unloads ITSU clusters from the RAM
+  //--------------------------------------------------------------------
+  for (Int_t i=0; i<kNLayers; i++) fgLayers[i].DeleteClusters();
+}
+
+AliCluster *AliITSUTrackerCooked::GetCluster(Int_t index) const {
+  //--------------------------------------------------------------------
+  //       Return pointer to a given cluster
+  //--------------------------------------------------------------------
+    Int_t l=(index & 0xf0000000) >> 28;
+    Int_t c=(index & 0x0fffffff) >> 00;
+    return fgLayers[l].GetCluster(c);
+}

diff --git a/ITS/UPGRADE/AliITSUTrackerCooked.h b/ITS/UPGRADE/AliITSUTrackerCooked.h
new file mode 100644 (file)
index 0000000..1a64687
--- /dev/null
+++ b/ITS/UPGRADE/AliITSUTrackerCooked.h
@@ -0,0 +1,89 @@
+#ifndef ALIITSUTRACKERCooked_H
+#define ALIITSUTRACKERCooked_H
+
+//-------------------------------------------------------------------------
+//                   The stand-alone ITSU tracker
+//    The pattern recongintion based on the "cooked covariance" approach
+//-------------------------------------------------------------------------
+
+#include "AliTracker.h"
+
+class TTree;
+class TClonesArray;
+class TObjArray;
+
+class AliESDEvent;
+class AliCluster;
+class AliITSUClusterPix;
+
+//-------------------------------------------------------------------------
+class AliITSUTrackerCooked : public AliTracker {
+public:
+  enum{kNLayers=7, kMaxClusterPerLayer=70000};
+  AliITSUTrackerCooked();
+  virtual ~AliITSUTrackerCooked();
+
+  // These functions must be implemented 
+  Int_t Clusters2Tracks(AliESDEvent *event);
+  Int_t PropagateBack(AliESDEvent *event);
+  Int_t RefitInward(AliESDEvent *event);
+  Int_t LoadClusters(TTree *ct);
+  void UnloadClusters();
+  AliCluster *GetCluster(Int_t index) const;
+
+  // internal helper classes
+  class AliITSUlayer;
+
+protected:
+  AliITSUTrackerCooked(const AliITSUTrackerCooked&);
+  // Other protected functions
+  Int_t MakeSeeds(Int_t layer1, Int_t layer2);
+  Bool_t AddCookedSeed(const Float_t r1[3], Int_t l1, Int_t i1,
+                       const Float_t r2[3], Int_t l2, Int_t i2,
+                       const AliCluster *c3,Int_t l3, Int_t i3);
+
+private:
+  AliITSUTrackerCooked &operator=(const AliITSUTrackerCooked &tr);
+
+  // Data members
+  // Internal tracker arrays, layers, modules, etc
+  static AliITSUlayer fgLayers[kNLayers];// Layers
+    
+  TObjArray *fSeeds; // Track seeds
+  
+  ClassDef(AliITSUTrackerCooked,1)   //ITSU stand-alone tracker
+};
+
+
+
+// The helper classes
+class AliITSUTrackerCooked::AliITSUlayer {
+  public:
+    AliITSUlayer();
+   ~AliITSUlayer(){DeleteClusters();}
+
+    void InsertClusters(TClonesArray *clusters);
+    void SetR(Double_t r) {fR=r;}
+    void DeleteClusters();
+
+    Int_t FindClusterIndex(Double_t z) const;
+    Double_t GetR() const {return fR;}
+    AliCluster *GetCluster(Int_t i) const { return fClusters[i]; } 
+    Int_t GetNumberOfClusters() const {return fN;}
+
+  protected:
+    AliITSUlayer(const AliITSUlayer&);
+    AliITSUlayer &operator=(const AliITSUlayer &tr);  
+    Int_t InsertCluster(AliCluster *c);
+
+    Double_t fR;                // mean radius of this layer
+
+    AliCluster *fClusters[kMaxClusterPerLayer]; // array of clusters
+    Int_t fN; //number of clusters 
+
+    Int_t fIndex[kMaxClusterPerLayer]; // Indexes of selected clusters 
+    Int_t fNsel;                       // number of selected clusters
+  };
+
+
+#endif

diff --git a/ITS/UPGRADE/CMakelibITSUpgradeRec.pkg b/ITS/UPGRADE/CMakelibITSUpgradeRec.pkg
old mode 100755 (executable)
new mode 100644 (file)
index 012bd7b..03bd904
--- a/ITS/UPGRADE/CMakelibITSUpgradeRec.pkg
+++ b/ITS/UPGRADE/CMakelibITSUpgradeRec.pkg
@@ -38,6 +38,8 @@ set ( SRCS
     AliITSUSeed.cxx
     AliITSUTrackerGlo.cxx
     AliITSUTrackerSA.cxx
+    AliITSUTrackerCooked.cxx
+    AliITSUTrackCooked.cxx
     AliITSUTrackCond.cxx
     AliITSUTrackHyp.cxx
     AliITSUMatLUT.cxx

diff --git a/ITS/UPGRADE/ITSUpgradeRecLinkDef.h b/ITS/UPGRADE/ITSUpgradeRecLinkDef.h
old mode 100755 (executable)
new mode 100644 (file)
index 6bb7fa3..4ce519c
--- a/ITS/UPGRADE/ITSUpgradeRecLinkDef.h
+++ b/ITS/UPGRADE/ITSUpgradeRecLinkDef.h
@@ -24,6 +24,8 @@
 #pragma link C++ class AliITSUSeed+;
 #pragma link C++ class AliITSUTrackerGlo+;
 #pragma link C++ class AliITSUTrackerSA+;
+#pragma link C++ class AliITSUTrackerCooked+;
+#pragma link C++ class AliITSUTrackCooked+;
 #pragma link C++ class AliITSUTrackCond+;
 #pragma link C++ class AliITSUTrackHyp+;
 #pragma link C++ class AliITSUMatLUT+;

diff --git a/ITS/UPGRADE/macros/testTrackerCooked.C b/ITS/UPGRADE/macros/testTrackerCooked.C
new file mode 100644 (file)
index 0000000..83db268
--- /dev/null
+++ b/ITS/UPGRADE/macros/testTrackerCooked.C
@@ -0,0 +1,91 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+   #include <Riostream.h>
+   #include <TSystem.h>
+   #include <TFile.h>
+   #include <TTree.h>
+   #include <TGeoManager.h>
+   #include <TGeoGlobalMagField.h>
+   #include <TStopwatch.h>
+
+   #include "AliRunLoader.h"
+   #include "AliHeader.h"
+   #include "AliGenEventHeader.h"
+   #include "AliMagF.h"
+   #include "AliESDEvent.h"
+   #include "AliITSUTrackerCooked.h"
+#endif
+
+extern TSystem *gSystem;
+
+const AliESDVertex *SetMCvertex(const AliRunLoader *rl, AliTracker *tr);
+
+void testTrackerCooked() {
+    gSystem->Load("libITSUpgradeBase");
+    gSystem->Load("libITSUpgradeSim");
+    gSystem->Load("libITSUpgradeRec");
+
+    TGeoGlobalMagField::Instance()->
+          SetField(new AliMagF("Maps","Maps", -1., -1., AliMagF::k5kG));
+    TGeoManager::Import("geometry.root");
+
+    TFile *esdFile=TFile::Open("AliESDs.root","recreate");
+    TTree *esdTree = new TTree("esdTree", "Tree with ESD objects");
+    AliESDEvent *esd=new AliESDEvent();
+    esd->CreateStdContent();
+    esd->WriteToTree(esdTree);
+    
+    TFile *clsFile=TFile::Open("ITS.RecPoints.root");
+    AliITSUTrackerCooked *tracker=new AliITSUTrackerCooked();
+    
+    AliRunLoader *rl = AliRunLoader::Open("galice.root","something");
+    rl->LoadHeader();
+
+    TStopwatch timer;    
+    Int_t nEvents=rl->GetNumberOfEvents();
+    for (Int_t i=0; i<nEvents; i++) {
+        cout<<"\nEvent number "<<i<<endl;
+        rl->GetEvent(i);
+
+        const AliESDVertex *vtx=SetMCvertex(rl,tracker);
+        esd->SetPrimaryVertexSPD(vtx);              
+
+        TTree *cTree=(TTree *)clsFile->Get(Form("Event%d/TreeR",i));
+        tracker->LoadClusters(cTree);
+        tracker->Clusters2Tracks(esd);
+        //tracker->PropagateBack(esd);
+        tracker->RefitInward(esd);
+        tracker->UnloadClusters();
+
+        Int_t n=esd->GetNumberOfTracks();
+        for (Int_t t=0; t<n; t++) {
+         AliESDtrack *track=esd->GetTrack(t);
+          if (!track->RelateToVertex(vtx, tracker->GetBz(), 33)) continue;
+          //Double_t r[3]; track->GetXYZ(r);
+         //cout<<r[0]<<' '<<r[1]<<' '<<r[2]-vtx->GetZ()<<endl;
+       }
+
+        esdTree->Fill();
+        esd->Reset();
+        delete vtx;
+    }
+    timer.Stop(); timer.Print();
+ 
+    delete tracker;
+    delete clsFile;
+    esdFile->cd();
+    esdTree->Write();
+    delete esd;
+    delete esdFile;
+}
+
+const AliESDVertex *SetMCvertex(const AliRunLoader *rl, AliTracker *tracker) {
+    AliGenEventHeader *h=rl->GetHeader()->GenEventHeader();
+    TArrayF vtx(3);
+    h->PrimaryVertex(vtx);
+    cout<<"Vertex "<<vtx[0]<<' '<<vtx[1]<<' '<<vtx[2]<<endl;
+    Double_t xyz[]={vtx[0],vtx[1],vtx[2]};
+    Double_t ers[]={2.,2.,2.};
+    tracker->SetVertex(xyz,ers);
+    AliESDVertex *vertex=new AliESDVertex(xyz,ers);
+    return vertex;
+}