CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MultiTrackMassKinematicConstraint.cc
Go to the documentation of this file.
3 
4 
5 AlgebraicVector MultiTrackMassKinematicConstraint::value(const std::vector<KinematicState> &states,
6  const GlobalPoint& point) const
7 {
8  if(states.size()<nPart) throw VertexException("MultiTrackMassKinematicConstraint::not enough states given");
9 
10  double sumEnergy = 0, sumPx=0, sumPy=0., sumPz=0.;
11 
12  double a;
13  for (unsigned int i=0;i<nPart;++i) {
14  a = -states[i].particleCharge() * states[i].magneticField()->inInverseGeV(states[i].globalPosition()).z();
15 
16  sumEnergy += states[i].kinematicParameters().energy();
17  sumPx += states[i].kinematicParameters()(3) - a*(point.y() - states[i].kinematicParameters()(1));
18  sumPy += states[i].kinematicParameters()(4) + a*(point.x() - states[i].kinematicParameters()(0));
19  sumPz += states[i].kinematicParameters()(5);
20  }
21 
22  double j_m = sumPx*sumPx + sumPy*sumPy + sumPz*sumPz;
23 
24  AlgebraicVector res(1,0);
25  res(1) = sumEnergy*sumEnergy - j_m - mass*mass;
26  return res;
27 }
28 
30  const GlobalPoint& point) const
31 {
32  if(states.size()<nPart) throw VertexException("MultiTrackMassKinematicConstraint::not enough states given");
33 
34  AlgebraicMatrix res(1,states.size()*7,0);
35 
36  double sumEnergy = 0, sumPx=0, sumPy=0., sumPz=0.;
37 
38  double a;
39  for (unsigned int i=0;i<nPart;++i) {
40  a = -states[i].particleCharge() * states[i].magneticField()->inInverseGeV(states[i].globalPosition()).z();
41 
42  sumEnergy += states[i].kinematicParameters().energy();
43  sumPx += states[i].kinematicParameters()(3) - a*(point.y() - states[i].kinematicParameters()(1));
44  sumPy += states[i].kinematicParameters()(4) + a*(point.x() - states[i].kinematicParameters()(0));
45  sumPz += states[i].kinematicParameters()(5);
46  }
47 
48  for (unsigned int i=0;i<nPart;++i) {
49  a = -states[i].particleCharge() * states[i].magneticField()->inInverseGeV(states[i].globalPosition()).z();
50 
51  //x derivatives:
52  res(1,1+i*7) = 2*a*sumPy;
53 
54  //y derivatives:
55  res(1,2+i*7) = -2*a*sumPx;
56 
57  //z components:
58  res(1,3+i*7) = 0.;
59 
60  //px components:
61  res(1,4+i*7) = 2*sumEnergy/states[i].kinematicParameters().energy()*states[i].kinematicParameters()(3) - 2*sumPx;
62 
63  //py components:
64  res(1,5+i*7) = 2*sumEnergy/states[i].kinematicParameters().energy()*states[i].kinematicParameters()(4) - 2*sumPy;
65 
66  //pz1 components:
67  res(1,6+i*7) = 2*sumEnergy/states[i].kinematicParameters().energy()*states[i].kinematicParameters()(5) - 2*sumPz;
68 
69  //mass components:
70  res(1,7+i*7) = 2*states[i].kinematicParameters().mass()*sumEnergy/states[i].kinematicParameters().energy();
71  }
72  return res;
73 }
74 
76  const GlobalPoint& point) const
77 {
78  AlgebraicMatrix res(1,3,0);
79  if(states.size()<nPart) throw VertexException("MultiTrackMassKinematicConstraint::not enough states given");
80 
81  double sumA = 0, sumPx=0, sumPy=0.;
82 
83  double a;
84  for (unsigned int i=0;i<nPart;++i) {
85  a = -states[i].particleCharge() * states[i].magneticField()->inInverseGeV(states[i].globalPosition()).z();
86  sumA += a;
87 
88  sumPx += states[i].kinematicParameters()(3) - a*(point.y() - states[i].kinematicParameters()(1));
89  sumPy += states[i].kinematicParameters()(4) + a*(point.x() - states[i].kinematicParameters()(0));
90  }
91 
92  //xv component
93  res(1,1) = - 2 * sumPy * sumA;
94 
95  //yv component
96  res(1,2) = 2 * sumPx * sumA;
97 
98  //zv component
99  res(1,3) = 0.;
100 
101  return res;
102 }
int i
Definition: DBlmapReader.cc:9
virtual AlgebraicVector value(const std::vector< KinematicState > &states, const GlobalPoint &point) const
Common base class.
T y() const
Definition: PV3DBase.h:63
virtual AlgebraicMatrix positionDerivative(const std::vector< KinematicState > &states, const GlobalPoint &point) const
CLHEP::HepMatrix AlgebraicMatrix
virtual AlgebraicMatrix parametersDerivative(const std::vector< KinematicState > &states, const GlobalPoint &point) const
CLHEP::HepVector AlgebraicVector
double a
Definition: hdecay.h:121
T x() const
Definition: PV3DBase.h:62
*vegas h *****************************************************used in the default bin number in original ***version of VEGAS is ***a higher bin number might help to derive a more precise ***grade subtle point
Definition: invegas.h:5