LCFIVertex  0.7.2
vertexfunctionclassic.cpp
1 
2 #include "../include/vertexfunctionclassic.h"
3 
4 #include "../include/gausstube.h"
5 #include "../include/gaussellipsoid.h"
6 #include "../include/interactionpoint.h"
7 #include "../../util/inc/vector3.h"
8 
9 namespace vertex_lcfi { namespace ZVTOP
10 {
11  VertexFunctionClassic::VertexFunctionClassic(std::vector<Track*> & Tracks, const double Kip, const double Kalpha, const Vector3 & JetAxis)
12  {
13  _Kip=Kip;
14  _Kalpha=Kalpha;
15  _JetAxis=JetAxis;
16 
17  for (std::vector<Track*>::iterator iTrack = Tracks.begin();iTrack != Tracks.end();++iTrack)
18  {
19  GaussTube* element= new GaussTube(*iTrack);
20  _AllElements.push_back(element);
21  _Tubes.push_back(element);
22  _ElementsNewedByThis.push_back(element);
23  }
24  _Ellipsoid=0;
25 
26  }
27 
28  VertexFunctionClassic::VertexFunctionClassic(std::vector<Track*> & Tracks, InteractionPoint* IP, const double Kip, const double Kalpha, const Vector3 & JetAxis)
29  {
30  _Kip=Kip;
31  _Kalpha=Kalpha;
32  _JetAxis=JetAxis;
33  for (std::vector<Track*>::iterator iTrack = Tracks.begin();iTrack != Tracks.end();++iTrack)
34  {
35  GaussTube* element= new GaussTube(*iTrack);
36  _AllElements.push_back(element);
37  _Tubes.push_back(element);
38  _ElementsNewedByThis.push_back(element);
39  }
40  if (IP)
41  {
42  GaussEllipsoid* element= new GaussEllipsoid(IP);
43  _AllElements.push_back(element);
44  _Ellipsoid = element;
45  _ElementsNewedByThis.push_back(element);
46  }
47  else
48  _Ellipsoid = 0;
49 
50  }
51 
52 
54  {
55  for (std::vector<VertexFunctionElement*>::iterator iElement = _ElementsNewedByThis.begin();iElement != _ElementsNewedByThis.end();++iElement)
56  delete *iElement;
57  }
58 
59 
60 
61  double VertexFunctionClassic::valueAt(const Vector3 & Point) const
62  {
63  double SumOfTubes = 0;
64  double SumOfSquaredTubes = 0;
65  double dlong = 0;
66  double dmag = 0;
67  //TODO make other constants parameters
68 
69  //Now add up the tubes
70  for (std::vector<GaussTube*>::const_iterator iTube = _Tubes.begin();iTube != _Tubes.end();++iTube)
71  {
72  double Tube = (*iTube)->valueAt(Point);
73  SumOfTubes += Tube;
74  SumOfSquaredTubes += (Tube*Tube);
75  }
76 
77  //And IP if we have one
78  double IPValue = 0;
79  if (_Ellipsoid)
80  IPValue = _Ellipsoid->valueAt(Point);
81 
82  double dtran = 0;
83  //Work out the distance to the jet axis and how far along the jet axis we are
84  if (_Ellipsoid)
85  {
86  dlong = Point.subtract(_Ellipsoid->ip()->position()).dot(_JetAxis) / _JetAxis.mag();
87 
88  if (dlong < -0.01) //100 Micron behind the ip
89  return -1.0;
90 
91  dmag = _Ellipsoid->ip()->distanceTo(Point);
92 
93  dtran = sqrt(dmag*dmag - dlong*dlong);
94  }
95  //Check if we are outside tube and add on kalpha adjustment if not
96  if (SumOfTubes > 0)
97  if (dtran > 0.005) //50 Micron
98  {
99  double alpha = acos((dlong + 0.01) / sqrt(((dlong + 0.01)*(dlong + 0.01)) + (dtran - 0.005)*(dtran - 0.005)));
100  return exp(-_Kalpha*alpha*alpha)*((_Kip*IPValue) + SumOfTubes - ( ((_Kip*IPValue*IPValue)+SumOfSquaredTubes) / ((_Kip*IPValue)+SumOfTubes)));
101  }
102  else
103  return (_Kip*IPValue) + SumOfTubes - ( ((_Kip*IPValue*IPValue)+SumOfSquaredTubes) / ((_Kip*IPValue)+SumOfTubes));
104  else
105  return 0;
106  }
107 
109  {
110  // TODO Implement if needed this funcion could be one class up
111  return Matrix3x3();
112  }
113 
114 
116  {
117  // TODO Implement if needed this funcion could be one class up
118  return Matrix3x3();
119  }
120 
121 
122 
123 }}
124 
double valueAt(const Vector3 &Point) const
Calculate the value of the ellipsoid at a point.
double valueAt(const Vector3 &Point) const
Find the value of the vertex function at Point.
const Vector3 & position() const
Return position of IP.
VertexFunctionClassic(std::vector< Track * > &Tracks, const double Kip, const double Kalpha, const Vector3 &JetAxis)
Constructor.
double distanceTo(const Vector3 &Point) const
Return distance from this interacion point to a point.
Matrix3x3 firstDervAt(const Vector3 &Point) const
Find the spacial derivative of the vertex function at Point (not implemented)
Gaussian tube component of the vertex function.
InteractionPoint * ip()
InteractionPoint object used.
Gaussian ellipsoid component of the vertex function.
Matrix3x3 secondDervAt(const Vector3 &Point) const
Find the 2nd spacial derivative of the vertex function at Point (not implemented) ...