VertexFitInterface/src/TrackParameters.cc

#include "TMath.h"
#include "MitEdm/DataFormats/interface/Types.h"
#include "MitEdm/VertexFitInterface/interface/TrackParameters.h"
using namespace std;
using namespace TMath;
using namespace reco;
using namespace mitedm;

TrackParameters::TrackParameters(const Track* trk, const TrackConvention cv, double bField) :
  iConvention_(cv),
  bField_     (bField),
  fCurv_      (0.5 * 0.0029979 * bField_ * -1.0)
{
  pars_.ResizeTo(5);
  cMat_.ResizeTo(5,5);
  if (cv == iCms) {
    pars_[0] = trk->qoverp();
    pars_[1] = trk->lambda();
    pars_[2] = trk->phi();
    pars_[3] = trk->dxy();
    pars_[4] = trk->dsz();
    for (int i=0; i<5; i++) {
      for (int j=0; j<5; j++)
        cMat_(i,j) = trk->covariance(i,j);
    }
  }
  else {
    cout << "TrackParameters::TrackParameters -- Requested transfer: CMS -> MVF convention"
         << " not yet implemented.\n";
  }
}

TrackParameters::TrackParameters(const TrackParameters &trk) :
  iConvention_(trk.iConvention_),
  bField_     (trk.bField_),
  fCurv_      (trk.fCurv_)
{
  pars_.ResizeTo(5);
  cMat_.ResizeTo(5,5);
  for (int i=0; i<5; i++) {
    pars_(i) = trk.pars(i);
    for (int j=0; j<5; j++)
      cMat_(i,j) = trk.cMat(i,j);
  }
}

TrackParameters TrackParameters::mvfTrack() const
{
  // Return a new set of track parameters in the required MVF convention
  TrackParameters outTk;
  outTk.pars_.ResizeTo(5);
  outTk.cMat_.ResizeTo(5,5);

  if (iConvention_ == iCms) {
    outTk.setPars(0,1.0/Tan(PiOver2() - pars_[1]));  // cotTheta
    outTk.setPars(1,fCurv_*pars_[0]/Cos(pars_[1]));  // curvature
    outTk.setPars(2,pars_[4]/Cos(pars_[1]));         // z0
    outTk.setPars(3,pars_[3]);                       // d0
    outTk.setPars(4,pars_[2]);                       // phi0
    if (outTk.pars()[4]<0.0)
      outTk.setPars(4,pars_[2] + 2*Pi());            // phi0 [0,2pi)
    // dPidQj(i,j) gives partial dP_i/dQ_j where Pi are the new parameters and Qj are the old ones
    TMatrixD dPidQj(5,5);
    dPidQj(0,1) = 1.0/(Sin(PiOver2() - pars_[1])*Sin(PiOver2() - pars_[1]));
    dPidQj(1,0) = fCurv_/Cos(pars_[1]);
    dPidQj(1,1) = fCurv_*pars_[0]*Tan(pars_[1])/Cos(pars_[1]);
    dPidQj(2,1) = pars_[4]*Tan(pars_[1])/Cos(pars_[1]);
    dPidQj(2,4) = 1.0/Cos(pars_[1]);
    dPidQj(3,3) = 1.0;
    dPidQj(4,2) = 1.0;
    for (int i=0; i<5; i++)
      for (int j=0; j<5; j++) {
        outTk.setCMat(i,j,0.0);
        for (int k=0; k<5; k++)
          for (int l=0; l<5; l++)
            outTk.addCMat(i,j,dPidQj(i,k)*dPidQj(j,l)*cMat_(k,l));
      }
  }
  else {
    outTk.setPars(0,pars_[0]);
    outTk.setPars(1,pars_[1]);
    outTk.setPars(2,pars_[2]);
    outTk.setPars(3,pars_[3]);
    outTk.setPars(4,pars_[4]);
    for (int i=0; i<5; i++)
      for (int j=0; j<5; j++) {
        outTk.setCMat(i,j,cMat_(i,j));
      }
  }
  return outTk;
}

