]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TPC/AliTPCTracklet.cxx
Some of the coding violations corrected
[u/mrichter/AliRoot.git] / TPC / AliTPCTracklet.cxx
index 7b3a523b7a1036990e4a5a2d894bc1ded978180d..9c670333ea61aa2782f67b37423702da6b8668de 100755 (executable)
@@ -37,6 +37,8 @@ using namespace std;
 ClassImp(AliTPCTracklet)
 
 const Double_t AliTPCTracklet::kB2C=0.299792458e-3;
+Float_t  AliTPCTracklet::fgEdgeCutY=3;
+Float_t  AliTPCTracklet::fgEdgeCutX=0;
 
 AliTPCTracklet::AliTPCTracklet() 
   : fNClusters(0),fNStoredClusters(0),fClusters(0),fSector(-1),fOuter(0),
@@ -59,6 +61,7 @@ AliTPCTracklet::AliTPCTracklet(const AliTPCseed *track,Int_t sector,
   
   for (Int_t i=0;i<160;++i) {
     AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (c && RejectCluster(c)) continue;
     if (c&&c->GetDetector()==sector)
       ++fNClusters;
   }
@@ -67,8 +70,9 @@ AliTPCTracklet::AliTPCTracklet(const AliTPCseed *track,Int_t sector,
     fClusters=new AliTPCclusterMI[fNClusters];
     for (Int_t i=0;i<160;++i) {
       AliTPCclusterMI *c=track->GetClusterPointer(i);
+      if (c && RejectCluster(c)) continue;
       if (c&&c->GetDetector()==sector)
-       fClusters[fNStoredClusters]=c;
+       fClusters[fNStoredClusters]=*c;
       ++fNStoredClusters;
     }
   }
@@ -88,13 +92,15 @@ AliTPCTracklet::AliTPCTracklet(const AliTPCseed *track,Int_t sector,
 
 }
 
-AliTPCTracklet::AliTPCTracklet(const TObjArray &clusters,Int_t sector,
-                              TrackType type,Bool_t storeClusters) {
+AliTPCTracklet::AliTPCTracklet(const TObjArray &/*clusters*/,Int_t sector,
+                              TrackType /*type*/,Bool_t /*storeClusters*/)
+  : fNClusters(0),fNStoredClusters(0),fClusters(0),fSector(sector),fOuter(0),
+    fInner(0),fPrimary(0) {
   //TODO: write it!
 }
 
 AliTPCTracklet::AliTPCTracklet(const AliTPCTracklet &t)
-  : fNClusters(t.fNClusters),fNStoredClusters(t.fNStoredClusters),fClusters(0),
+  : TObject(t),fNClusters(t.fNClusters),fNStoredClusters(t.fNStoredClusters),fClusters(0),
     fSector(t.fSector),fOuter(0),fInner(0),
     fPrimary(0) {
   ////
@@ -120,7 +126,6 @@ AliTPCTracklet& AliTPCTracklet::operator=(const AliTPCTracklet &t) {
   ////
   if (this!=&t) {
     fNClusters=t.fNClusters;
-    fNStoredClusters=fNStoredClusters;
     delete fClusters;
     if (t.fClusters) {
       fClusters=new AliTPCclusterMI[t.fNStoredClusters];
@@ -176,77 +181,121 @@ AliTPCTracklet::~AliTPCTracklet() {
   delete fPrimary;
 }
 
-void AliTPCTracklet::FitKalman(const AliTPCseed *track,Int_t sector) {
-  AliTPCseed *t=new AliTPCseed(*track);
-  if (!t->Rotate(TMath::DegToRad()*(sector%18*20.+10.)-t->GetAlpha())) {
-    delete t;
+
+
+
+
+void AliTPCTracklet::FitKalman(const AliTPCseed *seed,Int_t sector) {
+  //
+  // Fit using Kalman filter
+  //
+  AliTPCseed *track=new AliTPCseed(*seed);
+  if (!track->Rotate(TMath::DegToRad()*(sector%18*20.+10.)-track->GetAlpha())) {
+    delete track;
     return;
   }
   // fit from inner to outer row
-  AliTPCseed *outerSeed=new AliTPCseed(*t);
-  Int_t n=0;
-  for (Int_t i=0;i<160;++i) {
-    AliTPCclusterMI *c=t->GetClusterPointer(i);
-    if (c&&c->GetDetector()==sector) {
-      if (n==1)        {
-       outerSeed->ResetCovariance(100.);
-      }
-      ++n;
-      Double_t r[3]={c->GetX(),c->GetY(),c->GetZ()};
-      Double_t cov[3]={0.1,0.,0.1}; //TODO: correct error parametrisation
-      if (!outerSeed->PropagateTo(r[0]) ||
-         !static_cast<AliExternalTrackParam*>(outerSeed)->Update(&r[1],cov)) {
-       delete outerSeed;
-       outerSeed=0;
-       break;
-      }
+  Double_t covar[15];
+  for (Int_t i=0;i<15;i++) covar[i]=0;
+  covar[0]=1.*1.;
+  covar[2]=1.*1.;
+  covar[5]=1.*1./(64.*64.);
+  covar[9]=1.*1./(64.*64.);
+  covar[14]=0;  // keep pt
+  Float_t xmin=1000, xmax=-10000;
+  Int_t imin=158, imax=0;
+  for (Int_t i=0;i<160;i++) {
+    AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (!c) continue;
+    if (c->GetDetector()!=sector)  continue;
+    if (c->GetX()<xmin) xmin=c->GetX();
+    if (c->GetX()>xmax) xmax=c->GetX();
+    if (i<imin) imin=i;
+    if (i>imax) imax=i;
+  }
+  if(imax-imin<10) {
+    delete track;
+    return;
+  }
+
+  for (Float_t x=track->GetX(); x<xmin; x++) track->PropagateTo(x);
+  track->AddCovariance(covar);
+  //
+  AliExternalTrackParam paramIn;
+  AliExternalTrackParam paramOut;
+  Bool_t isOK=kTRUE;
+  //
+  //
+  //
+  for (Int_t i=imin; i<=imax; i++){
+    AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (!c) continue;
+    Double_t r[3]={c->GetX(),c->GetY(),c->GetZ()};
+    Double_t cov[3]={0.01,0.,0.01}; //TODO: correct error parametrisation
+    AliTPCseed::GetError(c, track,cov[0],cov[2]);
+    cov[0]*=cov[0];
+    cov[2]*=cov[2];
+    if (!track->PropagateTo(r[0])) {
+      isOK=kFALSE;
+      break;
     }
+    if (RejectCluster(c,track)) continue;
+    if ( !((static_cast<AliExternalTrackParam*>(track)->Update(&r[1],cov)))) isOK=kFALSE;
   }
-  if (outerSeed)
-    fOuter=new AliExternalTrackParam(*outerSeed);
-  delete outerSeed;
-  // fit from outer to inner rows
-  AliTPCseed *innerSeed=new AliTPCseed(*t);
-  n=0;
-  for (Int_t i=159;i>=0;--i) {
-    AliTPCclusterMI *c=t->GetClusterPointer(i);
-    if (c&&c->GetDetector()==sector) {
-      if (n==1)        {
-       innerSeed->ResetCovariance(100.);
-      }
-      ++n;
-      Double_t r[3]={c->GetX(),c->GetY(),c->GetZ()};
-      Double_t cov[3]={0.1,0.,0.1};
-      if (!innerSeed->PropagateTo(r[0]) ||
-         !static_cast<AliExternalTrackParam*>(innerSeed)->Update(&r[1],cov)) {
-       delete innerSeed;
-       innerSeed=0;
-       break;
-      }
+  if (!isOK) { delete track; return;}
+  track->AddCovariance(covar);
+  //
+  //
+  //
+  for (Int_t i=imax; i>=imin; i--){
+    AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (!c) continue;
+    Double_t r[3]={c->GetX(),c->GetY(),c->GetZ()};
+    Double_t cov[3]={0.01,0.,0.01}; 
+    AliTPCseed::GetError(c, track,cov[0],cov[2]);
+    cov[0]*=cov[0];
+    cov[2]*=cov[2];
+    if (!track->PropagateTo(r[0])) {
+      isOK=kFALSE;
+      break;
     }
+    if (RejectCluster(c,track)) continue;
+    if ( !((static_cast<AliExternalTrackParam*>(track)->Update(&r[1],cov)))) isOK=kFALSE;
   }
-  if (innerSeed)
-    fInner=new AliExternalTrackParam(*innerSeed);
-  // propagate to the primary vertex
-  if (innerSeed) {
-    AliTPCseed *primarySeed=new AliTPCseed(*innerSeed);
-    Double_t pos[]={0.,0.,0.};
-    Double_t sigma[]={.1,.1,.1}; //TODO: is this correct?
-    AliESDVertex vertex(pos,sigma);
-    if (primarySeed->PropagateToVertex(&vertex))
-      fPrimary=new AliExternalTrackParam(*primarySeed);
-    delete primarySeed;
-    // for better comparison one does not want to have alpha changed...
-    if (!fPrimary->Rotate(fInner->GetAlpha())) {
-      delete fPrimary;
-      fPrimary=0;
+  if (!isOK) { delete track; return;}
+  paramIn = *track;
+  track->AddCovariance(covar);
+  //
+  //
+  for (Int_t i=imin; i<=imax; i++){
+    AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (!c) continue;
+    Double_t r[3]={c->GetX(),c->GetY(),c->GetZ()};
+    Double_t cov[3]={0.01,0.,0.01}; 
+    AliTPCseed::GetError(c, track,cov[0],cov[2]);
+    cov[0]*=cov[0];
+    cov[2]*=cov[2];
+    if (!track->PropagateTo(r[0])) {
+      isOK=kFALSE;
+      break;
     }
+    if (RejectCluster(c,track)) continue;
+    if ( !((static_cast<AliExternalTrackParam*>(track)->Update(&r[1],cov)))) isOK=kFALSE;
   }
-  delete innerSeed;
-
-  delete t;
+  if (!isOK) { delete track; return;}
+  paramOut=*track;
+  //
+  //
+  //
+  fOuter=new AliExternalTrackParam(paramOut);
+  fInner=new AliExternalTrackParam(paramIn);
+  //
+  delete track;
 }
 
+
+
+
 void AliTPCTracklet::FitLinear(const AliTPCseed *track,Int_t sector,
                               TrackType type) {
   TLinearFitter fy(1);
@@ -262,11 +311,15 @@ void AliTPCTracklet::FitLinear(const AliTPCseed *track,Int_t sector,
     fy.SetFormula("1 ++ x ++ x*x");
     fz.SetFormula("1 ++ x");
     break;
+  case kKalman:
+  case kRiemann:
+    break;
   }
   Double_t xmax=-1.;
   Double_t xmin=1000.;
   for (Int_t i=0;i<160;++i) {
     AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (c && RejectCluster(c)) continue;
     if (c&&c->GetDetector()==sector) {
       Double_t x=c->GetX();
       fy.AddPoint(&x,c->GetY());
@@ -376,6 +429,7 @@ void AliTPCTracklet::FitRiemann(const AliTPCseed *track,Int_t sector) {
   Double_t xmin=1000.;
   for (Int_t i=0;i<160;++i) {
     AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (c && RejectCluster(c)) continue;
     if (c&&c->GetDetector()==sector) {
       Double_t x=c->GetX();
       Double_t y=c->GetY();
@@ -409,6 +463,7 @@ void AliTPCTracklet::FitRiemann(const AliTPCseed *track,Int_t sector) {
   Double_t phi=0.;
   for (Int_t i=0;i<160;++i) {
     AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (c && RejectCluster(c)) continue;
     if (c&&c->GetDetector()==sector) {
       Double_t x=c->GetX();
       Double_t y=c->GetY();
@@ -439,8 +494,8 @@ void AliTPCTracklet::FitRiemann(const AliTPCseed *track,Int_t sector) {
     fOuter=new AliExternalTrackParam(xmax,alpha,p,c);
 }
 
-Bool_t AliTPCTracklet::Riemann2Helix(Double_t *a,Double_t *ca,
-                                    Double_t *b,Double_t *cb,
+Bool_t AliTPCTracklet::Riemann2Helix(Double_t *a,Double_t */*ca*/,
+                                    Double_t *b,Double_t */*cb*/,
                                     Double_t x0,
                                     Double_t *p,Double_t *c) {
   //TODO: signs!
@@ -450,7 +505,7 @@ Bool_t AliTPCTracklet::Riemann2Helix(Double_t *a,Double_t *ca,
   Double_t R=.5*TMath::Sqrt(4.*a[0]+a[1]*a[1]+a[2]*a[2]);
   Double_t dx=x0-xr0;
   if (dx*dx>=R*R) return kFALSE;
-  Double_t dy=TMath::Sqrt(R*R-dx*dx); //sign!!
+  Double_t dy=TMath::Sqrt((R-dx)*(R+dx)); //sign!!
   if (TMath::Abs(yr0+dy)>TMath::Abs(yr0-dy))
     dy=-dy;
   Double_t y0=yr0+dy; 
@@ -499,6 +554,7 @@ TObjArray AliTPCTracklet::CreateTracklets(const AliTPCseed *track,
   Int_t sectors[72]={0};
   for (Int_t i=0;i<160;++i) {
     AliTPCclusterMI *c=track->GetClusterPointer(i);
+    if (c && RejectCluster(c)) continue;
     if (c)
       ++sectors[c->GetDetector()];
   }
@@ -512,12 +568,15 @@ TObjArray AliTPCTracklet::CreateTracklets(const AliTPCseed *track,
   return tracklets;
 }
 
-TObjArray AliTPCTracklet::CreateTracklets(const TObjArray &clusters,
-                                         TrackType type,
-                                         Bool_t storeClusters,
-                                         Int_t minClusters,
-                                         Int_t maxTracklets) {
+TObjArray AliTPCTracklet::CreateTracklets(const TObjArray &/*clusters*/,
+                                         TrackType /*type*/,
+                                         Bool_t /*storeClusters*/,
+                                         Int_t /*minClusters*/,
+                                         Int_t /*maxTracklets*/) {
+  // TODO!
 
+  TObjArray tracklets;
+  return tracklets;
 }
 
 Bool_t AliTPCTracklet::PropagateToMeanX(const AliTPCTracklet &t1,
@@ -548,19 +607,25 @@ Bool_t AliTPCTracklet::PropagateToMeanX(const AliTPCTracklet &t1,
       t2m=new AliExternalTrackParam(*t2.GetOuter());
     }
     Double_t mx=.5*(t1m->GetX()+t2m->GetX());
-    Double_t b1,b2;
+    //Double_t b1,b2;
     Double_t xyz[3];
     t1m->GetXYZ(xyz);
-    b1=GetBz(xyz);
+    //b1=GetBz(xyz);
+    Double_t b1[3]; AliTracker::GetBxByBz(xyz,b1);
     t2m->GetXYZ(xyz);
-    b2=GetBz(xyz);
+    //b2=GetBz(xyz);
+    Double_t b2[3]; AliTracker::GetBxByBz(xyz,b2);
     if (t1m->Rotate(t2m->GetAlpha()) 
-       && t1m->PropagateTo(mx,b1) 
-       && t2m->PropagateTo(mx,b2));
+       //&& t1m->PropagateTo(mx,b1) 
+       //&& t2m->PropagateTo(mx,b2));
+       && t1m->PropagateToBxByBz(mx,b1) 
+       && t2m->PropagateToBxByBz(mx,b2));
     else
       if (t2m->Rotate(t1m->GetAlpha())
-         && t1m->PropagateTo(mx,b1) 
-         && t2m->PropagateTo(mx,b2));
+         //&& t1m->PropagateTo(mx,b1) 
+         //&& t2m->PropagateTo(mx,b2));
+         && t1m->PropagateToBxByBz(mx,b1) 
+         && t2m->PropagateToBxByBz(mx,b2));
       else {
        delete t1m;
        delete t2m;
@@ -573,11 +638,9 @@ Bool_t AliTPCTracklet::PropagateToMeanX(const AliTPCTracklet &t1,
   return t1m&&t2m;
 }
 
-double AliTPCTracklet::GetBz(Double_t *xyz) {
-  if (AliTracker::UniformField())
-    return AliTracker::GetBz();
-  else
-    return AliTracker::GetBz(xyz);
+double AliTPCTracklet::GetBz(Double_t *xyz) 
+{
+  return AliTracker::GetBz(xyz);
 }
 
 void AliTPCTracklet::RandomND(Int_t ndim,const Double_t *p,const Double_t *c,
@@ -642,7 +705,7 @@ void AliTPCTracklet::Test(const char* filename) {
     TEllipse e=AliTPCTracklet::ErrorEllipse(0.,0.,4.,1.,1.8);
     e.Draw();
  */
-  TTreeSRedirector ds("AliTPCTrackletDebug.root");
+  TTreeSRedirector ds(filename);
   Double_t p[5]={0.};
   Double_t c[15]={4.,
                  0.,4.,
@@ -692,3 +755,18 @@ void AliTPCTracklet::Test(const char* filename) {
   }
   */
 }
+
+
+Bool_t AliTPCTracklet::RejectCluster(AliTPCclusterMI* cl, AliExternalTrackParam * param){
+  //
+  // check the acceptance of cluster
+  // Cut on edge effects
+  //
+  Bool_t isReject = kFALSE;
+  Float_t edgeY = cl->GetX()*TMath::Tan(TMath::Pi()/18);
+  Float_t dist  = edgeY - TMath::Abs(cl->GetY());
+  if (param)  dist  = edgeY - TMath::Abs(param->GetY());
+  if (dist<fgEdgeCutY) isReject=kTRUE;
+  if (cl->GetType()<0) isReject=kTRUE;
+  return isReject;
+}