Better formulas for the track seeding
authorbelikov <Iouri.Belikov@cern.ch>
Wed, 21 May 2014 16:17:38 +0000 (18:17 +0200)
committerbelikov <Iouri.Belikov@cern.ch>
Wed, 21 May 2014 16:17:38 +0000 (18:17 +0200)
ITS/UPGRADE/AliITSUTrackerCooked.cxx

index 96c2343..fc25920 100644 (file)
@@ -19,7 +19,8 @@ ClassImp(AliITSUTrackerCooked)
 // Constants hardcoded for the moment:
 //************************************************
 // seed "windows" in z and phi: MakeSeeds
-const Double_t kzWin=0.33, kpWin=3.14/4;
+const Double_t kzWin=0.33;
+const Double_t kminPt=0.05;
 // Maximal accepted impact parameters for the seeds 
 const Double_t kmaxDCAxy=3.;
 const Double_t kmaxDCAz= 3.;
@@ -112,26 +113,25 @@ f1(Double_t x1, Double_t y1, Double_t x2, Double_t y2, Double_t x3, Double_t y3)
     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);
+    Double_t xr=TMath::Abs(d/(d*x1-a)), yr=TMath::Abs(d/(d*y1-b));
     
-    return -xr*yr/sqrt(xr*xr+yr*yr);
+    Double_t crv=xr*yr/sqrt(xr*xr+yr*yr);
+    if (d>0) crv=-crv;
+
+    return crv;
 }
 
 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
+    // Initial approximation of the x-coordinate of the 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);
+
+  Double_t k1=(y2-y1)/(x2-x1), k2=(y3-y2)/(x3-x2);
+  Double_t x0=0.5*(k1*k2*(y1-y3) + k2*(x1+x2) - k1*(x2+x3))/(k2-k1);
+
+  return x0;
 }
 
 static Double_t 
@@ -166,11 +166,11 @@ AddCookedSeed(const Float_t r1[3], Int_t l1, Int_t i1,
     par[0]=y3;
     par[1]=z3;
     Double_t crv=f1(x1, y1, x2, y2, x3, y3); //curvature
-    Double_t cx0=f2(x1, y1, x2, y2, x3, y3); //curvature*x0
+    Double_t x0 =f2(x1, y1, x2, y2, x3, y3); //x-coordinate of the center
     Double_t tgl12=f3(x1, y1, x2, y2, z1, z2);
     Double_t tgl23=f3(x2, y2, x3, y3, z2, z3);
 
-    Double_t sf=x*crv - cx0;
+    Double_t sf=crv*(x-x0);
     if (TMath::Abs(sf) >= kAlmost1) return kFALSE;
     par[2]=sf;
 
@@ -211,7 +211,7 @@ Int_t AliITSUTrackerCooked::MakeSeeds() {
    if (fSeeds) {fSeeds->Delete(); delete fSeeds;}
    fSeeds=new TObjArray(77777);
 
-   Double_t zv=GetZ();
+   const Double_t zv=GetZ();
 
    AliITSUlayer &layer1=fgLayers[kSeedingLayer1];
    AliITSUlayer &layer2=fgLayers[kSeedingLayer2];
@@ -220,6 +220,9 @@ Int_t AliITSUTrackerCooked::MakeSeeds() {
    Double_t r2=layer2.GetR();
    Double_t r3=layer3.GetR();
 
+   const Double_t maxC  = TMath::Abs(GetBz()*kB2C/kminPt);
+   const Double_t kpWin = TMath::ASin(0.5*maxC*r1) - TMath::ASin(0.5*maxC*r2);
+
    Int_t nClusters1=layer1.GetNumberOfClusters();
    Int_t nClusters2=layer2.GetNumberOfClusters();
    Int_t nClusters3=layer3.GetNumberOfClusters();
@@ -246,8 +249,9 @@ Int_t AliITSUTrackerCooked::MakeSeeds() {
          if (TMath::Abs(phi2-phi1) > kpWin) continue;  //check in Phi
  
          Double_t zr3=z1 + (r3-r1)/(r2-r1)*(z2-z1);
-         Double_t d13=r1-r3, d32=r3-r2, d=r1-r2;
-         Double_t phir3=d32/d*phi1 + d13/d*phi2;  // FIXME
+         Double_t crv=f1(xyz1[0], xyz1[1], xyz2[0], xyz2[1], GetX(), GetY());
+         Double_t phir3 = phi1 + 0.5*crv*(r3 - r1); 
+
          Int_t start3=layer3.FindClusterIndex(zr3-kzWin/2);
          for (Int_t n3=start3; n3<nClusters3; n3++) {
              AliCluster *c3=layer3.GetCluster(n3);
@@ -255,8 +259,8 @@ Int_t AliITSUTrackerCooked::MakeSeeds() {
              //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(phir3-phi3) > kpWin/100) continue;  //check in Phi