#include <VtVertex.hh>

Inheritance diagram for VERTEX::Vertex:

[legend]

Collaboration diagram for VERTEX::Vertex:

[legend]

Public Member Functions
— Constructors —
	Vertex ()

	Vertex (const Vertex &rhs)

	Vertex (const Track_v &v)

	Vertex (Track &t1, Track &t2)

	~Vertex ()

— VertexIf access functions —
float	vx () const
	$x$ of vertex More...

float	vy () const
	$y$ of vertex More...

float	vz () const
	$z$ of vertex More...

float	chi2 () const
	$\chi^2$ of vertex fit More...

float	prob () const
	upper tail $\chi^2$ probability More...

unsigned short int	ndf () const
	degrees of freedom of vertex fit More...

unsigned short int	ntracks () const
	no of tracks in vertex $n$ More...

float	vtx_cov_x () const
	$\sigma_x^2$ of vertex More...

float	vtx_cov_y () const
	$\sigma_y^2$ of vertex More...

float	vtx_cov_z () const
	$\sigma_z^2$ of vertex More...

– Other access methods —
double	angle () const

double	dist () const
	available after call to rmsDistAngle() More...

double	vxerr () const
	$\sqrt{\sigma_{vx}^2}$ vertex $x$-error More...

double	vyerr () const
	$\sqrt{\sigma_{vy}^2}$ vertex $y$-error More...

double	vzerr () const
	$\sqrt{\sigma_{vz}^2}$ vertex $z$-error More...

double	pmaxfrac () const
	$\max\|p_i/p_j\|$ max. momentum fraction More...

const MATRIX::VtSymMatrix &	CS () const
	$3\times3$ Vertex covariance matrix More...

const MATRIX::VtSymMatrix &	covn () const
	$(3\cdot n+3)\times(3\cdot n+3)$ general covariance matrix More...

bool	valid () const
	is vertex valid? More...

void	set_invalid ()
	mark vertex as invalid More...

— Mass calculation routines —
@memo $m = \sqrt{\sum_i p_i^2 + m_i^2 - \sum_i \vec{p}_i^2}$
double	mass (const bool use=false) const
	rest-mass 0 or track rest mass More...

double	massCC (const bool use=false) const
	rest-mass 0 or conjugated track rest mass More...

double	mass (double m1) const
	single rest-mass calc., More...

double	mass (double m1, double m2) const
	2-track mass with particle rest-masses More...

double	mass (double m1, double m2, double m3) const
	3-track mass with rest-masses More...

double	mass (double m1, double m2, double m3, double m4) const
	4-track mass with rest-masses More...

double	mass (const Vector_d &m02, bool use=false, bool CC=false) const
	specify $n$ rest masses More...

— Vertexing routines —
const bool	findVertex2D (void)
	unweighted 2-D min. dist. More...

const bool	findVertex3D (void)
	unweighted 3-D min. dist. More...

const bool	findVertexVt ()
	Vt Kalman-filter. More...

const bool	VtEstimateVertex ()
	weighted 2-D min. dist. More...

const bool	VtEstimateVertexMath (double &x, double &y, double &z)
	estimate Vertex without changing Vertex object More...

const bool	VtEstimateVertexMathTA (double &x, double &y, double &z)
	estimate Vertex without changing Vertex object More...

const bool	VtMass ()
	Kalman filter with mass constraints. More...

— Adding tracks —
void	push_back (Track &track)

— Misc functions —
double	rmsDistAngle () const
	calc rms dist and rms angle More...

bool	calc_mother (bool use_kalman=true)
	calc track params of mother particle More...

bool	calc_mother_cov ()
	calc cov. matrix of mother particle More...

const bool	remove_last ()
	remove last track from vertex, do inverse Kalman filter More...

void	remove_worst ()
	remove track with biggest $\chi^2$ contribution, do inverse Kalman filter More...

double	track_chi2 (int i)
	get track $\chi^2$ contribution More...

int	track_worst ()
	get track with biggest $\chi^2$ contribution More...

void	use_momentum (const bool use)
	set for all tracks whether momentum should be used or not More...

void	use_kalman (const bool use)
	use refitted track parameters or not More...

bool	use_kalman () const

double	distance (double x, double y, double z) const
	$\chi^2$ distance to space point, $ndf = 3$ More...

double	distance (Track &t) const
	$\chi^2$ distance to track, $ndf = 2$ More...

double	distance (const Vertex &rhs) const
	$\chi^2$ distance to vertex, $ndf = 3$ More...

— Expert routines —
void	clear ()
	clear all track-vertex relations, makes vertex invalid More...

const bool	VtFilter ()

const bool	VtInverseFilter () const

void	VtSmoothX ()
	smooth vertex position More...

void	VtSmoothQ ()
	smooth momenta More...

const bool	VtRemoveTrack (Relation &track)
	remove track using inverse Kalman filter More...

const bool	VtFit ()
	Kalman filter + smoother. More...

const MATRIX::VtVector &	xv () const

const MATRIX::VtVector &	xvs () const
	vector of vertex-positions More...

std::ostream &	print (std::ostream &os) const
	called by operator<<() More...

double	chi2n () const
	$\chi^2$ after last filter step More...

double	chi2l () const
	$\chi^2$ contribution of last filter step More...

virtual void	remove (Relation *r)
	unbook relation More...

virtual const iterator	erase (const iterator &pos)
	erase a relation without refit More...

const MATRIX::VtSymMatrix &	bigcov ()
	compute general $(3\cdot n+3) \times (3\cdot n+3)$ covariance matrix More...

— Mass Constraint routines —
void	MassConstr (double m)
	add a mass constraint for all tracks, $m$ = mass of mother particle More...

MassC &	addMassConstr (double m=0.)
	add a mass constraint More...

void	clearMassConstr ()
	delete the mass constraints More...

unsigned int	nMassConstr () const
	number of mass constraints More...

Public Member Functions inherited from VERTEX::Track
	Track ()

	Track (const Track *const track)

	Track (const Track &rhs)

	Track (const MATRIX::VtVector &v, const MATRIX::CMatrix &c)

virtual	~Track ()

void	set (double x, double y, double z, double tx, double ty, double p, const MATRIX::CMatrix &c)

Track &	operator= (const Track &rhs)

float	x () const
	x at z=0 (in Vt: p(3,..)) More...

float	y () const
	y at z=0 (in Vt: p(4,..)) More...

float	x (float z) const
	Track $x$ position at $z$. More...

float	y (float z) const
	Track $y$ position at $z$. More...

float	z () const
	z = 0 in Vt More...

float	tx () const
	slope (in Vt: p(1,..)) More...

float	ty () const
	slope (in Vt: p(2,..)) More...

float	p () const
	momentum (in Vt: p(5,..)) More...

float	pt () const
	transv. momentum $\sqrt{p_x^2 + p_y^2}$ More...

float	pz () const
	$p \cdot e_z$ More...

float	chi2 () const
	dummy function: always return 0 More...

unsigned short int	ndf () const
	dummy function: always return 0 More...

float	phi () const
	azimuthal angle $\phi$ More...

float	theta () const
	polar angle $\theta = \cos^{-1}(e_z)$ More...

float	eta () const
	rapidity $-\log\tan(\theta/2.)$ More...

int	charge () const
	charge More...

float	energy (double rm=0.) const
	Energy with given rest-mass (in GeV) $E = \sqrt{p^2 + m^2}$. More...

float	xf (double rm=0.) const
	Feynman Variable $x_F = \frac{E+p_z}{(E+p_z)_{max}}$. More...

float	rap (double rm=0.) const
	Rapidity $y = \frac{1}{2}\ln(\frac{E+p_z}{E-p_z})$. More...

MATRIX::VtVector	evec () const
	$\vec{v} = (e_x,e_y,e_z)$ unit vector along refitted track More...

MATRIX::VtVector	tvec () const
	$\vec{v} = (t_x,t_y,1.)$ refitted slope vector More...

MATRIX::VtVector	pvec () const
	$\vec{v} = (p_x,p_y,p_z)$ refitted mom. vector More...

float	cov_x (double dz=0.) const
	get $\sigma_x^2\|_{z+dz}$ More...

float	cov_y (double dz=0.) const
	get $\sigma_y^2\|_{z+dz}$ More...

float	cov_tx () const
	get $\sigma_{t_x}^2$ More...

