bouchet's blog

void   StDcaGeometry::GetXYZ(Double_t xyzp[6], Double_t CovXyzp[21]) const {

  static const Float_t one = 1;

  Double_t sinP = TMath::Sin(mPsi);

  Double_t cosP = TMath::Cos(mPsi);

  Double_t pT   = pt();

  xyzp[0] = - mImp*sinP; // x 

  xyzp[1] =   mImp*cosP; // y

  xyzp[2] =   mZ;        // z

  xyzp[3] =   pT*cosP;   // px

  xyzp[4] =   pT*sinP;   // py

  xyzp[5] =   pT*mTan;   // pz

  Double_t dpTdPti = -pT*pT*TMath::Sign(one,mPti);

  Double_t f[30] = {

    //mImp,mZ,       mPsi,         mPti, mTan

    -sinP,  0, -mImp*cosP,            0,    0

    ,cosP,  0, -mImp*sinP,            0,    0

    ,   0,  1,          0,            0,    0

    ,   0,  0,   -pT*sinP, dpTdPti*cosP,    0

    ,   0,  0,    pT*cosP, dpTdPti*sinP,    0

    ,   0,  0,          0, dpTdPti*mTan,   pT

  };

  TRMatrix F(6,5,f);

  TRSymMatrix C(5,errMatrix());

  TRSymMatrix Cov(F,TRArray::kAxSxAT,C);

  TCL::ucopy(Cov.GetArray(),CovXyzp,21);

}

$STAR/StRoot/StarRoot/THelixTrack.cxx : DCA and error from THelixTrack

double THelixTrack::Dca(const double point[3]

                       ,double &dcaXY,double &dcaZ,double dcaEmx[3],int kind) const

/// Full 3d dca evaluation

/// point[3] - x,y,z of vertex

/// dcaXY - dca in xy plane

/// dcaZ  - dca in Z direction

/// dcaEmx[3] - err(dcaXY*dcaXY),err(dcaXY*dcaZ),err(dcaZ*dcaZ)

/// kind - 3=3d dca,2=2d dca

/// return distance to dca point

{

   double dif[3];

   double s = 0;

   assert(kind==2 || kind==3);

   if (kind==3) s = Path(point);

   else         s = Path(point[0],point[1]);

 

   THelixTrack th(*this);

   th.Move(s);

   const double *x=th.Pos();

   const double *d=th.Dir();

 

   for (int i=0;i<3;i++) {dif[i]=x[i]-point[i];}

   double nor = th.GetCos();

   double T[3][3]={{-d[1]/nor, d[0]/nor,    0}

                  ,{        0,        0,    1}

                  ,{ d[0]/nor, d[1]/nor,    0}};

 

   dcaXY = T[0][0]*dif[0]+T[0][1]*dif[1];

   dcaZ  = dif[2];

   THEmx_t *emx =th.Emx();

   dcaEmx[0] = emx->mHH;

   dcaEmx[1] = 0;

//      cos(Lambda) **4 to account that we are in the nearest point

   dcaEmx[2] = emx->mZZ*pow(th.GetCos(),4);

   return s;

}

$STAR/StRoot/StEvent/StDcaGeometry.h : signification of the track parameters

   /// signed impact parameter; Signed in such a way that:

    ///     x =  -impact*sin(Psi)

    ///     y =   impact*cos(Psi)

    Float_t  mImp;

    ///  Z-coordinate of this track (reference plane)

    Float_t  mZ;

    ///  Psi angle of the track

    Float_t  mPsi;

    /// signed invert pt [sign = sign(-qB)]

    Float_t  mPti;

    /// tangent of the track momentum dip angle

    Float_t  mTan;

    /// signed curvature

    Float_t  mCurv;

    

    /// pars errors

    Float_t  mImpImp;

    Float_t  mZImp, mZZ;

    Float_t  mPsiImp, mPsiZ, mPsiPsi;

    Float_t  mPtiImp, mPtiZ, mPtiPsi, mPtiPti;

    Float_t  mTanImp, mTanZ, mTanPsi, mTanPti, mTanTan;

Comparison of DCA distribution.

They were 3 methods to get the DCA :

by recalculation (method 1):

    TVector3 dirG;
    dirG = globP.Unit();
    Double_t cosl = dirG.Perp();
    Double_t cosp = dirG.x()/cosl;
    Double_t sinp = dirG.y()/cosl;
    DCA[0][0] = sinp * dcaG.x() - cosp * dcaG.y(); // switched sign
    DCA[1][0] = dcaG.z()/(cosl*cosl);

2. with the cov matrix of StDcaGeometry w/o moving the helix to the vertex (method 2):

CovGlobTrack_mTanPti[lgc],CovGlobTrack_mTanTan[lgc]};
    dcaGeometryC.set(parsKC, errsKC);
    THelixTrack     thelixKC = dcaGeometryC.thelix();
    Double_t ermx[3];
    thelixKC.Dca(vtx,parsKC[0],parsKC[1],ermx,2);
    DCA[0][1] = parsKC[0];
    DCA[1][1] = parsKC[1];

3. by calling double THelixTrack::Dca(const double point[3] at the vertex point (method 3):

    dcaGeometry.set(pars, errs);
    Double_t vtx[3] = {PrimVertexX[l], PrimVertexY[l], PrimVertexZ[l]};

    THelixTrack     thelix = dcaGeometry.thelix();
    StThreeVectorF Origin = dcaGeometry.origin();
    thelix.Move(thelix.Path(vtx));
    Double_t dcaEmx[3];
    Double_t dcaH[2];
    Double_t sigmaXY[2];
    Double_t sigmaZ[2];
    // dcaXY2 and dcaZ2 : via THelixTrack after moving the helix to the vertex
    DCA[0][2] = dcaGeometry.impact();
    DCA[1][2] = dcaGeometry.z();
    thelix.Dca(vtx,DCA[0][2],DCA[1][2],dcaEmx,2);

The plots :