CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | Private Attributes
KinematicConstrainedVertexUpdator Class Reference

#include <KinematicConstrainedVertexUpdator.h>

Public Member Functions

 KinematicConstrainedVertexUpdator ()
 
std::pair< std::pair
< std::vector< KinematicState >
, AlgebraicMatrix >
, RefCountedKinematicVertex
update (const AlgebraicVector &inState, const AlgebraicMatrix &inCov, const std::vector< KinematicState > &lStates, const GlobalPoint &lPoint, MultiTrackKinematicConstraint *cs) const
 
 ~KinematicConstrainedVertexUpdator ()
 

Private Attributes

VertexKinematicConstraintvConstraint
 
KinematicVertexFactoryvFactory
 

Detailed Description

Class caching the math part for KinematicConstrainedVertexFitter

Definition at line 13 of file KinematicConstrainedVertexUpdator.h.

Constructor & Destructor Documentation

KinematicConstrainedVertexUpdator::KinematicConstrainedVertexUpdator ( )
KinematicConstrainedVertexUpdator::~KinematicConstrainedVertexUpdator ( )

Definition at line 11 of file KinematicConstrainedVertexUpdator.cc.

References vConstraint, and vFactory.

12 {
13  delete vFactory;
14  delete vConstraint;
15 }

Member Function Documentation

std::pair< std::pair< std::vector< KinematicState >, AlgebraicMatrix >, RefCountedKinematicVertex > KinematicConstrainedVertexUpdator::update ( const AlgebraicVector inState,
const AlgebraicMatrix inCov,
const std::vector< KinematicState > &  lStates,
const GlobalPoint lPoint,
MultiTrackKinematicConstraint cs 
) const

Method updating the states. Takes a vector of full parameters: (x,y,z,particle_1,...,particle_n), corresponding linearization point: vector of states and GlobalPoint, and constraint to be applied during the vertex fit. Returns refitted vector of 7n+3 parameters and corresponding covariance matrix, where n - number of tracks.

Definition at line 18 of file KinematicConstrainedVertexUpdator.cc.

References g, i, edm::isNotFinite(), j, LogDebug, MultiTrackKinematicConstraint::numberOfEquations(), VertexKinematicConstraint::parametersDerivative(), MultiTrackKinematicConstraint::parametersDerivative(), VertexKinematicConstraint::positionDerivative(), MultiTrackKinematicConstraint::positionDerivative(), SpecificationBuilder_cfi::val(), VertexKinematicConstraint::value(), MultiTrackKinematicConstraint::value(), vConstraint, KinematicVertexFactory::vertex(), vFactory, PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by progressbar.ProgressBar::__next__(), MatrixUtil.Matrix::__setitem__(), MatrixUtil.Steps::__setitem__(), Vispa.Gui.VispaWidget.VispaWidget::autosize(), Vispa.Views.LineDecayView.LineDecayContainer::createObject(), Vispa.Views.LineDecayView.LineDecayContainer::deselectAllObjects(), Vispa.Gui.VispaWidgetOwner.VispaWidgetOwner::deselectAllWidgets(), Vispa.Gui.VispaWidget.VispaWidget::enableAutosizing(), progressbar.ProgressBar::finish(), KinematicConstrainedVertexFitter::fit(), Vispa.Gui.MenuWidget.MenuWidget::leaveEvent(), Vispa.Gui.VispaWidgetOwner.VispaWidgetOwner::mouseMoveEvent(), Vispa.Gui.MenuWidget.MenuWidget::mouseMoveEvent(), Vispa.Views.LineDecayView.LineDecayContainer::mouseMoveEvent(), Vispa.Gui.VispaWidgetOwner.VispaWidgetOwner::mouseReleaseEvent(), Vispa.Views.LineDecayView.LineDecayContainer::objectMoved(), MatrixUtil.Steps::overwrite(), Vispa.Views.LineDecayView.LineDecayContainer::removeObject(), Vispa.Gui.ConnectableWidget.ConnectableWidget::removePorts(), Vispa.Gui.FindDialog.FindDialog::reset(), Vispa.Gui.PortConnection.PointToPointConnection::select(), Vispa.Gui.VispaWidget.VispaWidget::select(), Vispa.Views.LineDecayView.LineDecayContainer::select(), Vispa.Gui.VispaWidget.VispaWidget::setText(), Vispa.Gui.VispaWidget.VispaWidget::setTitle(), Vispa.Gui.ZoomableWidget.ZoomableWidget::setZoom(), Vispa.Views.LineDecayView.LineDecayContainer::setZoom(), and Vispa.Gui.PortConnection.PointToPointConnection::updateConnection().

