Trk::JetFitterRoutines Class Reference

Februar 2007 (c) Atlas Detector Reconstruction Software. More...

#include <JetFitterRoutines.h>

Inheritance diagram for Trk::JetFitterRoutines:

Collaboration diagram for Trk::JetFitterRoutines:

Public Member Functions
virtual StatusCode	initialize () override
	JetFitterRoutines (const std::string &t, const std::string &n, const IInterface *p)
	Constructor. More...
	~JetFitterRoutines ()
	Destructor. More...
void	initializeToMinDistancesToJetAxis (VxJetCandidate *) const
	This method provides the initialization of all the tracks in the fit to the position of minimum distance to the initial flight axis. More...
void	performTheFit (VxJetCandidate *myJetCandidate, int num_maxiterations=30, bool treat_sign_flip=true, int num_signflip_maxiterations=30, double deltachi2_convergence=0.001) const
	This is the method where the fit is actually done. More...
void	smoothAllVertices (VxJetCandidate *myJetCandidate) const
	triggers the smoothing of all vertices (the tracks in the fit are updated using the constraint provided by the vertices along the jet axis they were fitted to) More...
void	updateAllVertices (VxJetCandidate *myJetCandidate) const
	One iteration of the Kalman Updated of all tracks to the actual fit is performed. More...
std::pair< double, bool >	fastProbabilityOfMerging (const VxVertexOnJetAxis *, const VxVertexOnJetAxis *, const VxJetCandidate *) const
	Calculates in a very fast way the probability two vertices along the jet axis to be compatible with each other, taking the constraint to be compatible with the previous fitted flight axis into account. More...
double	fullProbabilityOfMerging (const VxVertexOnJetAxis *firstVertex, const VxVertexOnJetAxis *secondVertex, const VxJetCandidate *myJetCandidate, int num_maxiterations=20, bool treat_sign_flip=true, int num_signflip_maxiterations=10, double deltachi2_convergence=1e-2) const
	Calculates in a complete way the probability two vertices along the jet axis to be compatible with each other, taking the constraint to be compatible with the previous fitted flight axis into account. More...
void	fillTableWithFullProbOfMerging (VxJetCandidate *myJetCandidate, int num_maxiterations=20, bool treat_sign_flip=true, int num_signflip_maxiterations=10, double deltachi2_convergence=1e-3, double threshold_probability=1e-5) const
	The VxClusteringTable of the VxJetCandidate provided in input is created, computing all the needed probabilities of merging between all the vertices two by two using a full probability estimation, with the complete fit performed with the options provided. More...
void	fillTableWithFastProbOfMerging (VxJetCandidate *myJetCandidate) const
	The VxClusteringTable of the VxJetCandidate provided in input is created, computing all the needed probabilities of merging between all the vertices two by two using the fast probability estimation (using for all combinations of vertices the fastProbabilityOfMerging method of the same class) More...
void	deleteVertexFromJetCandidate (VxVertexOnJetAxis *vertexToDelete, VxJetCandidate *myJetCandidate) const
	Deltes a vertex from the VxJetCandidate, doing everything is needed to preserve the internal coherence of the VxJetCandidate object (e.g. More...
ServiceHandle< StoreGateSvc > &	evtStore ()
	The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc. More...
const ServiceHandle< StoreGateSvc > &	evtStore () const
	The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc. More...
const ServiceHandle< StoreGateSvc > &	detStore () const
	The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc. More...
virtual StatusCode	sysInitialize () override
	Perform system initialization for an algorithm. More...
virtual StatusCode	sysStart () override
	Handle START transition. More...
virtual std::vector< Gaudi::DataHandle * >	inputHandles () const override
	Return this algorithm's input handles. More...
virtual std::vector< Gaudi::DataHandle * >	outputHandles () const override
	Return this algorithm's output handles. More...
Gaudi::Details::PropertyBase &	declareProperty (Gaudi::Property< T, V, H > &t)
Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, SG::VarHandleKey &hndl, const std::string &doc, const SG::VarHandleKeyType &)
	Declare a new Gaudi property. More...
Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, SG::VarHandleBase &hndl, const std::string &doc, const SG::VarHandleType &)
	Declare a new Gaudi property. More...
Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, SG::VarHandleKeyArray &hndArr, const std::string &doc, const SG::VarHandleKeyArrayType &)
Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, T &property, const std::string &doc, const SG::NotHandleType &)
	Declare a new Gaudi property. More...
Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, T &property, const std::string &doc="none")
	Declare a new Gaudi property. More...
void	updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream &	msg () const
MsgStream &	msg (const MSG::Level lvl) const
bool	msgLvl (const MSG::Level lvl) const

Static Public Member Functions
static const InterfaceID &	interfaceID ()

Protected Member Functions
void	renounceArray (SG::VarHandleKeyArray &handlesArray)
	remove all handles from I/O resolution More...
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void >	renounce (T &h)
void	extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
	Add StoreName to extra input/output deps as needed. More...

Private Types
typedef ServiceHandle< StoreGateSvc >	StoreGateSvc_t

Private Member Functions
std::pair< double, bool >	fastProbabilityOfMergingWithPrimary (const VxVertexOnJetAxis *otherVertex, const VxJetCandidate *myJetCandidate) const
	Internal method to calculate fast probability of merging, for merging with primary vertex. More...
std::pair< double, bool >	fastProbabilityOfMergingNoPrimary (const VxVertexOnJetAxis *, const VxVertexOnJetAxis *, const VxJetCandidate *myJetCandidate) const
	Internal method to calculate fast probability of merging, for merging with non primary vertex. More...
void	fillTableWithProbOfMerging (VxJetCandidate *myJetCandidate, bool fullcomputation, int num_maxiterations=20, bool treat_sign_flip=true, int num_signflip_maxiterations=10, double deltachi2_convergence=1e-2, double threshold_probability=1e-3) const
	Internal method to fill the VxClusteringTable of the VxJetCandidate object, independently on if the fullcomputation is requested or not. More...
bool	checkJetCandidate (VxJetCandidate *) const
	Internal method to provide a check if the VxJetCandidate has been initialized in a consistent way. More...
void	copyRecoPositionsToLinearizationPositions (VxJetCandidate &myJetCandidate) const
	Method to copy new reco positions to linearization positions after checking new positions are inside the ID. More...
Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
	specialization for handling Gaudi::Property<SG::VarHandleKey> More...
Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyArrayType &)
	specialization for handling Gaudi::Property<SG::VarHandleKeyArray> More...
Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleType &)
	specialization for handling Gaudi::Property<SG::VarHandleBase> More...
Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T, V, H > &t, const SG::NotHandleType &)
	specialization for handling everything that's not a Gaudi::Property<SG::VarHandleKey> or a <SG::VarHandleKeyArray> More...

Private Attributes
SG::ReadCondHandleKey< AtlasFieldCacheCondObj >	m_fieldCacheCondObjInputKey {this, "AtlasFieldCacheCondObj", "fieldCondObj", "Name of the Magnetic Field conditions object key"}
ToolHandle< JetFitterInitializationHelper >	m_initializationHelper
ToolHandle< JetFitterHelper >	m_helper
ToolHandle< KalmanVertexOnJetAxisUpdator >	m_updator
ToolHandle< KalmanVertexOnJetAxisSmoother >	m_smoother
ToolHandle< TrkDistanceFinderNeutralCharged >	m_minDistanceFinder
ToolHandle< TrkDistanceFinderNeutralNeutral >	m_minDistanceFinderNeutral
bool	m_fast
double	m_maxDRshift
bool	m_noPrimaryVertexRefit
double	m_maxR
double	m_maxZ
StoreGateSvc_t	m_evtStore
	Pointer to StoreGate (event store by default) More...
StoreGateSvc_t	m_detStore
	Pointer to StoreGate (detector store by default) More...
std::vector< SG::VarHandleKeyArray * >	m_vhka
bool	m_varHandleArraysDeclared

Detailed Description

Februar 2007 (c) Atlas Detector Reconstruction Software.

Author

Giacinto Piacquadio (University of Freiburg),

Christian Weiser (University of Freiburg)

This class contains the core functions for the "JetFitter" algorithm.

This is a new algorithm which uses a "Kalman filter" based approach to implement the "ghost track method" used for the first time by SLD (reference to be inserted here).

All the tracks inside a jet can be fitted with the constraint of lying on the same 'flight axis' from the primary vertex towards the jet. This class provides all the methods which are necessary to perform this fit and to evaluate the probability of two or more tracks to build one cluster on the jet axis.