float	cov_ty () const
	get $\sigma_{t_y}^2$ More...

float	cov_p () const
	get $\sigma_p^2$ More...

bool	isValid () const
	dummy function: always return true More...

void	valid ()
	dummy function: do nothing More...

void	invalid ()
	dummy function: do nothing More...

bool	propagate (const double zz)
	propagate track to $z$ More...

double	px () const
	$p \cdot e_x$ More...

double	py () const
	$p \cdot e_y$ More...

double	ex () const
	$t_x \cdot e_z$ More...

double	ey () const
	$t_y \cdot e_z$ More...

double	ez () const
	$1/\sqrt{1+t_x^2+t_y^2}$ More...

double	xerr (double dz=0) const
	$\sqrt{\sigma_x^2(z+dz)}$ More...

double	yerr (double dz=0) const
	$\sqrt{\sigma_y^2(z+dz)}$ More...

double	txerr () const
	$\sqrt{\sigma^2_{tx}}$ More...

double	tyerr () const
	$\sqrt{\sigma^2_{ty}}$ More...

double	perr () const
	$\sqrt{\sigma^2_p}$ More...

const MATRIX::CMatrix &	V () const
	covariance matrix More...

const MATRIX::CMatrix &	G () const
	inverse $5\times5$ covariance matrix More...

const MATRIX::CMatrix &	GM () const
	inverse $4\times4$ covariance matrix, without momentum More...

std::ostream &	print (std::ostream &os) const
	called by cout More...

double	rm () const
	rest-mass More...

void	rm (const double mass)
	set rest-mass (needed for mass constrained fits) More...

double	rmCC () const
	conjugated rest-mass More...

void	rmCC (const double mass)
	set conjugated rest-mass More...

void	delete_mom ()

Public Member Functions inherited from VERTEX::RelationList
	RelationList ()

virtual	~RelationList ()

void	push__back (Relation *rel)
	add a relation More...

Relation *	back () const
	return last relation More...

unsigned int	size () const
	return no of relations More...

iterator	begin ()

iterator	end ()

const_iterator	begin () const

const_iterator	end () const

reverse_iterator	rbegin ()

reverse_iterator	rend ()

const_reverse_iterator	rbegin () const

const_reverse_iterator	rend () const

Relation *const	operator[] (unsigned int i) const
	direct access More...

void	print (std::ostream &os) const

void	unbook (Relation *const rel)
	get rid of relation pointer More...

const RelationList &	operator= (const RelationList &rhs)

bool	operator== (const RelationList &rhs) const

— Operators —
bool	v_use_kalman

bool	v_valid

bool	v_Mvalid

bool	v_mother

bool	v_angdist

double	v_angle

double	v_dist

double	v_bk13

double	v_chi2

MATRIX::VtVector	kal_xv

MATRIX::VtVector	kal_xvs

MATRIX::VtSymMatrix	v_CS

MATRIX::VtSymMatrix	v_CINV

MATRIX::VtSymMatrix *	v_covn

Track_v	tracks

MassC_v	xmass

const Vertex &	operator= (const Vertex &rhs)

const bool	operator== (const Vertex &rhs) const

const double	operator- (const Vertex &rhs) const

void	add_track (Track &t)

const Track *	get_track (int i) const

const unsigned int	bigdim () const

void	invalid ()

	ClassDef (Vertex, 1)

Additional Inherited Members
Public Types inherited from VERTEX::RelationList
typedef ConstRelationIterator	const_iterator

typedef ConstReverseRelationIterator	const_reverse_iterator

typedef RelationIterator	iterator

typedef ReverseRelationIterator	reverse_iterator

Protected Member Functions inherited from VERTEX::Track
	ClassDef (Track, 1)

Protected Attributes inherited from VERTEX::Track
MATRIX::VtVector	t_p

MATRIX::CMatrix	t_V

MATRIX::CMatrix	t_G

MATRIX::CMatrix	t_GM

Detailed Description

Vertex class

Constructor & Destructor Documentation

◆ Vertex() [1/4]

Vertex::Vertex ( )

                  :
     Track       (),
     v_use_kalman(false),
     v_valid     (false),
     v_Mvalid    (false),
     v_mother    (false),
     v_angdist   (false),
     v_angle     (0.),
     v_dist      (0.),
     v_bk13      (0.),
     v_chi2      (0.),
     kal_xv      (VtVector(3)),
     kal_xvs     (VtVector(3)),
     v_CS        (VtSymMatrix(3,1.e5)),
     v_CINV      (CINV0),
     v_covn      (0),
     tracks      (),
     xmass       ()
   {}

◆ Vertex() [2/4]

Vertex::Vertex ( const Vertex & rhs )

                                  :
     Track       (rhs),
     v_use_kalman(rhs.v_use_kalman),
     v_valid     (rhs.v_valid),
     v_Mvalid    (rhs.v_Mvalid),
     v_mother    (rhs.v_mother),
     v_angdist   (rhs.v_angdist),
     v_angle     (rhs.v_angle),
     v_dist      (rhs.v_dist),
     v_bk13      (rhs.v_bk13),
     v_chi2      (rhs.v_chi2),
     kal_xv      (rhs.kal_xv),
     kal_xvs     (rhs.kal_xvs),
     v_CS        (rhs.v_CS),
     v_CINV      (v_CS.dsinv()),
     v_covn      (0)
     //    tracks  (rhs.tracks),
     //    xmass   (rhs.xmass)
   {
     // do not copy relations
   }

◆ Vertex() [3/4]

VERTEX::Vertex::Vertex ( const Track_v & v )

◆ Vertex() [4/4]

Vertex::Vertex	(	Track &	t1,
		Track &	t2
	)

                                     :
     Track       (),
     v_use_kalman(false),
     v_valid     (false),
     v_Mvalid    (false),
     v_mother    (false),
     v_angdist   (false),
     v_angle     (0.),
     v_dist      (0.),
     v_bk13      (0.),
     v_chi2      (0.),
     kal_xv      (VtVector(3)),
     kal_xvs     (VtVector(3)),
     v_CS        (VtSymMatrix(3,1.e5)),
     v_CINV      (CINV0),
     v_covn      (0),
     tracks      (),
     xmass       ()
  {
     createRelation(t1, *this);
     createRelation(t2, *this);
   }

◆ ~Vertex()

Vertex::~Vertex ( )

                   {
     // delete the mass constraints
     clearMassConstr();
     // delete tracks owned by vertex
     const Track_cit end = tracks.end();
     for(Track_cit it = tracks.begin(); it != end; ++it) {
       delete *it;
     }
     // delete big covariance matrix
     if(v_covn != 0) delete v_covn;
   }

Member Function Documentation

◆ add_track()

void Vertex::add_track ( Track & t )

                                  {
     tracks.push_back(&t);   // add track to the tracks list in a way that vertex become the owner
     push_back(t);           // set the relation
   }

◆ addMassConstr()

MassC & Vertex::addMassConstr ( double m = 0. )

add a mass constraint

                                        {
     xmass.push_back(new MassC(m));
     return *xmass.back();
   }

◆ angle()

double VERTEX::Vertex::angle ( ) const

$\sqrt{n\cdot\sum_it_{x,i}^2 - (\sum_i t_{x,i})^2 + n\cdot\sum_i t_{y,i}^2 - (\sum_i t_{y,i})^2} / n$

◆ bigcov()

const VtSymMatrix & Vertex::bigcov ( )

compute general $(3\cdot n+3) \times (3\cdot n+3)$ covariance matrix

                                     {
  
     // build up big cov. matrix (in Vt: covn)
     // covn consists of vertex + track parameter cov. matrices:
     //  vtx  tr1  tr2  tr3  ...
     //
     //  CS | ES | ES | ES | ...   CC = cor(q_k,q_l) = F_k*C_n*F_l^T
     //  ---+----+----+----+          = -ES_k * F_l
     //  ES | DS | CC | CC | ...
     //  ---+----+----+----+ 
     //  ES | CC | DS | CC | ...
     //  ---+----+----+----+
     //  ES | CC | CC | DS | ...
     //  ...
  
     if(v_covn == 0)
       v_covn = new VtSymMatrix(bigdim());
  
     VtSymMatrix& covn = *v_covn;
  
     //    VtSymMatrix covn(bigdim());
  
     covn.place_at(v_CS,0,0);
  
     unsigned int ipos = 3;
     const iterator last = end();
     for(iterator it = begin(); it != last; ++it) {
       const Kalman& kalman = it->kalman;
  
       covn.place_at(kalman.ES(), ipos, 0);
       covn.place_at(kalman.DS(), ipos, ipos);
  
       // fill cor(q_k^n,q_l^n) = F_k*C_n*F_l^T = -E_k^n * F_l^T
       unsigned int jpos = 6;
       iterator tmp=it;
       for(iterator it2=++tmp; it2 != last; ++it2) {
     const Kalman& kalman2 = it2->kalman;
     
     covn.place_at(-kalman.ES()*kalman2.F().T(), ipos,jpos);
     jpos+=3;
       }
  
       ipos+=3;
     }
  
     return covn;
   }

