14 #include "../include/vertexfitterlsm.h"
16 #include "../include/maxminfinder.h"
17 #include "../include/candidatevertex.h"
18 #include "../../inc/trackstate.h"
19 #include "../../inc/track.h"
20 #include "../include/interactionpoint.h"
21 #include "../../util/inc/matrix.h"
22 #include "../../util/inc/vector3.h"
24 namespace vertex_lcfi {
namespace ZVTOP
26 VertexFitterLSM::VertexFitterLSM()
27 :_UseManualSeed(false),_InitialStep(100.0/1000.0)
30 void VertexFitterLSM::fitVertex(
const std::vector<TrackState*> & Tracks, InteractionPoint* IP,
Vector3 & Result)
32 _trackStateList=Tracks;
48 int N = Tracks.size();
49 for (
int OuterIndex=0;OuterIndex < N-1;++OuterIndex)
51 for (
int InnerIndex=OuterIndex+1;InnerIndex < N;++InnerIndex)
53 Tracks[OuterIndex]->swimToStateNearest(Tracks[InnerIndex]);
54 Tracks[InnerIndex]->swimToStateNearest(Tracks[OuterIndex]->position());
55 Vector3 LeftPoint = Tracks[OuterIndex]->position();
56 Vector3 MidPoint = Tracks[OuterIndex]->position() + ((Tracks[InnerIndex]->position() - Tracks[OuterIndex]->position()) * 0.5);
57 Vector3 RightPoint = Tracks[InnerIndex]->position();
59 double p = (RightPoint.subtract(LeftPoint)).mag();
60 if (p > 0.0001/1000.0)
64 double ChiLeft = _chi2Contribution( MidPoint, Tracks[OuterIndex]);
66 double ChiRight = _chi2Contribution( MidPoint, Tracks[InnerIndex]);
68 double Sigma2Left = ((0.5*p)*(0.5*p))/ChiLeft;
69 double Sigma2Right = ((-0.5*p)*(-0.5*p))/ChiRight;
70 double Frac = Sigma2Left/(Sigma2Left+Sigma2Right);
71 Seed = Seed + (LeftPoint + ((RightPoint-LeftPoint)*Frac));
74 Seed = Seed + LeftPoint;
84 Tracks[0]->swimToStateNearest(IP->position());
85 Vector3 LeftPoint = IP->position();
86 Vector3 MidPoint = IP->position() + ((Tracks[0]->position() - IP->position()) * 0.5);
87 Vector3 RightPoint = Tracks[0]->position();
88 double p = (RightPoint.subtract(LeftPoint)).mag();
89 if (p > 0.0001/1000.0)
91 double ChiLeft = _chi2Contribution( MidPoint, IP);
92 double ChiRight = _chi2Contribution( MidPoint, Tracks[0]);
93 double Sigma2Left = ((0.5*p)*(0.5*p))/ChiLeft;
94 double Sigma2Right = ((-0.5*p)*(-0.5*p))/ChiRight;
95 double Frac = Sigma2Left/(Sigma2Left+Sigma2Right);
96 Seed = Seed + (LeftPoint + ((RightPoint-LeftPoint)*Frac));
99 Seed = Seed + LeftPoint;
102 if (Tracks.size()==0)
106 Seed = IP->position();
113 if (Tracks.size()>2 || (Tracks.size()>1 && IP))
123 ZVTOP::FunctionMinimiser<ZVTOP::VertexFitterLSM> minimiser(
this, _InitialStep, 6 );
124 std::vector<double> result = minimiser.Minimise(Seed.stlVector());
125 Result.x() = result[0];
126 Result.y() = result[1];
127 Result.z() = result[2];
135 void VertexFitterLSM::fitVertex(
const std::vector<TrackState*> & Tracks, InteractionPoint* IP,
Vector3 & Result,
double & ChiSquaredOfFit)
137 this->fitVertex(Tracks,IP,Result);
140 for( std::vector<TrackState*>::const_iterator i=Tracks.begin(); i<Tracks.end(); i++ )
142 ChiSquaredOfFit += _chi2Contribution( Result, (*i) );
145 ChiSquaredOfFit += _chi2Contribution( Result, IP );
148 void VertexFitterLSM::fitVertex(
const std::vector<TrackState*> & Tracks, InteractionPoint* IP,
Vector3 & Result,
double & ChiSquaredOfFit, std::map<TrackState*,double> & ChiSquaredOfTrack,
double & ChiSquaredOfIP)
151 this->fitVertex(Tracks,IP,Result);
154 ChiSquaredOfTrack.clear();
155 for( std::vector<TrackState*>::const_iterator i=Tracks.begin(); i<Tracks.end(); i++ )
157 double chi = _chi2Contribution( Result, (*i) );
158 ChiSquaredOfTrack.insert( std::pair<TrackState*,double>( (*i), chi ) );
159 ChiSquaredOfFit += chi;
162 ChiSquaredOfIP=_chi2Contribution( Result, IP );
163 ChiSquaredOfFit += ChiSquaredOfIP;
166 void VertexFitterLSM::fitVertex(
const std::vector<TrackState*> & Tracks, InteractionPoint* IP,
Vector3 & Result,
Matrix3x3 & ResultError,
double & ChiSquaredOfFit, std::map<TrackState*,double> & ChiSquaredOfTrack,
double & ChiSquaredOfIP)
168 this->fitVertex(Tracks,IP,Result,ChiSquaredOfFit,ChiSquaredOfTrack,ChiSquaredOfIP);
173 for( std::vector<TrackState*>::const_iterator i=Tracks.begin(); i<Tracks.end(); i++ )
175 ResultError += (*i)->vertexErrorContribution(Result);
177 ResultError = InvertMatrix(ResultError);
180 if (IP) ResultError = IP->errorMatrix();
183 double VertexFitterLSM::valueAt(
const Vector3 & point)
188 for( std::vector<TrackState*>::iterator i=_trackStateList.begin(); i<_trackStateList.end(); i++ )
190 chi2+=_chi2Contribution( point, (*i) );
195 if (_trackStateList.size()<2)
196 chi2+=_chi2Contribution( point, _ip );
201 double VertexFitterLSM::valueAt(
const std::vector<double> & point)
203 return this->valueAt(
Vector3(point[0],point[1],point[2]));
206 void VertexFitterLSM::setSeed(
Vector3 Seed)
208 _UseManualSeed =
true;
212 void VertexFitterLSM::setInitialStep(
double Step)
217 double VertexFitterLSM::_chi2Contribution(
const Vector3 & point, TrackState* pTrackState )
219 if( 0==pTrackState )
return 0;
220 else return pTrackState->chi2(point);
223 double VertexFitterLSM::_chi2Contribution(
const Vector3 & point, InteractionPoint* pIP )
225 if( 0==pIP )
return 0;
226 else return pIP->chi2(point);