The real finding algorithm is then implemented in the "SV" finder algorithm. In this class only the methods for 'fitting' are implemented.

There is no clear input and clear output in a similar fit, so this works more like a continous interaction between the finder and the fitter, which is provided in the form in a set of routines to use.

Definition at line 56 of file JetFitterRoutines.h.

Member Typedef Documentation

StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t

privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

JetFitterRoutines()

Trk::JetFitterRoutines::JetFitterRoutines	(	const std::string &	t,
		const std::string &	n,
		const IInterface *	p
	)

Constructor.

Definition at line 74 of file JetFitterRoutines.cxx.

                                           {
    
    AthAlgTool::initialize().ignore();
    ATH_CHECK( m_fieldCacheCondObjInputKey.initialize() );
 
    //retrieving the udator itself   
    ATH_CHECK( m_helper.retrieve()  );    
 
    ATH_CHECK( m_initializationHelper.retrieve() );

~JetFitterRoutines()

Trk::JetFitterRoutines::~JetFitterRoutines ( )

default

Destructor.

Member Function Documentation

checkJetCandidate()

bool Trk::JetFitterRoutines::checkJetCandidate ( VxJetCandidate *  myJetCandidate ) const

private

Internal method to provide a check if the VxJetCandidate has been initialized in a consistent way.

Definition at line 556 of file JetFitterRoutines.cxx.

             {
      return false;
    }
    // end if else rimary==0
 
    // check std::vector<VxVertexOnJetAxis*> (if pointers are not empty and if
    // all associated tracks are not empty)
    const std::vector<VxVertexOnJetAxis*>& tracksOfVertex =
      myJetCandidate->getVerticesOnJetAxis();
 
    auto badVertex = [](VxVertexOnJetAxis* pVertex) {
      return (pVertex == nullptr) or (pVertex->getNumVertex() < 0);
    };
    ok =
      (std::find_if(tracksOfVertex.begin(), tracksOfVertex.end(), badVertex) ==
       tracksOfVertex.end());
    if (not ok) {
      ATH_MSG_WARNING(
        "One of the VxTrackAtVertex is a null pointer or uninitialized. Not "
        "proceeding with the fit!"); // Two error messages combined into one
    }
    if (not ok) {
      return false;
    }
    
    //now check if there is some track at least to do the fit...
 
    if (tracksOfVertex.empty()&&sizeprimary==0) {
      ATH_MSG_DEBUG( "No tracks at primary, no tracks on jet axis. Not proceeding with the fit!" );
      return false;
    }
    
    //now check if the number of tracks corrisponds to the number of components of the recVertexPositions
    const Trk::RecVertexPositions & myRecVertexPositions=myJetCandidate->getRecVertexPositions();
 
    const Amg::VectorX& myPosition=myRecVertexPositions.position();
 
    if (static_cast<unsigned int>(tracksOfVertex.size()+5)!=static_cast<unsigned int>(myPosition.rows())) { 
        ATH_MSG_WARNING ( "The position matrix has " <<  myPosition.rows() 
                          << " components while " <<  tracksOfVertex.size()+5 
                          << " are expected. Not proceeding with the fit " );
    }
 
    const Amg::MatrixX & myErrorMatrix=myRecVertexPositions.covariancePosition();
 
    //check if symMatrix and myPosition are compatible (size matches)
    if (myPosition.rows()!=myErrorMatrix.rows()) {
      ATH_MSG_WARNING ("The dimension of the position vector and the covariance matrix does not match. Not performing fit...");
      return false;
    }
 
    //check if all the diagonal values of the covariance matrix are not zero
    for (int i=0;i<myPosition.rows();i++) {
      if (std::abs(myErrorMatrix(i,i))<1e-20) {
        ATH_MSG_WARNING ("Value of cov matrix component n. " << i << " has a value smaller than 1e-8. Not considered as possible. Not performing fit...");
        return false;
      }
    }
 
    //well, more checks in the future!!!
    
    return true;
 
  }//end method
 
  
  std::pair<double,bool> JetFitterRoutines::fastProbabilityOfMerging(const VxVertexOnJetAxis* firstVertex,
                                     const VxVertexOnJetAxis* secondVertex,
                                     const VxJetCandidate* myJetCandidate) const {
    //this method is for evaluating the probability of merging in a very fast (and rough) way
    //if above a threshold you should think about evaluating the "fullProbabilityOfMerging" :-)
 
    const VxVertexOnJetAxis * PrimaryVertex=myJetCandidate->getPrimaryVertex();
 
    if (firstVertex==PrimaryVertex) {
      return fastProbabilityOfMergingWithPrimary(secondVertex,myJetCandidate);
    } 
    if (secondVertex==PrimaryVertex) {
      return fastProbabilityOfMergingWithPrimary(firstVertex,myJetCandidate);

copyRecoPositionsToLinearizationPositions()

void Trk::JetFitterRoutines::copyRecoPositionsToLinearizationPositions ( VxJetCandidate &  myJetCandidate ) const

private

Method to copy new reco positions to linearization positions after checking new positions are inside the ID.

Definition at line 1141 of file JetFitterRoutines.cxx.

        {
          if (absLinJetTheta>1e-8)
          {
            ATH_MSG_DEBUG (" Closest distance of track to jet axis is outside ID envelope, R=" << R << ", setting to R= " << m_maxR);
            distOnAxis=m_maxR / abssinLinJetTheta;
          }
        }
        
        Z=distOnAxis*cosLinJetTheta;
        if (std::abs(Z)>m_maxZ)
        {
          if (abscosLinJetTheta>1e-8)
          {
            ATH_MSG_DEBUG( " Closest distance of track to jet axis is outside ID envelope, Z=" << Z << ", setting to Z= " << m_maxZ );
            distOnAxis=m_maxZ / cosLinJetTheta;
          }
        }
        
        linPositions[numRow(myVertex->getNumVertex())]=distOnAxis;
      }
    }
    
    VertexPositions & linVertexPositions=myJetCandidate.getLinearizationVertexPositions();
    linVertexPositions.setPosition(linPositions);
    //now set the linearization position for the next step to the actual fitted vertex
    //OLD CODE BEFORE CHECKS
    //myJetCandidate->setLinearizationVertexPositions(myJetCandidate->getRecVertexPositions());
  }
  
}//end namespace

declareGaudiProperty() [1/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > &  hndl, const SG::VarHandleKeyArrayType &    )

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKeyArray>

Definition at line 170 of file AthCommonDataStore.h.

  {
    return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                       hndl.value(), 
                                                       hndl.documentation());
 
  }

declareGaudiProperty() [2/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > &  hndl, const SG::VarHandleKeyType &    )

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

  {
    return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                       hndl.value(), 
                                                       hndl.documentation());
 
  }

declareGaudiProperty() [3/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > &  hndl, const SG::VarHandleType &    )

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleBase>

Definition at line 184 of file AthCommonDataStore.h.

  {
    return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                       hndl.value(), 
                                                       hndl.documentation());
  }

declareGaudiProperty() [4/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > &  t, const SG::NotHandleType &    )

inlineprivateinherited

specialization for handling everything that's not a Gaudi::Property<SG::VarHandleKey> or a <SG::VarHandleKeyArray>

Definition at line 199 of file AthCommonDataStore.h.

  {
    return PBASE::declareProperty(t);
  }

declareProperty() [1/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name, SG::VarHandleBase &  hndl, const std::string &  doc, const SG::VarHandleType &    )

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
hndl	Object holding the property value.
doc	Documentation string for the property.

This is the version for types that derive from SG::VarHandleBase. The property value object is put on the input and output lists as appropriate; then we forward to the base class.

Definition at line 245 of file AthCommonDataStore.h.

  {
    this->declare(hndl.vhKey());
    hndl.vhKey().setOwner(this);
 
    return PBASE::declareProperty(name,hndl,doc);
  }

declareProperty() [2/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name, SG::VarHandleKey &  hndl, const std::string &  doc, const SG::VarHandleKeyType &    )

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
hndl	Object holding the property value.
doc	Documentation string for the property.

This is the version for types that derive from SG::VarHandleKey. The property value object is put on the input and output lists as appropriate; then we forward to the base class.

Definition at line 221 of file AthCommonDataStore.h.

  {
    this->declare(hndl);
    hndl.setOwner(this);
 
    return PBASE::declareProperty(name,hndl,doc);
  }