TrackParameters TrackParameters::cmsTrack() const
{
  // Return a new set of track parameters in the required CMS convention
  TrackParameters outTk;
  outTk.pars_.ResizeTo(5);
  outTk.cMat_.ResizeTo(5,5);

  if (iConvention_ == iCms) {
    outTk.setPars(0,pars_[0]);
    outTk.setPars(1,pars_[1]);
    outTk.setPars(2,pars_[2]);
    outTk.setPars(3,pars_[3]);
    outTk.setPars(4,pars_[4]);
    for (int i=0; i<5; i++)
      for (int j=0; j<5; j++) {
        outTk.setCMat(i,j,cMat_(i,j));
      }
  }
  else {
    // Apply MVF to CMS conversion
    outTk.setPars(1,ATan(pars_[0]));                        // lambda
    outTk.setPars(0,pars_[1]*Cos(outTk.pars(1))/fCurv_);    // qoverp
    outTk.setPars(2,pars_[4]);                              // phi
    outTk.setPars(3,pars_[3]);                              // dxy
    outTk.setPars(4,Cos(outTk.pars(1))*pars_[2]);           // dsz
    if (outTk.pars(2)>Pi())
      outTk.setPars(2,pars_[4]-2*Pi());                     // phi0 [-pi,pi)
    // dPidQj(i,j) gives partial dP_i/dQ_j where Pi are the new parameters and Qj are the old ones
    TMatrixD dPidQj(5,5);
    dPidQj(1,0) = 1.0/(1.0+pars_[0]*pars_[0]);
    dPidQj(0,0) = pars_[1]*Sin(outTk.pars(1))*dPidQj(1,0)/fCurv_;
    dPidQj(0,1) = Cos(outTk.pars(1))/fCurv_;
    dPidQj(2,4) = 1.0;
    dPidQj(3,3) = 1.0;
    dPidQj(4,0) = -pars_[2]*Sin(outTk.pars(1))*dPidQj(1,0);
    dPidQj(4,2) = Cos(outTk.pars(1));
    for (int i=0; i<5; i++)
      for (int j=0; j<5; j++) {
        outTk.setCMat(i,j,0.0);
        for (int k=0; k<5; k++)
          for (int l=0; l<5; l++)
            outTk.addCMat(i,j,dPidQj(i,k)*dPidQj(j,l)*cMat_(k,l));
      }
  }
  return outTk;
}