◆ bigdim()

const unsigned int VERTEX::Vertex::bigdim ( ) const

private

◆ calc_mother()

bool Vertex::calc_mother ( bool use_kalman = true )

calc track params of mother particle

                                           {
     //    if(v_valid == false) return v_mother = false;
  
     // add up track momenta
     t_p[5] = 0;
 #if defined USE_ROOT
     SVector<double,3> psum;
 #else
     VtVector psum(0,0,0);
 #endif
  
     const iterator last = end();
     if(use_kalman == true) {
       for(iterator it = begin(); it != last; ++it) {
     const Kalman& t = it->kalman;
     t_p[5] += t.p();
     psum   += t.pvec();
       }
  
       if(fabs(psum[2])>1e-20) {
     t_p[3] = psum[0]/psum[2]; // tx = pxsum/pzsum
     t_p[4] = psum[1]/psum[2]; // ty = pysum/pzsum
       } else {
     return v_mother = false;
       }
  
       // x,y at z = 0
       t_p[0] = vx() - vz() * tx();
       t_p[1] = vy() - vz() * ty();
       t_p[2] = 0.;
  
     } else {
       double x, y, z;
       VtEstimateVertexMath(x,y,z);
       for(iterator it = begin(); it != last; ++it) {
     const Track& t = it->track;
     t_p[5] += t.p();
     psum   += t.pvec();
       }
  
       if(fabs(psum[2])>1e-20) {
     t_p[3] = psum[0]/psum[2]; // tx = pxsum/pzsum
     t_p[4] = psum[1]/psum[2]; // ty = pysum/pzsum
       } else {
     return v_mother = false;
       }
  
       // x,y at z = 0
       t_p[0] = x - z * tx();
       t_p[1] = y - z * ty();
       t_p[2] = 0.;
     }
  
     //    t_p[5] = sqrt(sqr(pxsum) + sqr(pysum) + sqr(pzsum));
  
     // what is the difference betwenn sum of track.p() and
     // calculation sqrt(sqr(pxsum)+sqr(pysum)+sqr(pzsum))?????
     //      cout << " reco-p: " << sqrt(sqr(pxsum) + sqr(pysum) + sqr(pzsum))
     //       << " orig reco p: " << p() << endl;
  
     return v_mother = true;
   }

◆ calc_mother_cov()

bool Vertex::calc_mother_cov ( )

calc cov. matrix of mother particle

                                {
     // no track parameters computed
     if(v_mother == false || v_valid == false) return false;
     // get big covariance matrix
     const VtSymMatrix& covn = (v_Mvalid==false)? bigcov() : *v_covn;
  
     // fill jacobian
     VtMatrix    jac(bigdim(),5);
  
     // unit vectors
     double epx = ex();
     double epy = ey();
     double epz = ez();
  
     // vertex part:
     jac(0,0) = jac(1,1) = 1.;
     jac(2,0) = -tx();
     jac(2,1) = -ty();
  
     double fact3 = 1./pz();
     //    double fact3 = 1./(p() * epz);
  
     unsigned int i = 2;
     const iterator last = end();
     for(iterator it = begin(); it != last; ++it) {
       const Kalman& t = it->kalman;
       
       // unit vectors along track:
       double ex = t.ex();
       double ey = t.ey();
       double ez = t.ez();
  
       double fact2 = t.pz();
       double fact  = fact2 * fact3;
  
       ++i;
       jac(i,2) = fact * (1. - sqr(ex) + tx() * ex * ez);
       jac(i,3) = fact * ex * (ty()*ez - ey);
       jac(i,0) = -vz() * jac(i,2);
       jac(i,1) = -vz() * jac(i,3);
  
       jac(i,4) = fact2 * ( epx*(1.-sqr(ex)) - epy*ex*ey - epz*ex*ez );
       
       ++i;
       jac(i,2) = fact*ey*(tx()*ez-ex);
       jac(i,3) = fact*(1. - sqr(ey) + ty()*ey*ez);
       jac(i,0) = -vz() * jac(i,2);
       jac(i,1) = -vz() * jac(i,3);
  
       jac(i,4) = fact2 * ( epy*(1.-sqr(ey)) - epx*ex*ey - epz*ey*ez );
  
       ++i;
       jac(i,2) = fact3 * (ex - tx()*ez);
       jac(i,3) = fact3 * (ey - ty()*ez);
       jac(i,0) = -vz() * jac(i,2);
       jac(i,1) = -vz() * jac(i,3);
       
       jac(i,4) = epx*ex + epy*ey + epz*ez;
     }
  
 #ifdef VtDEBUG
     cout << " jac: " << jac << endl;
 #endif
  
     // compute Track covariance matrix:
     t_V = covn.product(jac);
     // compute inverse
     t_G = t_V;
     t_G.VtDsinv();
     t_GM.copy(t_V.dsinv(4)); // inverse without momentum
  
     return true;
   }

◆ chi2()

float Vertex::chi2 ( ) const

$\chi^2$ of vertex fit

236 { return v_chi2; }

◆ chi2l()

double Vertex::chi2l ( ) const

$\chi^2$ contribution of last filter step

                              { 
     const_reverse_iterator it1;
     const_reverse_iterator it2 = it1 = rbegin();
     ++it2;
     return it1->kalman.chi2() - it2->kalman.chi2();
   }

◆ chi2n()

double Vertex::chi2n ( ) const

$\chi^2$ after last filter step

                              { 
     return back()->kalman.chi2();
   }

◆ ClassDef()

VERTEX::Vertex::ClassDef	(	Vertex	,
		1
	)

private

◆ clear()

void Vertex::clear ( void )

virtual

clear all track-vertex relations, makes vertex invalid

Reimplemented from VERTEX::RelationList.

                      {
     RelationList::clear();
     invalid();
     v_mother = v_use_kalman = false;
     // delete the mass constraints
     clearMassConstr();
     // delete tracks owned by vertex
     const Track_cit end = tracks.end();
     for(Track_cit it = tracks.begin(); it != end; ++it) {
       delete *it;
     }
     tracks.clear();
     // delete big covariance matrix
     if(v_covn != 0) { 
       delete v_covn; 
       v_covn = 0;
     }
   }

◆ clearMassConstr()

void Vertex::clearMassConstr ( )

delete the mass constraints

                                {
     // delete mass constraints owned by vertex
     const MassC_cit end = xmass.end();
     for(MassC_cit it = xmass.begin(); it != end; ++it) {
       delete *it;
     }
     xmass.clear();
     return;
   }

◆ covn()

const MATRIX::VtSymMatrix& VERTEX::Vertex::covn ( ) const

$(3\cdot n+3)\times(3\cdot n+3)$ general covariance matrix

◆ CS()

const MATRIX::VtSymMatrix& VERTEX::Vertex::CS ( ) const

$3\times3$ Vertex covariance matrix

◆ dist()

double VERTEX::Vertex::dist ( ) const

available after call to rmsDistAngle()

◆ distance() [1/3]

double Vertex::distance ( const Vertex & rhs ) const

$\chi^2$ distance to vertex, $ndf = 3$

                                                  {
     const VtSymMatrix& cov1 = back()->kalman.CINV();
     const VtSymMatrix& cov2 = rhs.back()->kalman.CINV();
     VtSymMatrix Cmat = cov1+cov2;
  
     if(Cmat.VtDsinv() == false) return 1000.;
  
     const VtVector x1  = cov1 * kal_xvs;
     const VtVector x2  = cov2 * rhs.kal_xvs;
     const VtVector xav = Cmat * (x1+x2);
     const VtVector dx1 = xav - kal_xvs;
     const VtVector dx2 = xav - rhs.kal_xvs;
  
     // get a chi2 value
     double sum1 = cov1.product(dx1);
     double sum2 = cov2.product(dx2);
     return (sum1 + sum2);
   }