declareProperty() [3/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name, SG::VarHandleKeyArray &  hndArr, const std::string &  doc, const SG::VarHandleKeyArrayType &    )

inlineinherited

Definition at line 259 of file AthCommonDataStore.h.

  {
 
    // std::ostringstream ost;
    // ost << Algorithm::name() << " VHKA declareProp: " << name 
    //     << " size: " << hndArr.keys().size() 
    //     << " mode: " << hndArr.mode() 
    //     << "  vhka size: " << m_vhka.size()
    //     << "\n";
    // debug() << ost.str() << endmsg;
 
    hndArr.setOwner(this);
    m_vhka.push_back(&hndArr);
 
    Gaudi::Details::PropertyBase* p =  PBASE::declareProperty(name, hndArr, doc);
    if (p != 0) {
      p->declareUpdateHandler(&AthCommonDataStore<PBASE>::updateVHKA, this);
    } else {
      ATH_MSG_ERROR("unable to call declareProperty on VarHandleKeyArray " 
                    << name);
    }
 
    return p;
 
  }

declareProperty() [4/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name, T &  property, const std::string &  doc, const SG::NotHandleType &    )

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
property	Object holding the property value.
doc	Documentation string for the property.

This is the generic version, for types that do not derive from SG::VarHandleKey. It just forwards to the base class version of declareProperty.

Definition at line 333 of file AthCommonDataStore.h.

  {
    return PBASE::declareProperty(name, property, doc);
  }

declareProperty() [5/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( const std::string &  name, T &  property, const std::string &  doc = "none"  )

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
property	Object holding the property value.
doc	Documentation string for the property.

This dispatches to either the generic declareProperty or the one for VarHandle/Key/KeyArray.

Definition at line 352 of file AthCommonDataStore.h.

  {
    typedef typename SG::HandleClassifier<T>::type htype;
    return declareProperty (name, property, doc, htype());
  }

declareProperty() [6/6]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > &  t )

inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

                                                                     {
    typedef typename SG::HandleClassifier<T>::type htype;
    return AthCommonDataStore<PBASE>::declareGaudiProperty(t, htype());
  }

deleteVertexFromJetCandidate()

void Trk::JetFitterRoutines::deleteVertexFromJetCandidate	(	VxVertexOnJetAxis *	vertexToDelete,
		VxJetCandidate *	myJetCandidate
	)		const

Deltes a vertex from the VxJetCandidate, doing everything is needed to preserve the internal coherence of the VxJetCandidate object (e.g.

track numbering, list of tracks at vertices along the jet axis and so on)