21 {
22  const MagneticField* field=lStates.front().magneticField();
23  AlgebraicMatrix d_matrix = vConstraint->parametersDerivative(lStates, lPoint);
24  AlgebraicMatrix e_matrix = vConstraint->positionDerivative(lStates, lPoint);
25  AlgebraicVector val_s = vConstraint->value(lStates, lPoint);
26  int vSize = lStates.size();
27 
28 //delta alpha
29  AlgebraicVector d_a(7*vSize + 3);
30  d_a(1) = lPoint.x();
31  d_a(2) = lPoint.y();
32  d_a(3) = lPoint.z();
33 
34  int cst = 0;
35  for(std::vector<KinematicState>::const_iterator i = lStates.begin(); i != lStates.end(); i++)
36  {
37  AlgebraicVector lst_par = asHepVector<7>(i->kinematicParameters().vector());
38  for(int j = 1; j<lst_par.num_row()+1; j++)
39  {d_a(3+7*cst+j) = lst_par(j);}
40  cst++;
41  }
42 
43 
44  AlgebraicVector delta_alpha = inPar - d_a;
45 
46 // cout<<"delta_alpha"<<delta_alpha<<endl;
47 //resulting matrix of derivatives and vector of values.
48 //their size depends of number of tracks to analyze and number of
49 //additional constraints to apply
52  if(cs == 0)
53  {
54 
55 //unconstrained vertex fitter case
56  g = AlgebraicMatrix(2*vSize,7*vSize+3,0);
57  val = AlgebraicVector(2*vSize);
58 
59 //filling the derivative matrix
60  g.sub(1,1,e_matrix);
61  g.sub(1,4,d_matrix);
62 
63 //filling the vector of values
64  val = val_s;
65  }else{
66 
67 //constrained vertex fitter case
68  int n_eq =cs->numberOfEquations();
69  g = AlgebraicMatrix(n_eq + 2*vSize,7*vSize + 3,0);
70  val = AlgebraicVector(n_eq + 2*vSize);
71  AlgebraicMatrix n_x = cs->positionDerivative(lStates, lPoint);
72  AlgebraicMatrix n_alpha = cs->parametersDerivative(lStates, lPoint);
73  AlgebraicVector c_val = cs->value(lStates, lPoint);
74 
75 //filling the derivative matrix
76 // cout<<"n_x:"<<n_x<<endl;
77 // cout<<"n_alpha"<<n_alpha<<endl;
78 
79  g.sub(1,1,n_x);
80  g.sub(1,4,n_alpha);
81  g.sub(n_eq+1, 1, e_matrix);
82  g.sub(n_eq+1, 4, d_matrix);
83 
84 //filling the vector of values
85  for(int i = 1;i< n_eq+1; i++)
86  {val(i) = c_val(i);}
87  for(int i = 1; i<(2*vSize+1); i++)
88  {val(i+n_eq) = val_s(i);}
89  }
90 
91  //check for NaN
92  for(int i = 1; i<=val.num_row();++i) {
93  if (edm::isNotFinite(val(i))) {
94  LogDebug("KinematicConstrainedVertexUpdator")
95  << "catched NaN.\n";
96  return std::pair<std::pair<std::vector<KinematicState>, AlgebraicMatrix>, RefCountedKinematicVertex >(
97  std::pair<std::vector<KinematicState>, AlgebraicMatrix>(std::vector<KinematicState>(), AlgebraicMatrix(1,0)),
99  }
100  }
101 
102 //debug feature
103  AlgebraicSymMatrix in_cov_sym(7*vSize + 3,0);
104 
105  for(int i = 1; i<7*vSize+4; ++i)
106  {
107  for(int j = 1; j<7*vSize+4; ++j)
108  {if(i<=j) in_cov_sym(i,j) = inCov(i,j);}
109  }
110 
111 //debug code
112  AlgebraicSymMatrix v_g_sym = in_cov_sym.similarity(g);
113 
114  int ifl1 = 0;
115  v_g_sym.invert(ifl1);
116  if(ifl1 !=0) {
117  LogDebug("KinematicConstrainedVertexFitter")
118  << "Fit failed: unable to invert SYM gain matrix\n";
119  return std::pair<std::pair<std::vector<KinematicState>, AlgebraicMatrix>, RefCountedKinematicVertex >(
120  std::pair<std::vector<KinematicState>, AlgebraicMatrix>(std::vector<KinematicState>(), AlgebraicMatrix(1,0)),
122  }
123 
124 // delta alpha is now valid!
125 //full math case now!
126  AlgebraicVector lambda = v_g_sym *(g*delta_alpha + val);
127 
128 //final parameters
129  AlgebraicVector finPar = inPar - in_cov_sym * g.T() * lambda;
130 
131 //covariance matrix business:
132  AlgebraicMatrix mFactor = in_cov_sym *(v_g_sym.similarityT(g))* in_cov_sym;
133 
134 //refitted covariance
135  AlgebraicMatrix rCov = in_cov_sym - mFactor;
136 
137 //symmetric covariance:
138  AlgebraicSymMatrix r_cov_sym(7*vSize+3,0);
139  for(int i = 1; i<7*vSize+4; ++i)
140  {
141  for(int j = 1; j<7*vSize+4; ++j)
142  {if(i<=j)r_cov_sym(i,j) = rCov(i,j);}
143  }
144 
145  AlgebraicSymMatrix pCov = r_cov_sym.sub(1,3);
146 
147 // chi2
148  AlgebraicVector chi = lambda.T()*(g*delta_alpha + val);
149 
150 //this is ndf without significant prior
151 //vertex so -3 factor exists here
152  float ndf = 2*vSize - 3;
153  if(cs != 0){ndf += cs->numberOfEquations();}
154 
155 
156 //making resulting vertex
157  GlobalPoint vPos (finPar(1),finPar(2),finPar(3));
158  VertexState st(vPos,GlobalError( asSMatrix<3>(pCov)));
159  RefCountedKinematicVertex rVtx = vFactory->vertex(st,chi(1),ndf);
160 
161 //making refitted states of Kinematic Particles
162  int i_int = 0;
163  std::vector<KinematicState> ns;
164  for(std::vector<KinematicState>::const_iterator i_st=lStates.begin(); i_st != lStates.end(); i_st++)
165  {
166  AlgebraicVector7 newPar;
167  for(int i =0; i<7; i++)
168  {newPar(i) = finPar(4 + i_int*7 + i);}
169 
170  AlgebraicSymMatrix nCovariance = r_cov_sym.sub(4 + i_int*7, 10 + i_int*7);
171  TrackCharge chl = i_st->particleCharge();
172  KinematicParameters nrPar(newPar);
173  KinematicParametersError nrEr(asSMatrix<7>(nCovariance));
174  KinematicState newState(nrPar,nrEr,chl, field);
175  ns.push_back(newState);
176  i_int++;
177  }
178  std::pair<std::vector<KinematicState>, AlgebraicMatrix> ns_m(ns,rCov);
179  return std::pair<std::pair<std::vector<KinematicState>, AlgebraicMatrix>, RefCountedKinematicVertex >(ns_m,rVtx);
180 }
#define LogDebug(id)
virtual AlgebraicMatrix positionDerivative(const std::vector< KinematicState > &, const GlobalPoint &) const =0
int i
Definition: DBlmapReader.cc:9
T y() const
Definition: PV3DBase.h:63
virtual AlgebraicVector value(const std::vector< KinematicState > &, const GlobalPoint &) const =0
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
GlobalErrorBase< double, ErrorMatrixTag > GlobalError
Definition: GlobalError.h:13
int TrackCharge
Definition: TrackCharge.h:4
bool isNotFinite(T x)
Definition: isFinite.h:10
ROOT::Math::SVector< double, 7 > AlgebraicVector7
Definition: Matrices.h:8
CLHEP::HepMatrix AlgebraicMatrix
virtual AlgebraicVector value(const std::vector< KinematicState > &states, const GlobalPoint &point) const
T z() const
Definition: PV3DBase.h:64
virtual int numberOfEquations() const =0
int j
Definition: DBlmapReader.cc:9
CLHEP::HepVector AlgebraicVector
ReferenceCountingPointer< KinematicVertex > RefCountedKinematicVertex
CLHEP::HepSymMatrix AlgebraicSymMatrix
virtual AlgebraicMatrix parametersDerivative(const std::vector< KinematicState > &, const GlobalPoint &) const =0
static RefCountedKinematicVertex vertex(const VertexState &state, float totalChiSq, float degreesOfFr)
T x() const
Definition: PV3DBase.h:62
virtual AlgebraicMatrix parametersDerivative(const std::vector< KinematicState > &states, const GlobalPoint &point) const
virtual AlgebraicMatrix positionDerivative(const std::vector< KinematicState > &states, const GlobalPoint &point) const

Member Data Documentation

VertexKinematicConstraint* KinematicConstrainedVertexUpdator::vConstraint
private
KinematicVertexFactory* KinematicConstrainedVertexUpdator::vFactory
private