#include <KinematicConstrainedVertexFitter.h>

Public Member Functions
RefCountedKinematicTree	fit (const std::vector< RefCountedKinematicParticle > &part)

RefCountedKinematicTree	fit (const std::vector< RefCountedKinematicParticle > &part, MultiTrackKinematicConstraint *cs)

RefCountedKinematicTree	fit (const std::vector< RefCountedKinematicParticle > &part, MultiTrackKinematicConstraint cs, GlobalPoint pt)

float	getCSum () const

int	getNit () const

	KinematicConstrainedVertexFitter ()

	KinematicConstrainedVertexFitter (const LinearizationPointFinder &fnd)

void	setParameters (const edm::ParameterSet &pSet)

	~KinematicConstrainedVertexFitter ()

Private Member Functions
void	defaultParameters ()

Private Attributes
float	csum

LinearizationPointFinder *	finder

int	iterations

ConstrainedTreeBuilder *	tBuilder

float	theMaxDelta

float	theMaxReducedChiSq

int	theMaxStep

float	theMinChiSqImprovement

KinematicConstrainedVertexUpdator *	updator

VertexKinematicConstraint *	vCons

Detailed Description

Class fitting the veretx out of set of tracks via usual LMS with Lagrange multipliers. Additional constraints can be applyed to the tracks during the vertex fit (solves non-factorizabele cases). Since the vertex constraint is included by default, do not add a separate VertexKinematicConstraint! Example: Vertex fit with collinear tracks..

Definition at line 21 of file KinematicConstrainedVertexFitter.h.

Constructor & Destructor Documentation

KinematicConstrainedVertexFitter::KinematicConstrainedVertexFitter ( )

Default constructor using LMSLinearizationPointFinder

Definition at line 7 of file KinematicConstrainedVertexFitter.cc.

References csum, defaultParameters(), finder, iterations, tBuilder, updator, and vCons.

                                                                    {
   finder = new DefaultLinearizationPointFinder();
   vCons = new VertexKinematicConstraint();
   updator = new KinematicConstrainedVertexUpdator();
   tBuilder = new ConstrainedTreeBuilder;
   defaultParameters();
   iterations = -1;
   csum = -1000.0;
 }

KinematicConstrainedVertexFitter::KinematicConstrainedVertexFitter ( const LinearizationPointFinder & fnd )

Constructor with user-provided LinearizationPointFinder

Definition at line 17 of file KinematicConstrainedVertexFitter.cc.

References LinearizationPointFinder::clone(), csum, defaultParameters(), finder, iterations, tBuilder, updator, and vCons.

                                                                                                       {
   finder = fnd.clone();
   vCons = new VertexKinematicConstraint();
   updator = new KinematicConstrainedVertexUpdator();
   tBuilder = new ConstrainedTreeBuilder;
   defaultParameters();
   iterations = -1;
   csum = -1000.0;
 }

KinematicConstrainedVertexFitter::~KinematicConstrainedVertexFitter ( )

Definition at line 27 of file KinematicConstrainedVertexFitter.cc.

References finder, tBuilder, updator, and vCons.

                                                                     {
   delete finder;
   delete vCons;
   delete updator;
   delete tBuilder;
 }

Member Function Documentation

void KinematicConstrainedVertexFitter::defaultParameters ( )

private

Definition at line 41 of file KinematicConstrainedVertexFitter.cc.

References theMaxDelta, theMaxReducedChiSq, theMaxStep, and theMinChiSqImprovement.

Referenced by fit(), and KinematicConstrainedVertexFitter().

                                                          {
   theMaxDelta = 0.01;
   theMaxStep = 1000;
   theMaxReducedChiSq = 225.;
   theMinChiSqImprovement = 50.;
 }

RefCountedKinematicTree KinematicConstrainedVertexFitter::fit ( const std::vector< RefCountedKinematicParticle > & part )

inline

Without additional constraint, this will perform a simple vertex fit using LMS with Lagrange multipliers method.

Definition at line 46 of file KinematicConstrainedVertexFitter.h.

Referenced by KineExample::analyze(), fit(), BPHKinematicFit::kinematicTree(), trackingPlots.Iteration::modules(), and ConversionVertexFinder::run().

                                                                                   {
     return fit(part, nullptr, nullptr);
   }

RefCountedKinematicTree KinematicConstrainedVertexFitter::fit	(	const std::vector< RefCountedKinematicParticle > &	part,
		MultiTrackKinematicConstraint *	cs
	)

inline