Definition at line 1027 of file JetFitterRoutines.cxx.

                        {
 
      // RS 19.04.2011 attempt to fix coverity defect 22750
      std::vector<VxTrackAtVertex*>::iterator TracksBegin=tracksAtJetCandidate->begin();
      std::vector<VxTrackAtVertex*>::iterator TracksEnd=tracksAtJetCandidate->end();
      for (std::vector<VxTrackAtVertex*>::iterator TracksIter=TracksBegin;TracksIter!=TracksEnd;) {
    
    if (*TracksIter==*TracksAtVtxIter) {
      Trk::VxTrackAtVertex* oldPointer=*TracksIter;
      TracksIter=tracksAtJetCandidate->erase(TracksIter);
      delete oldPointer;
      oldPointer=nullptr;
      TracksEnd=tracksAtJetCandidate->end();
      break; 
    } 
      ++TracksIter;
    
      }
    }
 
    if (numberOfTracksBefore-numberOfTracksToDelete!=(int)myJetCandidate->vxTrackAtVertex()->size()) {
      ATH_MSG_DEBUG( " MISMATCH in JetFitterRoutines: the jetcandidate had: " << numberOfTracksBefore << " tracks " << 
    " and " << numberOfTracksToDelete << " to delete = " << myJetCandidate->vxTrackAtVertex()->size() << " tracks left! " );
    }
 
 
    //now as a last step you need to delete the vertex you're not using anymore...
    
    std::vector<VxVertexOnJetAxis*> copyOfVerticesAtJetCandidate=myJetCandidate->getVerticesOnJetAxis();
 
    const std::vector<VxVertexOnJetAxis*>::iterator VerticesBegin=copyOfVerticesAtJetCandidate.begin();
    std::vector<VxVertexOnJetAxis*>::iterator VerticesEnd=copyOfVerticesAtJetCandidate.end();
     
    bool found=false;
    
    //have a variable to store the vector of tracks which then need to be removed (copy)...
    //    std::vector<VxTrackAtVertex*> tracksToRemove=vertexToDelete->getTracksAtVertex();
 
    for (std::vector<VxVertexOnJetAxis*>::iterator VerticesIter=VerticesBegin;VerticesIter!=VerticesEnd;) {
      if ((*VerticesIter)==vertexToDelete) {
    delete *VerticesIter;
    VerticesIter=copyOfVerticesAtJetCandidate.erase(VerticesIter);
    VerticesEnd=copyOfVerticesAtJetCandidate.end();
    found=true;
    break;
      } 
    ++VerticesIter;
      
    }
    
    if (!found) {
      ATH_MSG_WARNING( "Could not find the vertex to delete... Very strange... Check!!! " );
    }
 
    //update myJetCandidate with the new vector of tracks
    myJetCandidate->setVerticesOnJetAxis(copyOfVerticesAtJetCandidate);
   
    //update the numbering scheme
    m_initializationHelper->updateTrackNumbering(myJetCandidate);
 
  }
 
  void JetFitterRoutines::copyRecoPositionsToLinearizationPositions(VxJetCandidate & myJetCandidate) const
  {
 
    const VertexPositions & newLinVertexPositions=myJetCandidate.getRecVertexPositions();
    Amg::VectorX linPositions=newLinVertexPositions.position();
    
    const std::vector<VxVertexOnJetAxis*> & associatedVertices=myJetCandidate.getVerticesOnJetAxis();
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxBegin=associatedVertices.begin();
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxEnd=associatedVertices.end();
    
    for (std::vector<VxVertexOnJetAxis*>::const_iterator VtxIter=VtxBegin;VtxIter!=VtxEnd;++VtxIter) {
      VxVertexOnJetAxis* myVertex=(*VtxIter);
      if (myVertex!=nullptr) {
        

detStore()

const ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const

inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

{ return m_detStore; }

evtStore() [1/2]

ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )

inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

{ return m_evtStore; }

evtStore() [2/2]

const ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( ) const

inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 90 of file AthCommonDataStore.h.

{ return m_evtStore; }

extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase &  ExtraDeps )

protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

fastProbabilityOfMerging()

std::pair< double, bool > Trk::JetFitterRoutines::fastProbabilityOfMerging	(	const VxVertexOnJetAxis *	firstVertex,
		const VxVertexOnJetAxis *	secondVertex,
		const VxJetCandidate *	myJetCandidate
	)		const

Calculates in a very fast way the probability two vertices along the jet axis to be compatible with each other, taking the constraint to be compatible with the previous fitted flight axis into account.

As a result a pair is provided:

• first member: probability of the two vertices to be compatible with a single vertex
• second vertex: if true, constraining the two vertices to be a single vertex moves the flight axis direction in DeltaR by more than the variable m_maxDRshift, which is set by the jobOption property maxDRshift. (the idea is that the bigger this shift is, the worse the approximation given by the previous linearization of these two vertices is and so the fast estimation of the compatibility between the two vertices).

If not needed the second member can just be ignored.

Description of the algorithmic procedure: a single Kalman Update step is performed on top of the previous fit, just constraining the two vertices to lie on the same position on the flight axis (this is done using a measurement equation with zero associated error). This is very fast.

Definition at line 652 of file JetFitterRoutines.cxx.

                                                                                    {
    
 
    //first get a copy of all vertex positions (this you can't avoid I fear...)
    RecVertexPositions copyOfRecVertexPositions(myJetCandidate->getRecVertexPositions());
 
 
    double oldchi2=copyOfRecVertexPositions.fitQuality().chiSquared();
    double oldndf=copyOfRecVertexPositions.fitQuality().numberDoF();
 
    //#ifdef JetFitterRoutines_DEBUG2
    const Amg::VectorX & positionVector=copyOfRecVertexPositions.position();

fastProbabilityOfMergingNoPrimary()

std::pair< double, bool > Trk::JetFitterRoutines::fastProbabilityOfMergingNoPrimary ( const VxVertexOnJetAxis *  firstVertex, const VxVertexOnJetAxis *  secondVertex, const VxJetCandidate *  myJetCandidate  ) const

private

Internal method to calculate fast probability of merging, for merging with non primary vertex.

Definition at line 717 of file JetFitterRoutines.cxx.

                       {
      ATH_MSG_WARNING ("In the compatibility estimation chi2: " << chi2 << " ndf " << ndf << " giving back 0 prob");
      return std::pair<double,bool>(0,isshifted);
    }
 
    return std::pair<double,bool>(TMath::Prob(chi2,(int)std::floor(ndf+0.5)),isshifted);
 
 
  }
 
  double JetFitterRoutines::fullProbabilityOfMerging(const VxVertexOnJetAxis* firstVertex,
                            const VxVertexOnJetAxis* secondVertex,
                            const VxJetCandidate* myJetCandidate,
                            int num_maxiterations,
                            bool treat_sign_flip,
                            int num_signflip_maxiterations,
                            double deltachi2_convergence) const {
 
    
    if (firstVertex==nullptr||secondVertex==nullptr||myJetCandidate==nullptr) {
      ATH_MSG_WARNING ("zero pointer given to the full probability estimation. No estimation performed, zero prob returned ");
      return 0;
    }
 

fastProbabilityOfMergingWithPrimary()

std::pair< double, bool > Trk::JetFitterRoutines::fastProbabilityOfMergingWithPrimary ( const VxVertexOnJetAxis *  otherVertex, const VxJetCandidate *  myJetCandidate  ) const

private

Internal method to calculate fast probability of merging, for merging with primary vertex.

Definition at line 672 of file JetFitterRoutines.cxx.

                       {
      ATH_MSG_WARNING (" In the compatibility estimation chi2: " << chi2 << " ndf " << ndf << " giving back 0 prob ");
      return std::pair<double,bool>(0,isshifted);
    }
 
    return std::pair<double,bool>(TMath::Prob(chi2,(int)std::floor(ndf+0.5)),isshifted);
 
  }
 
  std::pair<double,bool> JetFitterRoutines::fastProbabilityOfMergingNoPrimary(const VxVertexOnJetAxis* firstVertex,
                                          const VxVertexOnJetAxis* secondVertex,
                                          const VxJetCandidate* myJetCandidate) const {
    
    //first get a copy of all vertex positions (this you can't avoid I fear...)
    RecVertexPositions copyOfRecVertexPositions(myJetCandidate->getRecVertexPositions());
    const FitQuality & copyOfRecVertexQuality=copyOfRecVertexPositions.fitQuality();
 
    double oldchi2=copyOfRecVertexQuality.chiSquared();
    double oldndf=copyOfRecVertexQuality.numberDoF();
 
    //#ifdef JetFitterRoutines_DEBUG2
    const Amg::VectorX & positionVector=copyOfRecVertexPositions.position();
    const Amg::MatrixX & positionCov=copyOfRecVertexPositions.covariancePosition();

fillTableWithFastProbOfMerging()

void Trk::JetFitterRoutines::fillTableWithFastProbOfMerging

(

VxJetCandidate *

myJetCandidate

)

const

The VxClusteringTable of the VxJetCandidate provided in input is created, computing all the needed probabilities of merging between all the vertices two by two using the fast probability estimation (using for all combinations of vertices the fastProbabilityOfMerging method of the same class)

This is very fast and gives reasonably good results in all the cases where the fit is not too much affected by non-linearities.

Definition at line 864 of file JetFitterRoutines.cxx.

fillTableWithFullProbOfMerging()

void Trk::JetFitterRoutines::fillTableWithFullProbOfMerging	(	VxJetCandidate *	myJetCandidate,
		int	num_maxiterations = 20,
		bool	treat_sign_flip = true,
		int	num_signflip_maxiterations = 10,
		double	deltachi2_convergence = 1e-3,
		double	threshold_probability = 1e-5
	)		const

The VxClusteringTable of the VxJetCandidate provided in input is created, computing all the needed probabilities of merging between all the vertices two by two using a full probability estimation, with the complete fit performed with the options provided.

To increase the speed, before using the full estimation of merging probability a fast merging is done (fastProbabilityOfMerging method) and only if the fastprobability of merging is higher than the threshold_probability provided, the full probability is computed.

Definition at line 849 of file JetFitterRoutines.cxx.

                                                                                             {
    fillTableWithProbOfMerging(myJetCandidate,false);
  }
 
 
  void JetFitterRoutines::fillTableWithProbOfMerging(VxJetCandidate* myJetCandidate,
                             bool fullcomputation,
                             int num_maxiterations,
                             bool treat_sign_flip,
                             int num_signflip_maxiterations,
                             double deltachi2_convergence,
                             double threshold_probability) const {
    
    if (myJetCandidate==nullptr) {

fillTableWithProbOfMerging()

void Trk::JetFitterRoutines::fillTableWithProbOfMerging ( VxJetCandidate *  myJetCandidate, bool  fullcomputation, int  num_maxiterations = 20, bool  treat_sign_flip = true, int  num_signflip_maxiterations = 10, double  deltachi2_convergence = 1e-2, double  threshold_probability = 1e-3  ) const

private

Internal method to fill the VxClusteringTable of the VxJetCandidate object, independently on if the fullcomputation is requested or not.

Definition at line 869 of file JetFitterRoutines.cxx.

                                     {
      delete clusteringTablePtr;
    }
    clusteringTablePtr=new Trk::VxClusteringTable();
 
    double highestprobability(0.);
 
    VxVertexOnJetAxis* primaryVertex=myJetCandidate->getPrimaryVertex();
 
    if (primaryVertex==nullptr) {
    ATH_MSG_WARNING( "VxJetCandidate provided has no primary vertex. No compatibility table calculated." );
    return;
    } 
 
    primaryVertex->setCompatibilityToPrimaryVtx(1);//stupid but assign prob 1 to primary vtx for consistency
 
    const std::vector<VxVertexOnJetAxis*> & tracksOnJetAxis=myJetCandidate->getVerticesOnJetAxis();
 
    //now evaluate probability of cluster forming with primary vertex
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxBegin=tracksOnJetAxis.begin();
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxEnd=tracksOnJetAxis.end();
    
    for (std::vector<VxVertexOnJetAxis*>::const_iterator VtxIter=VtxBegin;
     VtxIter!=VtxEnd;++VtxIter) {
 
      std::pair<double,bool> fastProbabilityAndNonLinearity=fastProbabilityOfMerging(primaryVertex,
                                             *VtxIter,
                                             myJetCandidate);
      
      
      ATH_MSG_VERBOSE ("Fast probability of merging between primary and " << 
                   (*VtxIter)->getNumVertex() << " is " << fastProbabilityAndNonLinearity.first);
      
#ifdef JetFitterRoutines_DEBUG2  
      ATH_MSG_DEBUG( "Fast probability of merging between primary and " << 
    (*VtxIter)->getNumVertex() << " is " << fastProbabilityAndNonLinearity.first <<  " and is max DR " << 
    fastProbabilityAndNonLinearity.second );
#endif  
 
      if (fullcomputation) {
    if (fastProbabilityAndNonLinearity.first>threshold_probability) {
      if (fastProbabilityAndNonLinearity.first>highestprobability/100.&&fastProbabilityAndNonLinearity.second) {
        
        double fullProbability=fullProbabilityOfMerging(primaryVertex,*VtxIter,
                               myJetCandidate,num_maxiterations,
                               treat_sign_flip,
                               num_signflip_maxiterations,
                               deltachi2_convergence);
        
 
#ifdef JetFitterRoutines_DEBUG2     
        ATH_MSG_DEBUG( "Full probability of merging with primary is " << fullProbability );
#endif
 
        ATH_MSG_DEBUG ("Full probability of merging with primary is " << fullProbability);
        
        //store this fullProbability into the VxClusteringTable of myJetCandidate
        
        clusteringTablePtr->setCompatibilityOfTo(PairOfVxVertexOnJetAxis(primaryVertex,*VtxIter),fullProbability);
        if (fullProbability>highestprobability) {
          highestprobability=fullProbability;
        }
        (*VtxIter)->setCompatibilityToPrimaryVtx(fullProbability);
      } else {
        clusteringTablePtr->setCompatibilityOfTo(PairOfVxVertexOnJetAxis(primaryVertex,*VtxIter),fastProbabilityAndNonLinearity.first);
        (*VtxIter)->setCompatibilityToPrimaryVtx(fastProbabilityAndNonLinearity.first);     
      }
    } else {
      (*VtxIter)->setCompatibilityToPrimaryVtx(fastProbabilityAndNonLinearity.first);
    }
      } else {
    if (fastProbabilityAndNonLinearity.first>threshold_probability) {
      clusteringTablePtr->setCompatibilityOfTo(PairOfVxVertexOnJetAxis(primaryVertex,*VtxIter),fastProbabilityAndNonLinearity.first);
      if (fastProbabilityAndNonLinearity.first>highestprobability) {
        highestprobability=fastProbabilityAndNonLinearity.first;
      }
    }
    (*VtxIter)->setCompatibilityToPrimaryVtx(fastProbabilityAndNonLinearity.first);
      }
    }//end first iteration over tracks for compatibility with primary vertex
    
    //now all the remaining combination have to be tested...
 
    for (std::vector<Trk::VxVertexOnJetAxis*>::const_iterator VtxIter2=VtxBegin;
     VtxIter2!=VtxEnd;++VtxIter2) {
      for (std::vector<Trk::VxVertexOnJetAxis*>::const_iterator VtxIter1=VtxBegin;
       VtxIter1!=VtxIter2;++VtxIter1) {
    
    std::pair<double,bool> fastProbabilityAndNonLinearity=fastProbabilityOfMerging(*VtxIter1,
                                               *VtxIter2,
                                               myJetCandidate);
    
    
    
       ATH_MSG_VERBOSE ("Fast probability of merging between vtx n " <<  (*VtxIter1)->getNumVertex()<<
                        " and " << (*VtxIter2)->getNumVertex() << " is " << 
                        fastProbabilityAndNonLinearity.first << " and is max DR " <<
                    fastProbabilityAndNonLinearity.second);
 
#ifdef JetFitterRoutines_DEBUG2  
    ATH_MSG_DEBUG( "Fast probability of merging between vtx n " <<  (*VtxIter1)->getNumVertex() << " and " << 
      (*VtxIter2)->getNumVertex() << " is " << fastProbabilityAndNonLinearity.first << " and is max DR " << 
      fastProbabilityAndNonLinearity.second );
#endif  
    if (fullcomputation) {
      if (fastProbabilityAndNonLinearity.first>threshold_probability) {
        if (fastProbabilityAndNonLinearity.first>highestprobability/100.&&fastProbabilityAndNonLinearity.second) {
          double fullProbability=fullProbabilityOfMerging(*VtxIter1,*VtxIter2,
                                 myJetCandidate,num_maxiterations,
                                 treat_sign_flip,
                                 num_signflip_maxiterations,
                                 deltachi2_convergence);
          
 
#ifdef JetFitterRoutines_DEBUG2     
          ATH_MSG_DEBUG( "Full probability of merging is " << fullProbability );
#endif
 
          ATH_MSG_VERBOSE ("Full probability of merging is " << fullProbability);
          
          //store this fullProbability into the VxClusteringTable of myJetCandidate
          
          clusteringTablePtr->setCompatibilityOfTo(PairOfVxVertexOnJetAxis(*VtxIter1,*VtxIter2),fullProbability);
          if (fullProbability>highestprobability) {
        highestprobability=fullProbability;
          }
        } else {
          clusteringTablePtr->setCompatibilityOfTo(PairOfVxVertexOnJetAxis(*VtxIter1,*VtxIter2),fastProbabilityAndNonLinearity.first);
        }
      }
    } else {
      if (fastProbabilityAndNonLinearity.first>threshold_probability) {
        clusteringTablePtr->setCompatibilityOfTo(PairOfVxVertexOnJetAxis(*VtxIter1,*VtxIter2),fastProbabilityAndNonLinearity.first);
        if (fastProbabilityAndNonLinearity.first>highestprobability) {
          highestprobability=fastProbabilityAndNonLinearity.first;
        }
      }
    }
      }//end first vtx loop
    }//end second vtx loop
    
  }//end fillTableWithFullProbOfMerging() method
 
  void JetFitterRoutines::deleteVertexFromJetCandidate(VxVertexOnJetAxis* vertexToDelete,
                               VxJetCandidate* myJetCandidate) const {
    if (vertexToDelete==myJetCandidate->getPrimaryVertex()) {
      ATH_MSG_WARNING ("YOU ARE deleting the primary vertex. This is not possible... ");
      return;
    }
 
    // now you need to *delete* the second vertex in copyOfRecVertexPositions and 
    // in copyOfLinearizationPositions
    const Amg::VectorX & recPosition=myJetCandidate->getRecVertexPositions().position();
    const Amg::VectorX & linPosition=myJetCandidate->getLinearizationVertexPositions().position();
    const Amg::VectorX & constraintPosition=myJetCandidate->getConstraintVertexPositions().position();
    const Amg::MatrixX & covPosition=myJetCandidate->getRecVertexPositions().covariancePosition();
    const Amg::MatrixX & covConstraintPosition=myJetCandidate->getConstraintVertexPositions().covariancePosition();

fullProbabilityOfMerging()

double Trk::JetFitterRoutines::fullProbabilityOfMerging	(	const VxVertexOnJetAxis *	firstVertex,
		const VxVertexOnJetAxis *	secondVertex,
		const VxJetCandidate *	myJetCandidate,
		int	num_maxiterations = 20,
		bool	treat_sign_flip = true,
		int	num_signflip_maxiterations = 10,
		double	deltachi2_convergence = 1e-2
	)		const

Calculates in a complete way the probability two vertices along the jet axis to be compatible with each other, taking the constraint to be compatible with the previous fitted flight axis into account.

Contrary to the fast estimation method, here the complete fit is repeated and so the flight axis direction is corrected according to the new constraint taking all non linearities into account (because relinearization of the tracks partecipating in the fit is performed).

The following steps are done: 1) a temporary copy of the full VxJetCandidate has to be created 2) the copies of the vertices to join have to be found 3) the two copied vertices are merged in the new VxJetCandidate, using the mergeVerticesInJetCandidate function of the JetFitterHelper class 4) the track Number is updated, using the updateTrackNumbering function of the JetFitterInitializationHelper class 5) the full fit is performed on the new VxJetCandidate object