void TrackParameters::print() const
{
  printf("\n==== Show track ====\n Track parameters with convention %1d.\n",int(iConvention_));
  for (int i=0; i<5; i++)
    printf(" Parameter [%1d]: %8.4f\n",i,pars_[i]);
  printf("Track covariance matrix");
  for (int i=0; i<5; i++) {
    printf("\n CV[%1d,0-4] ",i);
    for (int j=0; j<5; j++)
      printf(" %8.4f",cMat_(i,j));
  }
  printf("\n");
}
Revision:	1.1
Committed:	Thu Aug 28 22:09:20 2008 UTC (16 years, 8 months ago) by paus
Content type:	text/plain
Branch:	MAIN
Log Message:	Adding HelixIntersector interface.
#	User	Rev	Content
1	paus	1.1	#include "TMath.h"
2			#include "MitEdm/DataFormats/interface/Types.h"
3			#include "MitEdm/VertexFitInterface/interface/TrackParameters.h"
4			using namespace std;
5			using namespace TMath;
6			using namespace reco;
7			using namespace mitedm;
8
9			TrackParameters::TrackParameters(const Track* trk, const TrackConvention cv, double bField) :
10			iConvention_(cv),
11			bField_ (bField),
12			fCurv_ (0.5 * 0.0029979 * bField_ * -1.0)
13			{
14			pars_.ResizeTo(5);
15			cMat_.ResizeTo(5,5);
16			if (cv == iCms) {
17			pars_[0] = trk->qoverp();
18			pars_[1] = trk->lambda();
19			pars_[2] = trk->phi();
20			pars_[3] = trk->dxy();
21			pars_[4] = trk->dsz();
22			for (int i=0; i<5; i++) {
23			for (int j=0; j<5; j++)
24			cMat_(i,j) = trk->covariance(i,j);
25			}
26			}
27			else {
28			cout << "TrackParameters::TrackParameters -- Requested transfer: CMS -> MVF convention"
29			<< " not yet implemented.\n";
30			}
31			}
32
33			TrackParameters::TrackParameters(const TrackParameters &trk) :
34			iConvention_(trk.iConvention_),
35			bField_ (trk.bField_),
36			fCurv_ (trk.fCurv_)
37			{
38			pars_.ResizeTo(5);
39			cMat_.ResizeTo(5,5);
40			for (int i=0; i<5; i++) {
41			pars_(i) = trk.pars(i);
42			for (int j=0; j<5; j++)
43			cMat_(i,j) = trk.cMat(i,j);
44			}
45			}
46
47			TrackParameters TrackParameters::mvfTrack() const
48			{
49			// Return a new set of track parameters in the required MVF convention
50			TrackParameters outTk;
51			outTk.pars_.ResizeTo(5);
52			outTk.cMat_.ResizeTo(5,5);
53
54			if (iConvention_ == iCms) {
55			outTk.setPars(0,1.0/Tan(PiOver2() - pars_[1])); // cotTheta
56			outTk.setPars(1,fCurv_*pars_[0]/Cos(pars_[1])); // curvature
57			outTk.setPars(2,pars_[4]/Cos(pars_[1])); // z0
58			outTk.setPars(3,pars_[3]); // d0
59			outTk.setPars(4,pars_[2]); // phi0
60			if (outTk.pars()[4]<0.0)
61			outTk.setPars(4,pars_[2] + 2*Pi()); // phi0 [0,2pi)
62			// dPidQj(i,j) gives partial dP_i/dQ_j where Pi are the new parameters and Qj are the old ones
63			TMatrixD dPidQj(5,5);
64			dPidQj(0,1) = 1.0/(Sin(PiOver2() - pars_[1])*Sin(PiOver2() - pars_[1]));
65			dPidQj(1,0) = fCurv_/Cos(pars_[1]);
66			dPidQj(1,1) = fCurv_pars_[0]Tan(pars_[1])/Cos(pars_[1]);
67			dPidQj(2,1) = pars_[4]*Tan(pars_[1])/Cos(pars_[1]);
68			dPidQj(2,4) = 1.0/Cos(pars_[1]);
69			dPidQj(3,3) = 1.0;
70			dPidQj(4,2) = 1.0;
71			for (int i=0; i<5; i++)
72			for (int j=0; j<5; j++) {
73			outTk.setCMat(i,j,0.0);
74			for (int k=0; k<5; k++)
75			for (int l=0; l<5; l++)
76			outTk.addCMat(i,j,dPidQj(i,k)dPidQj(j,l)cMat_(k,l));
77			}
78			}
79			else {
80			outTk.setPars(0,pars_[0]);
81			outTk.setPars(1,pars_[1]);
82			outTk.setPars(2,pars_[2]);
83			outTk.setPars(3,pars_[3]);
84			outTk.setPars(4,pars_[4]);
85			for (int i=0; i<5; i++)
86			for (int j=0; j<5; j++) {
87			outTk.setCMat(i,j,cMat_(i,j));
88			}
89			}
90			return outTk;
91			}
92
93			TrackParameters TrackParameters::cmsTrack() const
94			{
95			// Return a new set of track parameters in the required CMS convention
96			TrackParameters outTk;
97			outTk.pars_.ResizeTo(5);
98			outTk.cMat_.ResizeTo(5,5);
99
100			if (iConvention_ == iCms) {
101			outTk.setPars(0,pars_[0]);
102			outTk.setPars(1,pars_[1]);
103			outTk.setPars(2,pars_[2]);
104			outTk.setPars(3,pars_[3]);
105			outTk.setPars(4,pars_[4]);
106			for (int i=0; i<5; i++)
107			for (int j=0; j<5; j++) {
108			outTk.setCMat(i,j,cMat_(i,j));
109			}
110			}
111			else {
112			// Apply MVF to CMS conversion
113			outTk.setPars(1,ATan(pars_[0])); // lambda
114			outTk.setPars(0,pars_[1]*Cos(outTk.pars(1))/fCurv_); // qoverp
115			outTk.setPars(2,pars_[4]); // phi
116			outTk.setPars(3,pars_[3]); // dxy
117			outTk.setPars(4,Cos(outTk.pars(1))*pars_[2]); // dsz
118			if (outTk.pars(2)>Pi())
119			outTk.setPars(2,pars_[4]-2*Pi()); // phi0 [-pi,pi)
120			// dPidQj(i,j) gives partial dP_i/dQ_j where Pi are the new parameters and Qj are the old ones
121			TMatrixD dPidQj(5,5);
122			dPidQj(1,0) = 1.0/(1.0+pars_[0]*pars_[0]);
123			dPidQj(0,0) = pars_[1]Sin(outTk.pars(1))dPidQj(1,0)/fCurv_;
124			dPidQj(0,1) = Cos(outTk.pars(1))/fCurv_;
125			dPidQj(2,4) = 1.0;
126			dPidQj(3,3) = 1.0;
127			dPidQj(4,0) = -pars_[2]Sin(outTk.pars(1))dPidQj(1,0);
128			dPidQj(4,2) = Cos(outTk.pars(1));
129			for (int i=0; i<5; i++)
130			for (int j=0; j<5; j++) {
131			outTk.setCMat(i,j,0.0);
132			for (int k=0; k<5; k++)
133			for (int l=0; l<5; l++)
134			outTk.addCMat(i,j,dPidQj(i,k)dPidQj(j,l)cMat_(k,l));
135			}
136			}
137			return outTk;
138			}
139
140			void TrackParameters::print() const
141			{
142			printf("\n==== Show track ====\n Track parameters with convention %1d.\n",int(iConvention_));
143			for (int i=0; i<5; i++)
144			printf(" Parameter [%1d]: %8.4f\n",i,pars_[i]);
145			printf("Track covariance matrix");
146			for (int i=0; i<5; i++) {
147			printf("\n CV[%1d,0-4] ",i);
148			for (int j=0; j<5; j++)
149			printf(" %8.4f",cMat_(i,j));
150			}
151			printf("\n");
152			}