◆ distance() [2/3]

double Vertex::distance	(	double	x,
		double	y,
		double	z
	)		const

$\chi^2$ distance to space point, $ndf = 3$

                                                                {
     VtVector dx(xx-vx(),yy-vy(),zz-vz());
     return back()->kalman.CINV().product(dx);
   }

◆ distance() [3/3]

double Vertex::distance ( Track & t ) const

$\chi^2$ distance to track, $ndf = 2$

                                         {
     Vertex tmp;
     const const_iterator last = end();
     for(RelationList::const_iterator it = begin(); it != last; ++it) {
         Track& tp = it->track;
     tmp.push_back(tp);
     }
     tmp.push_back(t);
     if (tmp.findVertexVt() == false) return -10.;
     return tmp.track_chi2(tmp.ntracks()-1);
   }

◆ erase()

const iterator Vertex::erase ( const iterator & pos )

virtual

erase a relation without refit

Reimplemented from VERTEX::RelationList.

                                                   {
     invalid();
     return RelationList::erase(pos);
   }

◆ findVertex2D()

const bool Vertex::findVertex2D ( void )

unweighted 2-D min. dist.

                                       {  
  
     unsigned int v_nrseg = size();
     if(v_nrseg < 2) return false;
     
     double sumtx=0, sumty=0, sumrx=0, sumry=0, sumrt=0, sumtt=0;
  
     double tx, ty, rx, ry;
     const RelationList::iterator last = end();
     for(RelationList::iterator pos = begin(); pos != last; ++pos) {
       const Track& track = pos->track;
  
       // slopes and position at z=0
       tx = track.tx();  rx = track.x() - tx*track.z();
       ty = track.ty();  ry = track.y() - ty*track.z();
       
       // update cumulants
       sumtx += tx;   sumrx += rx;
       sumty += ty;   sumry += ry;
       sumrt += rx*tx   + ry*ty;
       sumtt += sqr(tx) + sqr(ty);
     }
  
 #ifdef VtDEBUG
     cout << " sumtx: " << sumtx << " sumty: " << sumty
      << " sumrx: " << sumrx << " sumry: " << sumry
      << " sumrt: " << sumrt << " sumtt: " << sumtt
      << endl;
 #endif
  
     // get vertex z-position
  
     double z = -(v_nrseg*sumrt - sumrx*sumtx - sumry*sumty) / 
                 (v_nrseg*sumtt - sqr<double>(sumtx) - sqr<double>(sumty));
   
     // and the transverse coordinates
     double x = (sumrx + z*sumtx)/static_cast<double>(v_nrseg);
     double y = (sumry + z*sumty)/static_cast<double>(v_nrseg);
  
     // fill vertex parameters
     kal_xvs[0]  = x;         // x-vertex
     kal_xvs[1]  = y;         // y-vertex
     kal_xvs[2]  = z;         // z-vertex
  
     v_valid = true;
     return true;
   }

◆ findVertex3D()

const bool Vertex::findVertex3D ( void )

unweighted 3-D min. dist.

                                       {
  
     if(size() < 2) return false;
  
     double c[3], v[3], tt[9];
     for(int i=0; i<9; i++) tt[i] = 0.;
     for(int i=0; i<3; i++) c[i]  = v[i] = 0.;
  
     double t[3], r[3], dt;
     const RelationList::iterator last = end();
     for(RelationList::iterator pos = begin(); pos != last; ++pos) {
       const Track& track = pos->track;
  
       t[0] = track.tx();  r[0] = track.x();
       t[1] = track.ty();  r[1] = track.y();
       t[2] = 1;           r[2] = track.z();
       
       double tsq = sqr(t[0]) + sqr(t[1]) + sqr(t[2]); 
       
       for(int i=0; i<3; i++) {
     for(int j=0; j<3; j++) {
       dt = t[i]*t[j]/tsq; if(i==j) dt -= 1; 
       c[i]      += dt*r[j];
       tt[3*i+j] += dt;
     }
       }
     }
     
     int ifail = 0, ndim = 3;
     double work[3];
     Dinv1(&ndim, tt, &ndim, work, &ifail);
     
     if(ifail) {
       cout << "findVertex: matrix inversion failed" << endl;
       return false;
     }
     
     // get the vertex position
     for(int i=0; i<3; i++) {
       for(int j=0; j<3; j++) v[i] += tt[3*i+j]*c[j];
     }
     printf("findVertex3D v: %f %f %f\n",v[0],v[1],v[2]);
  
     // fill vertex parameters
     kal_xvs[0]  = v[0];      // x-vertex
     kal_xvs[1]  = v[1];      // y-vertex
     kal_xvs[2]  = v[2];      // z-vertex
  
     // require that vertex-z is smaller than smallest hit-z of rsegs
     // to avoid vertexing the downstream end of rseg's
     //    for(RelationList::iterator it = begin(); it != last; ++it) {
     //      if(vz() > it->track.z() ) return false;
     //    }
  
     v_valid = true;
     return true;
   }

◆ findVertexVt()

const bool Vertex::findVertexVt ( )

Vt Kalman-filter.

                                   {  
     v_valid = v_Mvalid = false;
     if(ntracks() < 2) { return false; }
     
     // 1st step: find approximate Vertex from average of all given tracks.
     // Initialize Kalman-filter
     VtEstimateVertex();
     // kal_xvs[0] = vest[0];
     // kal_xvs[1] = vest[1];
     // kal_xvs[2] = vest[2];
     
     // 2nd step: initialize kalman filter
     v_CS   = VtSymMatrix(3,1.e5);
     v_CINV = CINV0;
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
       tit->kalman.init();
     }
     
 #ifdef VtDEBUG
     cout << "-------------------------------" << endl;
     cout << "Estimated Vertex: (" << vx() << "," << vy() << "," << vz() << ")" << endl;
     cout << "Number of Tracks: " << size();
     //copy(begin(), end(), ostream_iterator<Relation>(cout," "));
     cout << endl;
 #endif
  
     double prChi2 = 0.;
     const Kalman& trackn = back()->kalman; // last element
  
     // 3rd step: iterate the kalman filter
     unsigned int iter = 0;
     do {
  
       if(VtFit() == false) return false;
  
       double chi2 = trackn.chi2();
  
 #ifdef VtDEBUG
       cout << "----> chi2n: " << chi2 << " prChi2: " << prChi2 << " fabs(chi2-prChi2): " 
        << fabs(chi2 - prChi2) << endl;
 #endif
  
       if(fabs(chi2 - prChi2) > 1.) {
     prChi2 = chi2;
       } else
     if(iter > 0)
       break;
     } while(++iter < 20);
  
 #ifdef VtDEBUG
     cout << "-------------------------------" << endl;
     cout << "No of iterations: " << iter << endl;
     cout << *this;
 #endif
  
     if(iter<20) { //  && chi2()<10.
       v_valid = true;
     } else {
       //      cout << "-----> VERTEX FIT FAILED!! <--------" << endl;
     }
     return v_valid;
   }

◆ get_track()

const Track * Vertex::get_track ( int i ) const

266 {return tracks.at(i);}

◆ invalid()

void Vertex::invalid ( )

private

                        {
     v_valid = v_Mvalid = false;
  
     if(v_covn != 0) {
       delete v_covn;
       v_covn = 0;
     }
  
     return;
   }

◆ mass() [1/6]

double Vertex::mass ( const bool use = false ) const

rest-mass 0 or track rest mass

                                           {
     return mass(vector<double>(ntracks(),0), use, false);
   }

◆ mass() [2/6]

double Vertex::mass	(	const Vector_d &	m02,
		bool	use = `false`,
		bool	CC = `false`
	)		const