As a result the only full probability of merging is provided as a double.

Definition at line 765 of file JetFitterRoutines.cxx.

                                                           {
      ATH_MSG_WARNING ("Empty primary vertex found when estimating fullProbOfMerging. 0 prob returned.");
      return 0;
    }
 
    oldToNewVtxPointers[primaryOfFirst]=primaryOfSecond;
 
    unsigned int sizeOfVertices=vectorOfOldJetCand.size();
    if (vectorOfNewJetCand.size()!=sizeOfVertices) {
      ATH_MSG_WARNING ("Old and new track of vertices do not match during fullProbOfMerging. 0 prob returned.");
      return 0;
    }
 
    for (unsigned int s=0;s<sizeOfVertices;s++) {
      const VxVertexOnJetAxis* pointer1=vectorOfOldJetCand[s];
      VxVertexOnJetAxis* pointer2=vectorOfNewJetCand[s];
      if (pointer1==nullptr||pointer2==nullptr) {
    ATH_MSG_WARNING ("One of the pointers of the original or copied vector of vertices is empty during fullProbOfMerging. Skipping it...");
      } else {
    oldToNewVtxPointers[pointer1]=pointer2;
      }
    }
 
    //now merge firstVertex and secondVertex
    VxVertexOnJetAxis* newFirstVertex=oldToNewVtxPointers[firstVertex];
    VxVertexOnJetAxis* newSecondVertex=oldToNewVtxPointers[secondVertex];
 
    if (newFirstVertex==nullptr||newSecondVertex==nullptr) {
      ATH_MSG_WARNING ("No correspondence to the given firstVertex or secondVertex in fullProbOfMerging. Returning 0 prob.");
      return 0.;
    }
 
    VxVertexOnJetAxis & commonVertex=m_helper->mergeVerticesInJetCandidate(*newFirstVertex,
                                       *newSecondVertex,
                                       newJetCandidate);
 
    //now you need to update the numbering scheme
    m_initializationHelper->updateTrackNumbering(&newJetCandidate);
    
    performTheFit(&newJetCandidate,num_maxiterations,
          treat_sign_flip,
          num_signflip_maxiterations,
          deltachi2_convergence);
 
    //    const FitQuality & qualityOfMergedVertex=commonVertex.fitQuality();
 
#ifdef JetFitterRoutines_DEBUG
    ATH_MSG_DEBUG( " End estimating full prob of merging... chi2 " << commonVertex.fitQuality().chiSquared() << " ndf " << commonVertex.fitQuality().numberDoF() );