LMS with Lagrange multipliers fit of vertex constraint and user-specified constraint.

Definition at line 53 of file KinematicConstrainedVertexFitter.h.

References fwrapper::cs, defaultParameters(), fit(), getCSum(), getNit(), and DiDispStaMuonMonitor_cfi::pt.

Referenced by trackingPlots.Iteration::modules().

                                                                                                                      {
     return fit(part, cs, nullptr);
   };

RefCountedKinematicTree KinematicConstrainedVertexFitter::fit	(	const std::vector< RefCountedKinematicParticle > &	part,
		MultiTrackKinematicConstraint *	cs,
		GlobalPoint *	pt
	)

LMS with Lagrange multipliers fit of vertex constraint, user-specified constraint and user-specified starting point.

Definition at line 48 of file KinematicConstrainedVertexFitter.cc.

References ConstrainedTreeBuilder::buildTree(), csum, dumpMFGeometry_cfg::delta, Vispa.Plugins.EdmBrowser.EdmDataAccessor::eq(), finder, LinearizationPointFinder::getLinearizationPoint(), mps_fire::i, input, KinematicState::isValid(), iterations, dqmiolumiharvest::j, KinematicState::kinematicParameters(), KinematicState::kinematicParametersError(), LogDebug, KinematicParametersError::matrix(), allConversions_cfi::maxDelta, ecalTrigSettings_cff::particles, DiDispStaMuonMonitor_cfi::pt, InputSort::sort(), tBuilder, theMaxDelta, theMaxReducedChiSq, theMaxStep, theMinChiSqImprovement, KinematicConstrainedVertexUpdator::update(), updator, KinematicParameters::vector(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by trackingPlots.Iteration::modules().

                                                                                {
   if (part.size() < 2)
     throw VertexException("KinematicConstrainedVertexFitter::input states are less than 2");
 
   //sorting out the input particles
   InputSort iSort;
   std::pair<std::vector<RefCountedKinematicParticle>, std::vector<FreeTrajectoryState> > input = iSort.sort(part);
   const std::vector<RefCountedKinematicParticle> &particles = input.first;
   const std::vector<FreeTrajectoryState> &fStates = input.second;
 
   // linearization point
   // (only compute it using the linearization point finder if no point was passed to the fit function):
   GlobalPoint linPoint;
   if (pt != nullptr) {
     linPoint = *pt;
   } else {
     linPoint = finder->getLinearizationPoint(fStates);
   }
 
   //initial parameters:
   int vSize = particles.size();
   AlgebraicVector inPar(3 + 7 * vSize, 0);
 
   //final parameters
   AlgebraicVector finPar(3 + 7 * vSize, 0);
 
   //initial covariance
   AlgebraicMatrix inCov(3 + 7 * vSize, 3 + 7 * vSize, 0);
 
   //making initial vector of parameters and initial particle-related covariance
   int nSt = 0;
   std::vector<KinematicState> inStates;
   for (std::vector<RefCountedKinematicParticle>::const_iterator i = particles.begin(); i != particles.end(); i++) {
     KinematicState state = (*i)->stateAtPoint(linPoint);
     if (!state.isValid()) {
       LogDebug("KinematicConstrainedVertexFitter") << "State is invalid at point: " << linPoint << std::endl;
       return ReferenceCountingPointer<KinematicTree>(new KinematicTree());
     }
     AlgebraicVector prPar = asHepVector<7>(state.kinematicParameters().vector());
     for (int j = 1; j < 8; j++) {
       inPar(3 + 7 * nSt + j) = prPar(j);
     }
     AlgebraicSymMatrix l_cov = asHepMatrix<7>(state.kinematicParametersError().matrix());
     inCov.sub(4 + 7 * nSt, 4 + 7 * nSt, l_cov);
     inStates.push_back(state);
     ++nSt;
   }
 
   //initial vertex error matrix components (huge error method)
   //and vertex related initial vector components
   double in_er = 100.;
   inCov(1, 1) = in_er;
   inCov(2, 2) = in_er;
   inCov(3, 3) = in_er;
 
   inPar(1) = linPoint.x();
   inPar(2) = linPoint.y();
   inPar(3) = linPoint.z();
 
   //constraint equations value and number of iterations
   double eq;
   int nit = 0;
   iterations = 0;
   csum = 0.0;
 
   std::vector<KinematicState> lStates = inStates;
   GlobalPoint lPoint = linPoint;
   RefCountedKinematicVertex rVtx;
   AlgebraicMatrix refCCov;
 
   double chisq = 1e6;
   bool convergence = false;
   //iterarions over the updator: each time updated parameters
   //are taken as new linearization point
   do {
     eq = 0.;
     std::pair<std::pair<std::vector<KinematicState>, AlgebraicMatrix>, RefCountedKinematicVertex> lRes =
         updator->update(inPar, inCov, lStates, lPoint, cs);
 
     const std::vector<KinematicState> &newStates = lRes.first.first;
 
     if (particles.size() != newStates.size()) {
       LogDebug("KinematicConstrainedVertexFitter") << "updator failure\n";
       return ReferenceCountingPointer<KinematicTree>(new KinematicTree());
     }
 
     rVtx = lRes.second;
 
     double newchisq = rVtx->chiSquared();
     if (nit > 2 && newchisq > theMaxReducedChiSq * rVtx->degreesOfFreedom() &&
         (newchisq - chisq) > (-theMinChiSqImprovement)) {
       LogDebug("KinematicConstrainedVertexFitter") << "bad chisq and insufficient improvement, bailing\n";
       return ReferenceCountingPointer<KinematicTree>(new KinematicTree());
     }
     chisq = newchisq;
 
     const GlobalPoint &newPoint = rVtx->position();
 
     double maxDelta = 0.0;
 
     double deltapos[3];
     deltapos[0] = newPoint.x() - lPoint.x();
     deltapos[1] = newPoint.y() - lPoint.y();
     deltapos[2] = newPoint.z() - lPoint.z();
     for (int i = 0; i < 3; ++i) {
       double delta = deltapos[i] * deltapos[i] / rVtx->error().matrix()(i, i);
       if (delta > maxDelta)
         maxDelta = delta;
     }
 
     for (std::vector<KinematicState>::const_iterator itold = lStates.begin(), itnew = newStates.begin();
          itnew != newStates.end();
          ++itold, ++itnew) {
       for (int i = 0; i < 7; ++i) {
         double deltapar = itnew->kinematicParameters()(i) - itold->kinematicParameters()(i);
         double delta = deltapar * deltapar / itnew->kinematicParametersError().matrix()(i, i);
         if (delta > maxDelta)
           maxDelta = delta;
       }
     }
 
     lStates = newStates;
     lPoint = newPoint;
 
     refCCov = lRes.first.second;
     nit++;
     convergence = maxDelta < theMaxDelta || (nit == theMaxStep && maxDelta < 4.0 * theMaxDelta);
 
   } while (nit < theMaxStep && !convergence);
 
   if (!convergence) {
     return ReferenceCountingPointer<KinematicTree>(new KinematicTree());
   }
 
   // std::cout << "old full cov matrix" << std::endl;
   // std::cout << refCCov << std::endl;
 
   // cout<<"number of relinearizations "<<nit<<endl;
   // cout<<"value obtained: "<<eq<<endl;
   iterations = nit;
   csum = eq;
 
   return tBuilder->buildTree(particles, lStates, rVtx, refCCov);
 }

float KinematicConstrainedVertexFitter::getCSum ( ) const

Definition at line 197 of file KinematicConstrainedVertexFitter.cc.

References csum.

Referenced by fit().

197 { return csum; }

KinematicConstrainedVertexFitter::csum

float csum

Definition: KinematicConstrainedVertexFitter.h:81

int KinematicConstrainedVertexFitter::getNit ( ) const

Definition at line 195 of file KinematicConstrainedVertexFitter.cc.

References iterations.

Referenced by fit(), and ConversionVertexFinder::run().

195 { return iterations; }

KinematicConstrainedVertexFitter::iterations

int iterations

Definition: KinematicConstrainedVertexFitter.h:80

void KinematicConstrainedVertexFitter::setParameters ( const edm::ParameterSet & pSet )

Configuration through PSet: number of iterations(maxDistance) and stopping condition (maxNbrOfIterations)

Definition at line 34 of file KinematicConstrainedVertexFitter.cc.

References edm::ParameterSet::getParameter(), theMaxDelta, theMaxReducedChiSq, theMaxStep, and theMinChiSqImprovement.

                                                                                 {
   theMaxDelta = pSet.getParameter<double>("maxDelta");
   theMaxStep = pSet.getParameter<int>("maxNbrOfIterations");
   theMaxReducedChiSq = pSet.getParameter<double>("maxReducedChiSq");
   theMinChiSqImprovement = pSet.getParameter<double>("minChiSqImprovement");
 }