specify $n$ rest masses

                                                                   {
  
     double esum  = 0.;
 #if defined USE_ROOT
     SVector<double,3> psum;
 #else
     VtVector psum(0,0,0);
 #endif
  
     unsigned int i=0;
     const const_iterator last = end();
     for(RelationList::const_iterator it = begin(); it != last; ++it) {
       const Track& t = it->track;
       double p = 0.;
       if(v_use_kalman == false && use == false) {
     p     = t.p();
     psum += t.pvec();
       } else {
     const Kalman& k = it->kalman;
     p     = k.p();
     psum += k.pvec();
       }
  
       // in case no rest-mass given, use permanent rest-mass
       double mm = m02[i++];
  
       // handle track rest masses
       if(mm == 0.) {
     double rm = t.rm();
     if(CC == true) rm = t.rmCC(); // use conjugated mass
     mm = rm * rm;
       }
       esum  += sqrt(p*p + mm);
     }
     //   cout << "unllllock mutex" << endl;
     //    pthread_mutex_unlock(&VTmutex);
  
 #if defined USE_ROOT
     double mass2 = esum*esum - mag2(psum);
 #else
     double mass2 = esum*esum - psum(0)*psum(0) - psum(1)*psum(1) -
       psum(2)*psum(2);
 #endif
     if(mass2 < 0) return 0.;
     else return sqrt(mass2);
   }

◆ mass() [3/6]

double Vertex::mass ( double m1 ) const

single rest-mass calc.,

                                      {
     return mass(vector<double>(ntracks(),m1*m1));
   }

◆ mass() [4/6]

double Vertex::mass	(	double	m1,
		double	m2
	)		const

2-track mass with particle rest-masses

                                                 {
     if(ntracks() != 2) return -1;
     vector<double> tmp(2);
     tmp[0] = m1*m1, tmp[1] = m2*m2;
     return mass(tmp);
   }

◆ mass() [5/6]

double Vertex::mass	(	double	m1,
		double	m2,
		double	m3
	)		const

3-track mass with rest-masses

                                                            {
     if(ntracks() != 3) return -1;
     vector<double> tmp(3);
     tmp[0] = m1*m1, tmp[1] = m2*m2, tmp[2] = m3*m3;
     return mass(tmp);
   }

◆ mass() [6/6]

double Vertex::mass	(	double	m1,
		double	m2,
		double	m3,
		double	m4
	)		const

4-track mass with rest-masses

                                                                       {
     if(ntracks() != 4) return -1;
     vector<double> tmp(4);
     tmp[0] = m1*m1, tmp[1] = m2*m2, tmp[2] = m3*m3, tmp[4] = m4*m4;
     return mass(tmp);
   }

◆ massCC()

double Vertex::massCC ( const bool use = false ) const

rest-mass 0 or conjugated track rest mass

                                             {
     return mass(vector<double>(ntracks(),0), use, true);
   }

◆ MassConstr()

void Vertex::MassConstr ( double m )

add a mass constraint for all tracks, $m$ = mass of mother particle

                                   {
     MassC* massc = new MassC(m);
     const RelationList::iterator last = end();
     for(RelationList::iterator it = begin(); it != last; ++it) {
       massc->push_back(*it);
     }
  
     xmass.push_back(massc);
     return;
   }

◆ ndf()

unsigned short int Vertex::ndf ( ) const

degrees of freedom of vertex fit

239 { return ((v_Mvalid == true) ? 2*ntracks()-3 + xmass.size() : 2*ntracks()-3); }

◆ nMassConstr()

unsigned int VERTEX::Vertex::nMassConstr ( ) const

number of mass constraints

◆ ntracks()

unsigned short int Vertex::ntracks ( ) const

no of tracks in vertex $n$

240 { return size(); }

◆ operator-()

const double Vertex::operator- ( const Vertex & rhs ) const

                                                         {
     return sqrt(sqr(vx()-rhs.vx()) + sqr(vy()-rhs.vy()) + sqr(vz()-rhs.vz()));
   }

◆ operator=()

const Vertex & Vertex::operator= ( const Vertex & rhs )

                                                    {
     if( this == &rhs ) { return *this; }
     
     Track::operator=(rhs);
  
     v_use_kalman   = rhs.v_use_kalman;
     v_valid        = rhs.v_valid;
     v_Mvalid       = rhs.v_Mvalid;
     v_mother       = rhs.v_mother;
     v_angdist      = rhs.v_angdist;
     v_angle        = rhs.v_angle;
     v_dist         = rhs.v_dist;
     v_bk13         = rhs.v_bk13;
     v_chi2         = rhs.v_chi2;
  
     kal_xv         = rhs.kal_xv;
     kal_xvs        = rhs.kal_xvs;
     v_CS           = rhs.v_CS;
  
     // possible if tracks & xmass contain SmartPointers to the objects
     //      tracks         = rhs.tracks;
     //      xmass          = rhs.xmass;
  
     return *this;
   }

◆ operator==()

const bool Vertex::operator== ( const Vertex & rhs ) const

                                                        {
     return RelationList::operator==(rhs);
   }

◆ pmaxfrac()

double Vertex::pmaxfrac ( ) const

$\max|p_i/p_j|$ max. momentum fraction

                                 {
     double maxfrac = 0.;
  
     const RelationList::const_iterator last = end();
     for(RelationList::const_iterator tit = begin(); tit != last; ++tit) {
       double p1 = fabs(tit->track.p());
  
       RelationList::const_iterator it=tit;
       for(++it; it != last; ++it) {
     double p2 = fabs(it->track.p());
     double frac = max(p1/p2, p2/p1);
     if( frac > maxfrac ) maxfrac = frac;
       }
     }
  
     return maxfrac;
   }

◆ print()

std::ostream & Vertex::print ( std::ostream & os ) const

called by operator<<()

                                                 {
     os.setf(ios::right,ios::adjustfield);
     os.setf(ios::fixed,ios::floatfield);
     if(v_mother == true) {
       os << " Mother Track: ";
       Track::print(os);
       os << endl;
     }
 //       if (1) return os;
     os << " Vertex:"
        << " vx: "       << setw(6) << setprecision(3) << vx()
        << " vy: "       << setw(6) << setprecision(3) << vy()
        << " vz: "       << setw(7) << setprecision(3) << vz()
 //       if (1) return os;
        << " chi2: "     << setw(6) << setprecision(2) << chi2()
        << " prob: "     << setw(6) << setprecision(2) << prob()
        << " angle: "    << setw(5) << setprecision(3) << angle() 
        << " dist: "     << setw(6) << setprecision(4) << dist()
        << " mass: "     << setw(6) << setprecision(4) << mass(v_use_kalman)
        << " no of tracks: " << ntracks() 
        << endl;
  
     if(xmass.size() > 0) {
       os << "Mass constraints: " << endl;
       const MassC_cit end = xmass.end();
       for(MassC_cit it = xmass.begin(); it != end; ++it) {
     os << "   " << **it << endl;
       }
     }
  
 #ifdef VtDEBUG
     os << " Relations: ";
  
     RelationList::print(os);
 #endif
  
     os << "Vertex cov. matrix: " << CS();
  
 #ifdef VtDEBUG
     if(v_covn != 0) {
       cout << " bigcov: " << *v_covn << endl;
     }
 #endif
  
     return os;
   }

◆ prob()

float Vertex::prob ( ) const

upper tail $\chi^2$ probability

237 { return CERN_prob(chi2(), ndf()); }

VERTEX::Vertex::ndf

unsigned short int ndf() const

degrees of freedom of vertex fit

Definition: VtVertex.C:238

◆ push_back()

void Vertex::push_back ( Track & track )

                                      {
     createRelation(track, *this);
     invalid();
   }

◆ remove()

void Vertex::remove ( Relation * r )

virtual

unbook relation

Reimplemented from VERTEX::RelationList.

                                 {
     RelationList::remove(r);
     invalid();
   }

◆ remove_last()

const bool Vertex::remove_last ( )

remove last track from vertex, do inverse Kalman filter

                                  {
     return VtRemoveTrack(*back());
   }

◆ remove_worst()

void Vertex::remove_worst ( )

remove track with biggest $\chi^2$ contribution, do inverse Kalman filter

                             {
     double prchi2 = 0.;
     double maxchi2 = 0.;
     Relation* worst = 0;
  
     // a sorted list would be better here
     const iterator last = end();
     for(iterator it=begin(); it!=last; ++it) {
       const Kalman& k = it->kalman;
       double dchi2 = k.chi2() - prchi2;
       prchi2 = k.chi2();
  
       if(dchi2 > maxchi2){
     maxchi2 = dchi2;
     worst = &(*it);
       }
     }
     
     if(worst != 0) VtRemoveTrack(*worst);
     //    if(worst != 0) remove(worst);
     return;
   }

◆ rmsDistAngle()