#endif
 
    return (double)TMath::Prob(commonVertex.fitQuality().chiSquared(),
                  (int)std::floor(commonVertex.fitQuality().numberDoF()+0.5));
 
  }
 
  void JetFitterRoutines::fillTableWithFullProbOfMerging(VxJetCandidate* myJetCandidate,
                             int num_maxiterations,
                             bool treat_sign_flip,
                             int num_signflip_maxiterations,
                             double deltachi2_convergence,
                             double threshold_probability) const {
    fillTableWithProbOfMerging(myJetCandidate,
                   true,
                   num_maxiterations,
                   treat_sign_flip,
                   num_signflip_maxiterations,
                   deltachi2_convergence,
                   threshold_probability);
  }

initialize()

StatusCode Trk::JetFitterRoutines::initialize ( )

overridevirtual

Definition at line 108 of file JetFitterRoutines.cxx.

                                                                                                {
    
    ATH_MSG_DEBUG ("initializingToMinDistancesToJetAxis is now implemented! Will converge faster!!! Neutrals are fully supported...");
 
    VertexPositions & linVertexPositions=myJetCandidate->getLinearizationVertexPositions();
    Amg::VectorX linPositions=linVertexPositions.position();
    
    
    const Trk::RecVertexPositions & recVtxPosition=myJetCandidate->getRecVertexPositions();
    const Amg::VectorX & recPosition=recVtxPosition.position();
    
    Amg::Vector3D a(recPosition[Trk::jet_xv],
                    recPosition[Trk::jet_yv],

initializeToMinDistancesToJetAxis()

void Trk::JetFitterRoutines::initializeToMinDistancesToJetAxis

(

VxJetCandidate *

myJetCandidate

)

const

This method provides the initialization of all the tracks in the fit to the position of minimum distance to the initial flight axis.

This should considerably improve the speed of convergence. It is still not implemented in the first version of the fit.

Definition at line 130 of file JetFitterRoutines.cxx.

                                    {//Was that your intention? to be checked... 15.03.2007
      SG::ReadCondHandle<AtlasFieldCacheCondObj> readHandle{m_fieldCacheCondObjInputKey};
      const AtlasFieldCacheCondObj* fieldCondObj{*readHandle};
      if (!readHandle.isValid()) {
         std::stringstream msg;
         msg << "Failed to retrieve magmnetic field conditions data " << m_fieldCacheCondObjInputKey.key() << ".";
         throw std::runtime_error(msg.str());
      }
      MagField::AtlasFieldCache fieldCache;
      fieldCondObj->getInitializedCache (fieldCache);
      for (std::vector<VxVertexOnJetAxis*>::const_iterator VtxIter=VtxBegin;VtxIter!=VtxEnd;++VtxIter) {
        VxVertexOnJetAxis* myVertex=(*VtxIter);
        if (myVertex!=nullptr) {
          
          const std::vector<VxTrackAtVertex*> & tracksAtVertex=myVertex->getTracksAtVertex();
          if (tracksAtVertex.size()>1) { 
              ATH_MSG_DEBUG( "Warning in JetFitterInitializationHelper.Number of tracks at vertex is bigger than one, " 
                             << "even during initialization phase. Skipping this vertex (already initialized)..." );
          } 
          else if (tracksAtVertex.empty())
          {
            ATH_MSG_WARNING( "No track at vertex. Internal fitter error. Contact author (GP) ... " );
          }
          else
          {
            ATH_MSG_VERBOSE( " TrackAtVertexSize is: " << tracksAtVertex.size() );
            ATH_MSG_VERBOSE( " pointer to first element: " << tracksAtVertex[0] );
            ATH_MSG_VERBOSE( " pointer to initialPerigee: " << tracksAtVertex[0]->initialPerigee() );
 
        // RS 19.04.2011 try to fix coverity defect 22750
        const Trk::Perigee* ptr = dynamic_cast<const Trk::Perigee*>((tracksAtVertex[0]->initialPerigee()));
            if (ptr) 
            {
 
              double distOnAxis=-999.;
              std::pair<Amg::Vector3D,double> result;
              try {
                result=m_minDistanceFinder->getPointAndDistance(myJetAxis,*ptr,distOnAxis, fieldCache);
 
                double R=distOnAxis*sinRecJetTheta;
                double Z=distOnAxis*cosRecJetTheta;
 
                if (std::abs(R)>m_maxR)
                {
 
                  if (absRecJetTheta>1e-8)
                  {
                    ATH_MSG_DEBUG( " Closest distance of track to jet axis is outside ID envelope, R=" << R << ", setting to R= " << m_maxR );
                    distOnAxis=m_maxR /abssinRecJetTheta;
                  }
                }
 
                Z=distOnAxis*cosRecJetTheta;
                if (std::abs(Z)>m_maxZ)
                {
                  if (abscosRecJetTheta>1e-8)
                  {
                    ATH_MSG_DEBUG( " Closest distance of track to jet axis is outside ID envelope, Z=" << Z << ", setting to Z= " << m_maxZ );
                    distOnAxis=m_maxZ / cosRecJetTheta;
                  }
                }
                
                linPositions[numRow(myVertex->getNumVertex())]=distOnAxis;
                  
 
              } catch (Error::NewtonProblem e) {
                ATH_MSG_WARNING( "Problem with Newton finder " << e.p );
              } catch (...) {
                ATH_MSG_ERROR( "Could not catch error " );
              }
              ATH_MSG_DEBUG ("initializingToMinDistancesToJetAxis for vertex number... " << 
                             myVertex->getNumVertex() << " to distance " << distOnAxis << 
                             " distance to axis " <<  result.second);
            }
            // FIXME   THIS PART IS DEAD IN NEW TRACKING EDM:
            else if (dynamic_cast<const Trk::NeutralPerigee*>((tracksAtVertex[0]->initialPerigee()))!=nullptr)
            {
              double distOnAxis=-999.;
              std::pair<Amg::Vector3D,double> result;
 
              const Trk::NeutralPerigee * neutralTrack = 
                dynamic_cast<const Trk::NeutralPerigee*>((tracksAtVertex[0]->initialPerigee()));
              
              NeutralTrack myNeutralTrack(neutralTrack->position(),neutralTrack->momentum());
              
              result=m_minDistanceFinderNeutral->getPointAndDistance(myJetAxis,myNeutralTrack,distOnAxis);
 
              double R=distOnAxis*sinRecJetTheta;
              double Z=distOnAxis*cosRecJetTheta;
 
              if (std::abs(R)>m_maxR)
              {
                
                if (absRecJetTheta>1e-8)
                {
                  ATH_MSG_DEBUG( " Closest distance of track to jet axis is outside ID envelope, R=" << R << ", setting to R= " << m_maxR );
                  distOnAxis=m_maxR / abssinRecJetTheta;
                }
              }
              
              Z=distOnAxis*cosRecJetTheta;
              if (std::abs(Z)>m_maxZ)
              {
                if (abscosRecJetTheta>1e-8)
                {
                  ATH_MSG_DEBUG( " Closest distance of track to jet axis is outside ID envelope, Z=" << Z << ", setting to Z= " << m_maxZ );
                  distOnAxis=m_maxZ / cosRecJetTheta;
                }
              }
              
              linPositions[numRow(myVertex->getNumVertex())]=distOnAxis; 
              ATH_MSG_DEBUG( "initializingToMinDistancesToJetAxis for vertex from NEUTRAL number... " << 
                               myVertex->getNumVertex() << " to distance " << 
                               distOnAxis << " distance to axis " <<  result.second );
              
            }
            else 
            {
              ATH_MSG_WARNING( "Could not cast to neither CHARGED or NEUTRAL! This error is not FATAL" );
            }
          }
        } else {
          ATH_MSG_WARNING( "Warning in JetFitterInitializationHelper.Inconsistency found. Pointer to VxVertexOnJetAxis should be different from zero. Skipping track..." );
          throw ("Warning in JetFitterInitializationHelper.Inconsistency found. Pointer to VxVertexOnJetAxis should be different from zero. Skipping track...");
        }
      }
 
      linVertexPositions.setPosition(linPositions);
      
    } else {
      ATH_MSG_DEBUG ("No Associated Vertices found! no initialization to minimum distance is possible.");
    }
  }
 
  void JetFitterRoutines::performTheFit(VxJetCandidate* myJetCandidate,
                    int num_maxiterations,//default 20
                    bool treat_sign_flip,//default true
                    int num_signflip_maxiterations,//default 10
                    double deltachi2_convergence) const {//default 0.001
    
    bool isInitialized=this->checkJetCandidate(myJetCandidate);
    
    if (!isInitialized) {
      ATH_MSG_DEBUG ("JetFitter found the candidate was not correctly initialized. Not proceeding with the fitt...");
      return;
    }
    
 

inputHandles()

virtual std::vector<Gaudi::DataHandle*> AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const

overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

interfaceID()

static const InterfaceID& Trk::JetFitterRoutines::interfaceID ( )

inlinestatic

Definition at line 61 of file JetFitterRoutines.h.

                                            {
      return IID_JetFitterRoutines;
    }

msg() [1/2]

MsgStream& AthCommonMsg< AlgTool >::msg ( ) const

inlineinherited

Definition at line 24 of file AthCommonMsg.h.

                                {
    return this->msgStream();
  }

msg() [2/2]

MsgStream& AthCommonMsg< AlgTool >::msg ( const MSG::Level  lvl ) const

inlineinherited

Definition at line 27 of file AthCommonMsg.h.

                                                  {
    return this->msgStream(lvl);
  }

msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level  lvl ) const

inlineinherited

Definition at line 30 of file AthCommonMsg.h.

                                               {
    return this->msgLevel(lvl);
  }

outputHandles()

virtual std::vector<Gaudi::DataHandle*> AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const

overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

performTheFit()

void Trk::JetFitterRoutines::performTheFit	(	VxJetCandidate *	myJetCandidate,
		int	num_maxiterations = 30,
		bool	treat_sign_flip = true,
		int	num_signflip_maxiterations = 30,
		double	deltachi2_convergence = 0.001
	)		const

This is the method where the fit is actually done.

The VxJetCandidate in input must be properly initialized (initial positions + proper initial covariance matrix for primary vertex, direction + all the tracks along the jet axis).

To initialized the VxJetCandidate properly please use the methods provided by the JetFitterInitializationHelper class.

The parameters to be provided are:

• num_maxiterations: the number of maximal iterations (in case of no convergence),
• treat_sign_flip: true if the special procedure for avoding flips to negative flightlengths has to be used
• deltachi2_convergence: if the difference in fit chi2 to the previous iteration is smaller than the value provided convergence is reached

Definition at line 300 of file JetFitterRoutines.cxx.

                         {
      do {
 
    myJetCandidate->setRecVertexPositions(RecVertexPositions(myJetCandidate->getConstraintVertexPositions()));
 
    ATH_MSG_DEBUG( myJetCandidate->getRecVertexPositions() );
    
    //linearization position is now stored in a different VertexPositions attached to the VxJetCandidate...
    m_initializationHelper->linearizeAllTracks(myJetCandidate,true);
 
    //now update all vertices with the new information
    updateAllVertices(myJetCandidate);
 
    //now do the chi2 update of all tracks after end of fit iterations
    //(this permits to invert the big weight matrix only once at the end!!)
    //  updateChi2AllVertices(myJetCandidate);
    
    num_iteration_signflip+=1;
    
    const RecVertexPositions & myRecPosition=myJetCandidate->getRecVertexPositions();
 
    const FitQuality & myFitQuality=myRecPosition.fitQuality();
 
    double actualchi2=myFitQuality.chiSquared();
 
    auto const absActLastChi2 = std::abs(actualchi2-lastchi2);
 
#ifdef JetFitterRoutines_DEBUG
    ATH_MSG_DEBUG( " last chi2 " << lastchi2 << " actual chi2 " << actualchi2 << " difference " << 
      absActLastChi2<< " < " << deltachi2_convergence << " ? " <<  " ndf " << myFitQuality.numberDoF() );
#endif
 
    if (absActLastChi2<deltachi2_convergence) {
      converged=true;
    } else {
 
          //GP 12.3.2012 Check that no vertex is outside the ID acceptance
          copyRecoPositionsToLinearizationPositions(*myJetCandidate);
 
    }
 
    lastchi2=actualchi2;
    
      } while ((!converged)&&num_iteration_signflip<num_signflip_maxiterations);
 
    
      if (converged) {
    ATH_MSG_VERBOSE( " For sign flip treatment there was convergence after " << num_iteration_signflip );
      } 
    }
    
    //    if (msgSvc()->outputLevel()== MSG::DEBUG||msgSvc()->outputLevel()== MSG::VERBOSE) {
#ifdef JetFitterRoutines_DEBUG
      ATH_MSG_DEBUG( "JetFitterRoutines: after convergence with sign flip treatment: " << myJetCandidate->getRecVertexPositions() );
#endif
      //    }
 
    int num_iteration=0;
    lastchi2=-99.;
    converged=false;
   
    do {
      
 
      myJetCandidate->setRecVertexPositions(RecVertexPositions(myJetCandidate->getConstraintVertexPositions()));
 
 
      //linearization position is now stored in a different VertexPositions attached to the VxJetCandidate...
      m_initializationHelper->linearizeAllTracks(myJetCandidate,false);
      
      //now update all vertices with the new information
      updateAllVertices(myJetCandidate);
 
      //now do the chi2 update of all tracks after end of fit iterations
      //(this permits to invert the big weight matrix only once at the end!!)
      //      updateChi2AllVertices(myJetCandidate);
 
      num_iteration+=1;
      
      const RecVertexPositions & myRecPosition=myJetCandidate->getRecVertexPositions();
 
      ATH_MSG_VERBOSE( " num_iteration (full fit): " << num_iteration << " det " << myRecPosition.covariancePosition().determinant() <<  " recVertex " << myJetCandidate->getRecVertexPositions() );
 
      const FitQuality & myFitQuality=myRecPosition.fitQuality();
      double actualchi2=myFitQuality.chiSquared(); 
      
      auto const absActLastChi2 = std::abs(actualchi2-lastchi2);      
 
      //    if (msgSvc()->outputLevel()== MSG::VERBOSE) {
#ifdef JetFitterRoutines_DEBUG
      ATH_MSG_DEBUG( " without sign flip: last chi2 " << lastchi2 << " actual chi2 " << actualchi2 << " difference " << 
    absActLastChi2 << " < " << deltachi2_convergence << " ? " <<  " ndf " << myFitQuality.numberDoF() );
#endif
      //    }
 
      if (absActLastChi2<deltachi2_convergence) {
    converged=true;
      } else {
    //now set the linearization position for the next step to the actual fitted vertex
          //GP 12.3.2012 Check that no vertex is outside the ID acceptance
          copyRecoPositionsToLinearizationPositions(*myJetCandidate);
          //OLD myJetCandidate->setLinearizationVertexPositions(myJetCandidate->getRecVertexPositions());
      }
      lastchi2=actualchi2;
      
    } while ((!converged)&&num_iteration<num_maxiterations);
    
    if (converged) {
      ATH_MSG_VERBOSE( " Fit without sign flip treatment there was convergence after " << num_iteration );
    } 
    
    //    if (msgSvc()->outputLevel()== MSG::VERBOSE) {
    #ifdef JetFitterRoutines_DEBUG
    ATH_MSG_DEBUG( "JetFitterRoutines: after convergence without sign flip treatment: " << myJetCandidate->getRecVertexPositions() );
    #endif
      //    }
    if (num_iteration>=num_maxiterations)
    {
      ATH_MSG_DEBUG( "There wasn't convergence in JetFitter after: " << num_maxiterations );
    }
 
    //now only the smoothing is missing as a last step... (updated momenta, chi2 + ndf of clusters,...)
    
    smoothAllVertices(myJetCandidate);
    
    Trk::VxJetFitterDebugInfo * & myDebugInfo=myJetCandidate->getDebugInfo();
    if (myDebugInfo!=nullptr) {
      delete myDebugInfo;
    }
    myDebugInfo=new VxJetFitterDebugInfo();
    myDebugInfo->setNumFitIterations(num_iteration);
    myDebugInfo->setSignFlipNumFitIterations(num_iteration_signflip);
 
  }
 
 
  void JetFitterRoutines::updateAllVertices(VxJetCandidate* myJetCandidate) const {
 
//    ATH_MSG_DEBUG( " Updating PV " );
 
    int n_iteration=0;
    
 
    if (!m_noPrimaryVertexRefit) {
      //new iteration
      VxVertexOnJetAxis* myPrimary=myJetCandidate->getPrimaryVertex();
      const std::vector<VxTrackAtVertex*> & primaryVectorTracks=myPrimary->getTracksAtVertex();
      
      const std::vector<VxTrackAtVertex*>::const_iterator primaryVectorTracksBegin=primaryVectorTracks.begin();

renounce()

std::enable_if_t<std::is_void_v<std::result_of_t<decltype(&T::renounce)(T)> > && !std::is_base_of_v<SG::VarHandleKeyArray, T> && std::is_base_of_v<Gaudi::DataHandle, T>, void> AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T &  h )

inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

  {
    h.renounce();
    PBASE::renounce (h);
  }

renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray &  handlesArray )

inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

                                                          {
    handlesArray.renounce();
  }

smoothAllVertices()

void Trk::JetFitterRoutines::smoothAllVertices

(

VxJetCandidate *

myJetCandidate

)

const

triggers the smoothing of all vertices (the tracks in the fit are updated using the constraint provided by the vertices along the jet axis they were fitted to)

Definition at line 529 of file JetFitterRoutines.cxx.

                                                                                                 {
      if (!m_fast) {
    m_smoother->update(*VtxIter,myJetCandidate);
      } else {
    m_smoother->fastUpdate(*VtxIter,myJetCandidate);
      }
    }
  }
 
 
 
 
  bool JetFitterRoutines::checkJetCandidate(VxJetCandidate* myJetCandidate) const {
    
    int sizeprimary(0);
 
    //check if primaryvertex exists
    const VxVertexOnJetAxis* myPrimary=myJetCandidate->getPrimaryVertex();
    if (myPrimary==nullptr) {
      ATH_MSG_WARNING( "No primary vertex found in VxJetCandidate class. Initialization was not done correctly..." );
      return false;
    }
 

sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )

overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in DerivationFramework::CfAthAlgTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and asg::AsgMetadataTool.

sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )

overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

updateAllVertices()

void Trk::JetFitterRoutines::updateAllVertices