double Vertex::rmsDistAngle ( void ) const

calc rms dist and rms angle

                                         {
  
     // and the rms distance to the vertex
     v_dist = 0.;
     double tsq, rsq, rt, t[3], r[3];
     double sumtx = 0., sumty = 0.;
     double sumtx2 = 0., sumty2 = 0.;
     const RelationList::const_iterator last = end();
     for(RelationList::const_iterator pos = begin(); pos != last; ++pos) {
       const Track& track = pos->track;
  
       t[0] = track.tx();  r[0] = track.x() - kal_xvs[0];
       t[1] = track.ty();  r[1] = track.y() - kal_xvs[1];
       t[2] = 1;           r[2] = track.z() - kal_xvs[2];
  
       sumtx += t[0]; sumtx2 += sqr(t[0]);
       sumty += t[1]; sumty2 += sqr(t[1]);
  
       tsq = sqr(t[0]) + sqr(t[1]) + sqr(t[2]);
       rsq = sqr(r[0]) + sqr(r[1]) + sqr(r[2]);
       rt  = r[0]*t[0] + r[1]*t[1] + r[2]*t[2];
  
       v_dist += rsq - sqr(rt)/tsq;
     } 
  
     unsigned int v_nrseg = size();
     v_angle = sqrt(sumtx2*v_nrseg - sqr(sumtx) + sumty2*v_nrseg -
            sqr(sumty))/static_cast<double>(v_nrseg); 
     
     v_dist /= static_cast<double>(v_nrseg-1);
     v_dist  = sqrt(v_dist);
     v_angdist = true;
     return v_dist;
   }

◆ set_invalid()

void VERTEX::Vertex::set_invalid ( )

mark vertex as invalid

◆ track_chi2()

double Vertex::track_chi2 ( int i )

get track $\chi^2$ contribution

                                    {
     double prchi2 = 0.;
     double dchi2  = 0.;
     // a sorted list would be better here
     const iterator last = end();
     int ntr = 0;
     for(iterator it=begin(); it!=last; ++it) {
       const Kalman& k = it->kalman;
       dchi2 = k.chi2() - prchi2;
       if (ntr == i) return dchi2;
       prchi2 = k.chi2();
       ntr++;
     }
     
     return dchi2;
   }

◆ track_worst()

int Vertex::track_worst ( )

get track with biggest $\chi^2$ contribution

                           {
     double prchi2 = 0.;
     double maxchi2 = 0.;
  
     // a sorted list would be better here
     const iterator last = end();
     int num = -1, ntr = 0;
     for(iterator it=begin(); it!=last; ++it) {
       const Kalman& k = it->kalman;
       double dchi2 = k.chi2() - prchi2;
       prchi2 = k.chi2();
  
       if(dchi2 > maxchi2){
     maxchi2 = dchi2;
     num = ntr;
       }
       ntr++;
     }
     
     return num;
   }

◆ use_kalman() [1/2]

bool VERTEX::Vertex::use_kalman ( ) const

◆ use_kalman() [2/2]

void VERTEX::Vertex::use_kalman ( const bool use )

use refitted track parameters or not

◆ use_momentum()

void Vertex::use_momentum ( const bool use )

set for all tracks whether momentum should be used or not

                                           {
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
       tit->kalman.use_momentum(use);
     }
     return;
   }

◆ valid()

bool VERTEX::Vertex::valid ( ) const

is vertex valid?

◆ VtEstimateVertex()

const bool Vertex::VtEstimateVertex ( )

weighted 2-D min. dist.

                                       {
     //return VtEstimateVertexMath(kal_xvs[0], kal_xvs[1], kal_xvs[2]);
     return VtEstimateVertexMathTA(kal_xvs[0], kal_xvs[1], kal_xvs[2]);
   }

◆ VtEstimateVertexMath()

const bool Vertex::VtEstimateVertexMath	(	double &	x,
		double &	y,
		double &	z
	)

estimate Vertex without changing Vertex object

                                                                          {
  
     double tx_sum  = 0.;
     double x_sum   = 0.;
     double xw_sum  = 0.;
     double xtx_sum = 0.;
     double tx2_sum = 0.;
     double ty_sum  = 0.;
     double y_sum   = 0.;
     double yw_sum  = 0.;
     double yty_sum = 0.;
     double ty2_sum = 0.;
  
     // fill cumulants
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
       const Track& track = tit->track;
  
       double x        = track.x();
       double y        = track.y();
       double tx       = track.tx();
       double ty       = track.ty();
       double xweight  = 1./track.V().x();
       double yweight  = 1./track.V().y();
       double xweight2 = xweight*xweight;
       double yweight2 = yweight*yweight;
  
       tx_sum  += tx * xweight;
       x_sum   += x  * xweight;
       xw_sum  += xweight;
       xtx_sum += x  * tx * xweight2;
       tx2_sum += tx * tx * xweight2;
  
       ty_sum  += ty * yweight;
       y_sum   += y  * yweight;
       yw_sum  += yweight;
       yty_sum += y  * ty * yweight2;
       ty2_sum += ty * ty * yweight2;
  
     } // for track
  
     double det = -tx2_sum - ty2_sum + tx_sum*tx_sum/xw_sum + ty_sum*ty_sum/yw_sum;
  
 //#ifdef VtDEBUG    
     cout << " txsum: " << tx_sum << " xsum: " << x_sum << " xwsum: " << xw_sum
      << " xtxsum: " << xtx_sum << " tx2sum: " << tx2_sum << endl;
     cout << " tysum: " << ty_sum << " ysum: " << y_sum << " ywsum: " << yw_sum
      << " ytysum: " << yty_sum << " ty2sum: " << ty2_sum << endl;
     cout << " det: " << det << endl;
 //#endif
  
     if(det == 0.) {
 //#ifdef VtDEBUG
       cout << "VtVertex::VtEstimateVertex error: determinant == 0.!" << endl;
 //#endif
       return false;
     }
  
 //      kal_xvs[2] = ( xtx_sum + yty_sum - tx_sum*x_sum/xw_sum - ty_sum*y_sum/yw_sum ) / det;
 //      kal_xvs[0] = ( x_sum + tx_sum * vz() ) / xw_sum;
 //      kal_xvs[1] = ( y_sum + ty_sum * vz() ) / yw_sum;
  
     z = ( xtx_sum + yty_sum - tx_sum*x_sum/xw_sum - ty_sum*y_sum/yw_sum ) / det;
     x = ( x_sum + tx_sum * z ) / xw_sum;
     y = ( y_sum + ty_sum * z ) / yw_sum;
  
     v_valid = true;
     return true;
   }

◆ VtEstimateVertexMathTA()

const bool Vertex::VtEstimateVertexMathTA	(	double &	x,
		double &	y,
		double &	z
	)

estimate Vertex without changing Vertex object

                                                                               {
   // Tomoko Ariga version inserted by VT 1/09/2011
   // calc vertex point with given segments. just topological calculation
  
     double Det,Ax,Ay,Az,Bx,By,Bz,Cx,Cy,Cz,Dx,Dy,Dz;
     Ax=0.;Ay=0.;Az=0.;Bx=0.;By=0.;Bz=0.;Cx=0.;Cy=0.;Cz=0.;Dx=0.;Dy=0.;Dz=0.;
  
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
         const Track& track = tit->track;
  
         double x        = track.x();
         double y        = track.y();
         double ax       = track.tx();
         double ay       = track.ty();
         double az = 1.0;
     double z  = track.z();
     double a = ax*ax+ay*ay+az*az;
     double c = -ax*x-ay*y-az*z;
     double b = (ax*ax+ay*ay);
 //      double w = 1.0/a/a; // weight for small angle tracks.
     double w = 1.0; // no weight
  
     Ax+=2.0*w/b*( a*(ay*ay+az*az) );
     Bx+=2.0*w/b*( -a*ax*ay );
     Cx+=2.0*w/b*( -a*ax*az );
     Dx+=2.0*w/b*( -(a*x+c*ax)*(ax*ax-a)-(a*y+c*ay)*ax*ay-(a*z+c*az)*az*ax );
  
     Ay+=2.0*w/b*( -a*ay*ax );
     By+=2.0*w/b*( a*(az*az+ax*ax) );
     Cy+=2.0*w/b*( -a*ay*az );
     Dy+=2.0*w/b*( -(a*y+c*ay)*(ay*ay-a)-(a*z+c*az)*ay*az-(a*x+c*ax)*ax*ay );
  
     Az+=2.0*w/b*( -a*az*ax );
     Bz+=2.0*w/b*( -a*az*ay );
     Cz+=2.0*w/b*( a*b );
     Dz+=2.0*w/b*( -(a*z+c*az)*(az*az-a)-(a*x+c*ax)*az*ax-(a*y+c*ay)*ay*az );
  
      }
  
     Det=fabs( Ax*(By*Cz-Cy*Bz)-Bx*(Ay*Cz-Cy*Az)+Cx*(Ay*Bz-By*Az) );
     if( Det < 0.00000001 ) {
       cout << "VtVertex::VtEstimateVertexMathTA error: determinant is too small! " << Det << endl;
       return false;
     }
     xx=( (By*Cz-Cy*Bz)*Dx-(Bx*Cz-Cx*Bz)*Dy+(Bx*Cy-Cx*By)*Dz)/Det;
     yy=(-(Ay*Cz-Cy*Az)*Dx+(Ax*Cz-Cx*Az)*Dy-(Ax*Cy-Cx*Ay)*Dz)/Det;
     zz=( (Ay*Bz-By*Az)*Dx-(Ax*Bz-Bx*Az)*Dy+(Ax*By-Bx*Ay)*Dz)/Det;
     
     v_valid = true;
     return true;
   }

◆ VtFilter()

const bool Vertex::VtFilter ( )

                               {
     
     // initialize by vertex position & cov. matrix
     //    const VtVector*    prkal_xv = &kal_xv0;
     const VtVector*    prkal_xv = &kal_xvs;
     const VtSymMatrix* prCINV   = &v_CINV;
     double             prChi2   = 0.;
     v_bk13 = vz();
  
     // cout << "bk13: " << v_bk13 << endl;
  
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
       Kalman& kalman = tit->kalman;
  
       //      cout << "----> FILTER TRACK " << track;
 #ifdef VtDEBUG
       cout << endl;
 #endif
  
       // compute contributions of track k
       if(kalman.filter(v_bk13, *prCINV, *prkal_xv) == false) {
     return false; // forget this vertex
       }
  
       kalman.filter_chi2(v_bk13, prChi2, *prCINV, *prkal_xv);
     
       //      cout << "Track-Kalman chi2: " << kalman.chi2() << endl;
  
       // store pointers to previous (k-1) values
       prkal_xv = &kalman.xv();
       prCINV   = &kalman.CINV();
       prChi2   = kalman.chi2();
     } // for tracks
  
     return true;
   }

◆ VtFit()

const bool Vertex::VtFit ( )

Kalman filter + smoother.

                            {
     if(VtFilter() == false) return false;
     else {
       VtSmoothX();
       VtSmoothQ();
     }
     return true;
   }

◆ VtInverseFilter()

const bool Vertex::VtInverseFilter ( ) const

                                            {
     
     // get last track (in Vt: x_n)
     const Kalman& trackn = back()->kalman;
  
     const VtSymMatrix* prCINV = &v_CINV;
  
     const RelationList::const_iterator last = end();
     for(RelationList::const_iterator tit = begin(); tit != last; ++tit) {
       Kalman& kalman = tit->kalman;
       
       if(kalman.inverse_filter(v_bk13, trackn.CINV(), *prCINV,
                    trackn.xv()) == false) {
     return false;
       }
       prCINV   = &kalman.CINV();
     }
  
     return true;
   }

◆ VtMass()

const bool Vertex::VtMass ( )

Kalman filter with mass constraints.

                             {
     v_Mvalid = false;
  
     // should be assured:
     use_momentum(true);
     if(valid() == false && findVertexVt() == false) return false;
  
     // check kinematic boundaries
     const MassC_cit xend = xmass.end();
     for(MassC_cit it = xmass.begin(); it != xend; ++it) {
       if((*it)->check_kinematic() == false)
     return false;
     }
  
     // fill general state vectors (in Vt: alp(1,...), alpc(1,...))
     VtVector alpha (bigdim());
     VtVector alphac(bigdim());
  
     // initialize general state vector alphac (in Vt: alpc)
     alphac.place_at(xvs(),0);
  
     unsigned int i=3;
     const iterator last = end();
     for(iterator it = begin(); it != last; ++it,i+=3) {
       it->kalman.alpc_init();
       alphac.place_at(it->kalman.alpc(),i);
     }
  
 #ifdef VtDEBUG
     cout << " alphac0: " << alphac << endl;
 #endif
  
     double chi2old = -100.;
     double dchi2 = 0.;
     double mc_chi2 = 0.;
     const unsigned int maxiter = 100 * nMassConstr();
  
     // declare some matrices which are needed outside do-loop
     VtSymMatrix gaVga(nMassConstr());
     VtMatrix    Vgam(bigdim(),nMassConstr());
  
     // iteration
     unsigned int iter = 0;
     do {
  
       if(iter > 0) {
     // redo filter and smoother
     kal_xvs.copy(alphac,0);
     unsigned int i=3;
     for(iterator it = begin(); it != last; ++it,i+=3) {
       Kalman& kalm = it->kalman;
       kalm.qvs_nc().copy(alphac,i);
       kalm.alpc_init();
     }
     if((v_valid = VtFit()) == false)
       return false;
       }
  
       // fill general state vector alpha (in Vt: alp)
       alpha.place_at(xvs(),0);
  
       i=3;
       for(iterator it = begin(); it != last; ++it,i+=3) {
     Kalman& kalm = it->kalman;
  
     kalm.alp_init();
     kalm.calc_ealpc(); // n in Vt (eq. 40)
     alpha.place_at(kalm.alp(),i);
       }
  
 #ifdef VtDEBUG
       cout << " alpha: " << alpha << endl;
 #endif
  
       // get general covariance matrix v_covn (in Vt: Valp)
       const VtSymMatrix& Valp = bigcov();
  
 #ifdef VtDEBUG
       cout << "Valp: " << Valp << endl;
 #endif
  
       // build up matrix gam (in Vt: Gamma^T, eq. 34, 44, 45, 46)
       VtMatrix gam(bigdim(),nMassConstr());
  
       // vector \vect{g}(\vect{\alpha}^(0))
       VtVector galp(nMassConstr());
  
       unsigned int j=0;
       for(MassC_cit mit = xmass.begin(); mit != xend; ++mit,++j) {
     MassC* mC = *mit;
     mC->calc_Sl();
     galp[j] = mC->galp();
  
     i=3;
     for(iterator it = begin(); it != last; ++it,i+=3) {
       Relation& r = *it;
       if(find(r.massc.begin(),r.massc.end(),mC) != r.massc.end()) {
         gam.place_at(mC->gam(r),i,j);
       }
     }
       } // end of for(MassC_cit)
  
 #ifdef VtDEBUG
       cout << " galp: " << galp << endl;
       cout << " gam: " << gam << endl;
 #endif
  
       // calc (Gamma * V_alpha * Gamma^T)^-1 (in Vt: eq. 36,37,39)
       gaVga = Valp.product(gam);
       if(gaVga.VtDsinv() == false)
     return false;
       
 #ifdef VtDEBUG
       cout << " gaVga^-1: " << gaVga << endl;
 #endif
  
       // calc g(alpha) + Gamma * (alpha - alpha_c)
       VtVector gnum = galp + gam.T() * (alpha - alphac);
  
       // calc \vect{\mu} (in Vt: xmu, eq. 36)
       VtVector xmu = gaVga * gnum;
  
 #ifdef VtDEBUG
       cout << " xmu: " << xmu << endl;
 #endif
  
       // calc V_alpha*Gamma^T
       VtMatrix Vgam = Valp * gam;
  
       // calc alpha_c (in Vt: eq. 37)
       alphac = alpha - Vgam * xmu;
       
 #ifdef VtDEBUG
       cout << " alphac: " << alphac << endl;
 #endif
  
       // calc chi^2 (in Vt: chi2m, eq. 39)
       // this is an incremental chi2!
       // (chi2 due to mass constraint, to be added to
       // chi2 of unconstraint vertex fit!)
       mc_chi2 = gaVga.product(gnum);
  
       double chi  = mc_chi2 + chi2n();
       dchi2 = fabs(chi - chi2old);
       chi2old = chi;
  
 #ifdef VtDEBUG
       cout << "chi2m: " << setprecision(7) << mc_chi2 
        << " chi: " << chi2n()
        << " fabs(chi2m + chi2() - chi2old: " << dchi2
        << endl;
 #endif
     
     } while(++iter < maxiter && dchi2 > 0.00001);
  
 #ifdef VtDEBUG
     cout << " iterations: " << iter << endl;
 #endif
     if(iter >= maxiter) return false;
  
     VtSymMatrix& cov = *v_covn;
  
     // v_covn = Valp => v_covn = Valp - gaVga.product(Vgam.T())
     // calc v_covn (in Vt: cov(alpha_c), eq. 38)
     cov -= gaVga.product(Vgam.T());
  
 #ifdef VtDEBUG
     cout << " v_covn: " << cov << endl;
 #endif
  
     // get results
     v_CS.copy(cov);
     // add mass constraint chi2 term to vertex fit term
     v_chi2 += mc_chi2;
  
     kal_xvs.copy(alphac,0);
     i=3;
     for(iterator it = begin(); it != last; ++it,i+=3) {
       it->kalman.qvs_nc().copy(alphac,i);
     }
     
     return v_Mvalid = true;
   } // end of VtMass()