(

VxJetCandidate *

myJetCandidate

)

const

One iteration of the Kalman Updated of all tracks to the actual fit is performed.

Difference between fitting to primary vertex or to a vertex along the flight axis is taken into account.

Definition at line 457 of file JetFitterRoutines.cxx.

                                                                      {
    
    ATH_MSG_VERBOSE( " RecVertexPositions before update " << myJetCandidate->getRecVertexPositions() );
    
    if ((!m_fast)) {
      m_updator->add(*primaryVectorIter,myPrimary,myJetCandidate);
    } else {
      m_updator->addWithFastUpdate(*primaryVectorIter,myPrimary,myJetCandidate);
    }
 
        const RecVertexPositions & myRecPosition=myJetCandidate->getRecVertexPositions();
 
        ATH_MSG_VERBOSE( " Determinant after PRIMARY VTX update: " << n_iteration <<  " det: "  << myRecPosition.covariancePosition().determinant() <<  " recVertex " << myJetCandidate->getRecVertexPositions() );
      }
    }
 
    n_iteration=0;
    
    // const RecVertexPositions & myRecPositionBeg=myJetCandidate->getRecVertexPositions();
 
    const std::vector<VxVertexOnJetAxis*> & associatedVertices=myJetCandidate->getVerticesOnJetAxis();
    
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxBegin=associatedVertices.begin();
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxEnd=associatedVertices.end();
    
    for (std::vector<VxVertexOnJetAxis*>::const_iterator VtxIter=VtxBegin;VtxIter!=VtxEnd;++VtxIter) {
 
//      ATH_MSG_DEBUG( " Updating an SV along jet axis " );
 
      const std::vector<VxTrackAtVertex*> & tracksAtVertex=(*VtxIter)->getTracksAtVertex();
 
      const std::vector<VxTrackAtVertex*>::const_iterator TracksBegin=tracksAtVertex.begin();
      const std::vector<VxTrackAtVertex*>::const_iterator TracksEnd=tracksAtVertex.end();
      
      for (std::vector<VxTrackAtVertex*>::const_iterator  TrackVectorIter=TracksBegin;
       TrackVectorIter!=TracksEnd;++TrackVectorIter) {
    
    ATH_MSG_VERBOSE( " RecVertexPositions before update " << myJetCandidate->getRecVertexPositions() );
    
    if (!m_fast) {
      m_updator->add(*TrackVectorIter,*VtxIter,myJetCandidate);
    } else {
      m_updator->addWithFastUpdate(*TrackVectorIter,*VtxIter,myJetCandidate);
    }
 
        n_iteration+=1;
 
        const RecVertexPositions & myRecPosition=myJetCandidate->getRecVertexPositions();
 
        ATH_MSG_VERBOSE( " Determinant after sec update: " << n_iteration <<  " det: "  << myRecPosition.covariancePosition().determinant() <<  " recVertex " << myJetCandidate->getRecVertexPositions() );
      }
    }
 
     myJetCandidate->getRecVertexPositions().finalizePosition();
  }
 
 void JetFitterRoutines::smoothAllVertices(VxJetCandidate* myJetCandidate) const {
    //new iteration
    VxVertexOnJetAxis* myPrimary=myJetCandidate->getPrimaryVertex();
 
    if (!m_fast) {
      m_smoother->update(myPrimary,myJetCandidate);    
    } else {
      m_smoother->fastUpdate(myPrimary,myJetCandidate); 
    }
 
    const std::vector<VxVertexOnJetAxis*> & associatedVertices=myJetCandidate->getVerticesOnJetAxis();
    
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxBegin=associatedVertices.begin();
    const std::vector<VxVertexOnJetAxis*>::const_iterator VtxEnd=associatedVertices.end();

updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase &  )

inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

                                               {
    // debug() << "updateVHKA for property " << p.name() << " " << p.toString() 
    //         << "  size: " << m_vhka.size() << endmsg;
    for (auto &a : m_vhka) {
      std::vector<SG::VarHandleKey*> keys = a->keys();
      for (auto k : keys) {
        k->setOwner(this);
      }
    }
  }

Member Data Documentation

m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore

privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore

privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

m_fast

bool Trk::JetFitterRoutines::m_fast

private

Definition at line 285 of file JetFitterRoutines.h.

m_fieldCacheCondObjInputKey

SG::ReadCondHandleKey<AtlasFieldCacheCondObj> Trk::JetFitterRoutines::m_fieldCacheCondObjInputKey {this, "AtlasFieldCacheCondObj", "fieldCondObj", "Name of the Magnetic Field conditions object key"}

private

Definition at line 269 of file JetFitterRoutines.h.

m_helper

ToolHandle<JetFitterHelper> Trk::JetFitterRoutines::m_helper

private

Definition at line 273 of file JetFitterRoutines.h.

m_initializationHelper

ToolHandle<JetFitterInitializationHelper> Trk::JetFitterRoutines::m_initializationHelper

private

Definition at line 272 of file JetFitterRoutines.h.

m_maxDRshift

double Trk::JetFitterRoutines::m_maxDRshift

private

Definition at line 294 of file JetFitterRoutines.h.

m_maxR

double Trk::JetFitterRoutines::m_maxR

private

Definition at line 317 of file JetFitterRoutines.h.

m_maxZ

double Trk::JetFitterRoutines::m_maxZ

private

Definition at line 318 of file JetFitterRoutines.h.

m_minDistanceFinder

ToolHandle<TrkDistanceFinderNeutralCharged> Trk::JetFitterRoutines::m_minDistanceFinder

private

Definition at line 276 of file JetFitterRoutines.h.

m_minDistanceFinderNeutral

ToolHandle<TrkDistanceFinderNeutralNeutral> Trk::JetFitterRoutines::m_minDistanceFinderNeutral

private

Definition at line 277 of file JetFitterRoutines.h.

m_noPrimaryVertexRefit

bool Trk::JetFitterRoutines::m_noPrimaryVertexRefit

private

Definition at line 311 of file JetFitterRoutines.h.

m_smoother

ToolHandle<KalmanVertexOnJetAxisSmoother> Trk::JetFitterRoutines::m_smoother

private

Definition at line 275 of file JetFitterRoutines.h.

m_updator

ToolHandle<KalmanVertexOnJetAxisUpdator> Trk::JetFitterRoutines::m_updator

private

Definition at line 274 of file JetFitterRoutines.h.

m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared

privateinherited

Definition at line 399 of file AthCommonDataStore.h.

m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka

privateinherited

Definition at line 398 of file AthCommonDataStore.h.

---

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

• JetFitterRoutines.h
• JetFitterRoutines.cxx