◆ VtRemoveTrack()

const bool Vertex::VtRemoveTrack ( Relation & track )

remove track using inverse Kalman filter

                                                    {
     if(ntracks() <= 2 || valid() == false) return false; // exclude weird cases
  
     if(VtInverseFilter() == false) return false;
  
     // get last track (in Vt: x_n)
     const Kalman& trackn = back()->kalman;
  
     // compute chi2_n -> chi2_n - chi2_k,S (eq. 22)
     double dchi2 = trackn.chi2() - trackk.kalman.chi2s();
     
     // new smoothed vertex position
     kal_xvs = trackk.kalman.xnk();
  
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
       // calculate qvs
       tit->kalman.calc_qvs(kal_xvs);
     }
  
     remove(&trackk);
  
     v_chi2 = back()->kalman.set_chi2(dchi2);
  
     // refit the vertex
     return v_valid = VtFit();
   }

◆ VtSmoothQ()

void Vertex::VtSmoothQ ( )

smooth momenta

                          {
    
 #ifdef VtDEBUG
     cout << "---> VtSmoothQ called! <---" << endl;
 #endif
  
     const Kalman& trackn = back()->kalman;
     const RelationList::iterator last = end();
     for(RelationList::iterator tit = begin(); tit != last; ++tit) {
       tit->kalman.smooth(v_bk13, kal_xvs, trackn.C());
     }
     
     return;
   }

◆ VtSmoothX()

void Vertex::VtSmoothX ( )

smooth vertex position

                          {
     // get last track (in Vt: x_n)
     const Kalman& kaln = back()->kalman;
     
 #ifdef VtDEBUG
     cout << "--> VtSmooth called!" << endl;
 #endif
  
     // update vertex position saved in xvs!!
     kal_xvs = kaln.xv(); 
  
 #ifdef VtDEBUG
     cout << "kal_xvs: "    << kal_xvs  << endl;
     cout << "track->C(): " << kaln.C() << endl;
 #endif
  
     // save smoothed covariance matrix of vertex
     v_CS = kaln.C();
  
     // save chi2
     v_chi2 = kaln.chi2();
  
 #ifdef VtDEBUG
     cout << "v_CS: " << v_CS << endl;
 #endif
  
     return;
   }

◆ vtx_cov_x()

float Vertex::vtx_cov_x ( ) const

$\sigma_x^2$ of vertex

254 { return v_CS(0,0); }

◆ vtx_cov_y()

float Vertex::vtx_cov_y ( ) const

$\sigma_y^2$ of vertex

255 { return v_CS(1,1); }

◆ vtx_cov_z()

float Vertex::vtx_cov_z ( ) const

$\sigma_z^2$ of vertex

256 { return v_CS(2,2); }

◆ vx()

float Vertex::vx ( ) const

$x$ of vertex

233 { return kal_xvs[0]; }

◆ vxerr()

double VERTEX::Vertex::vxerr ( ) const

$\sqrt{\sigma_{vx}^2}$ vertex $x$-error

◆ vy()

float Vertex::vy ( ) const

$y$ of vertex

234 { return kal_xvs[1]; }

◆ vyerr()

double VERTEX::Vertex::vyerr ( ) const

$\sqrt{\sigma_{vy}^2}$ vertex $y$-error

◆ vz()

float Vertex::vz ( ) const

$z$ of vertex

235 { return kal_xvs[2]; }

◆ vzerr()

double VERTEX::Vertex::vzerr ( ) const

$\sqrt{\sigma_{vz}^2}$ vertex $z$-error

◆ xv()

const MATRIX::VtVector& VERTEX::Vertex::xv ( ) const

◆ xvs()

const MATRIX::VtVector& VERTEX::Vertex::xvs ( ) const

vector of vertex-positions

Member Data Documentation

◆ kal_xv

MATRIX::VtVector VERTEX::Vertex::kal_xv

private

◆ kal_xvs

MATRIX::VtVector VERTEX::Vertex::kal_xvs

private

◆ tracks

Track_v VERTEX::Vertex::tracks

private

◆ v_angdist

bool VERTEX::Vertex::v_angdist

mutableprivate

◆ v_angle

double VERTEX::Vertex::v_angle

mutableprivate

◆ v_bk13

double VERTEX::Vertex::v_bk13

private

◆ v_chi2

double VERTEX::Vertex::v_chi2

private

◆ v_CINV

MATRIX::VtSymMatrix VERTEX::Vertex::v_CINV

private

◆ v_covn

MATRIX::VtSymMatrix* VERTEX::Vertex::v_covn

private

◆ v_CS

MATRIX::VtSymMatrix VERTEX::Vertex::v_CS

private

◆ v_dist

double VERTEX::Vertex::v_dist

mutableprivate

◆ v_mother

bool VERTEX::Vertex::v_mother

private

◆ v_Mvalid

bool VERTEX::Vertex::v_Mvalid

private

◆ v_use_kalman

bool VERTEX::Vertex::v_use_kalman

private

◆ v_valid

bool VERTEX::Vertex::v_valid

private

◆ xmass

MassC_v VERTEX::Vertex::xmass

private

The documentation for this class was generated from the following files:

src/libVt++/vt++/include/VtVertex.hh
src/libVt++/vt++/VtVertex.C

Public Member Functions

— Operators —

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ Vertex() [1/4]

◆ Vertex() [2/4]

◆ Vertex() [3/4]

◆ Vertex() [4/4]

◆ ~Vertex()

Member Function Documentation

◆ add_track()

◆ addMassConstr()

◆ angle()

◆ bigcov()

◆ bigdim()

◆ calc_mother()

◆ calc_mother_cov()

◆ chi2()

◆ chi2l()

◆ chi2n()

◆ ClassDef()

◆ clear()

◆ clearMassConstr()

◆ covn()

◆ CS()

◆ dist()

◆ distance() [1/3]

◆ distance() [2/3]

◆ distance() [3/3]

◆ erase()

◆ findVertex2D()

◆ findVertex3D()

◆ findVertexVt()

◆ get_track()

◆ invalid()

◆ mass() [1/6]

◆ mass() [2/6]

◆ mass() [3/6]

◆ mass() [4/6]

◆ mass() [5/6]

◆ mass() [6/6]

◆ massCC()

◆ MassConstr()

◆ ndf()

◆ nMassConstr()

◆ ntracks()

◆ operator-()

◆ operator=()

◆ operator==()

◆ pmaxfrac()

◆ print()

◆ prob()

◆ push_back()

◆ remove()

◆ remove_last()

◆ remove_worst()

◆ rmsDistAngle()

◆ set_invalid()

◆ track_chi2()

◆ track_worst()

◆ use_kalman() [1/2]

◆ use_kalman() [2/2]

◆ use_momentum()

◆ valid()

◆ VtEstimateVertex()

◆ VtEstimateVertexMath()

◆ VtEstimateVertexMathTA()

◆ VtFilter()

◆ VtFit()

◆ VtInverseFilter()

◆ VtMass()

◆ VtRemoveTrack()

◆ VtSmoothQ()

◆ VtSmoothX()

◆ vtx_cov_x()

◆ vtx_cov_y()

◆ vtx_cov_z()

◆ vx()

◆ vxerr()