From 2c2063fc39e986aab1c6fcef800da8d7b9cc6982 Mon Sep 17 00:00:00 2001
From: Christian Gumpert <Christian.Gumpert@cern.ch>
Date: Mon, 7 Dec 2015 18:51:59 +0100
Subject: [PATCH] fix compiler warnings, workaround const_correctness in
material update (TrkExEngine-00-00-38)
* fix warnings in StaticNavigationEngine
* use const_cast to update parameters with material effects
* tag as TrkExEngine-00-00-38
2015-11-04 Noemi Calace <noemi.calace@cern.ch>
* Adding direction in the energy loss calculation
* tag TrkExEngine-00-00-37
2015-10-26 Noemi Calace <noemi.calace@cern.ch>
* Updating configuration of the AtlasExtrapolationEngine
* Making changes for propation to the endlayer
* tag TrkExEngine-00-00-36
2015-10-16 Sharka Todorova <sarka.todorova@cern.ch>
* first import of StepEngine
* updates in ExtrapolationEngine/StaticNavigationEngine
* tag TrkExEngine-00-00-35
2015-08-04 Shaun Roe
* fix coverity defects
...
(Long ChangeLog diff - truncated)
Former-commit-id: f8f9e563d870a8e68aa0b061f9ec773151046e4f
---
.../TrkExEngine/CMakeLists.txt | 38 --
.../TrkExEngine/ExtrapolationEngine.h | 3 +
.../TrkExEngine/ExtrapolationEngine.icc | 34 +-
.../TrkExEngine/TrkExEngine/StaticEngine.icc | 4 +-
.../TrkExEngine/StaticNavigationEngine.h | 18 +-
.../TrkExEngine/StaticNavigationEngine.icc | 35 +-
.../TrkExEngine/TrkExEngine/StepEngine.h | 134 ++++++++
.../TrkExEngine/TrkExEngine/StepEngine.icc | 144 ++++++++
.../python/AtlasExtrapolationEngine.py | 2 +
.../TrkExEngine/src/ExtrapolationEngine.cxx | 8 +
.../TrkExEngine/src/MaterialEffectsEngine.cxx | 10 +-
.../TrkExEngine/src/PropagationEngine.cxx | 23 +-
.../src/StaticNavigationEngine.cxx | 8 +
.../TrkExEngine/src/StepEngine.cxx | 324 ++++++++++++++++++
.../src/components/TrkExEngine_entries.cxx | 3 +
15 files changed, 710 insertions(+), 78 deletions(-)
delete mode 100644 Tracking/TrkExtrapolation/TrkExEngine/CMakeLists.txt
create mode 100644 Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.h
create mode 100644 Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.icc
create mode 100644 Tracking/TrkExtrapolation/TrkExEngine/src/StepEngine.cxx
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/CMakeLists.txt b/Tracking/TrkExtrapolation/TrkExEngine/CMakeLists.txt
deleted file mode 100644
index 054d1e7e788..00000000000
--- a/Tracking/TrkExtrapolation/TrkExEngine/CMakeLists.txt
+++ /dev/null
@@ -1,38 +0,0 @@
-################################################################################
-# Package: TrkExEngine
-################################################################################
-
-# Declare the package name:
-atlas_subdir( TrkExEngine )
-
-# Declare the package's dependencies:
-atlas_depends_on_subdirs( PUBLIC
- Control/AthenaBaseComps
- DetectorDescription/GeoPrimitives
- GaudiKernel
- Tracking/TrkDetDescr/TrkDetDescrInterfaces
- Tracking/TrkDetDescr/TrkGeometry
- Tracking/TrkDetDescr/TrkSurfaces
- Tracking/TrkDetDescr/TrkVolumes
- Tracking/TrkEvent/TrkEventPrimitives
- Tracking/TrkEvent/TrkNeutralParameters
- Tracking/TrkEvent/TrkParameters
- Tracking/TrkExtrapolation/TrkExInterfaces
- Tracking/TrkExtrapolation/TrkExUtils
- PRIVATE
- Tracking/TrkDetDescr/TrkDetDescrUtils )
-
-# External dependencies:
-find_package( Eigen )
-
-# Component(s) in the package:
-atlas_add_component( TrkExEngine
- src/*.cxx
- src/components/*.cxx
- INCLUDE_DIRS ${EIGEN_INCLUDE_DIRS}
- LINK_LIBRARIES ${EIGEN_LIBRARIES} AthenaBaseComps GeoPrimitives GaudiKernel TrkDetDescrInterfaces TrkGeometry TrkSurfaces TrkVolumes TrkEventPrimitives TrkNeutralParameters TrkParameters TrkExInterfaces TrkExUtils TrkDetDescrUtils )
-
-# Install files from the package:
-atlas_install_headers( TrkExEngine )
-atlas_install_python_modules( python/*.py )
-
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.h b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.h
index 26e6b13c3dd..6fd4d265fad 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.h
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.h
@@ -27,6 +27,7 @@ namespace Trk {
class TrackingGeometry;
class IPropagationEngine;
+ class INavigationEngine;
/** @class ExtrapolationEngine
@@ -102,6 +103,8 @@ namespace Trk {
ToolHandle<IPropagationEngine> m_propagationEngine; //!< the used propagation engine for navigation initialization
std::vector<const IExtrapolationEngine*> m_eeAccessor; //!< the extrapolation engines for
+ ToolHandle<INavigationEngine> m_navigationEngine; //!< access to tracking geometry (unique?)
+
//!< forces a global search for the initialization, allows to switch TrackingGeometries in one job
bool m_forceSearchInit;
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.icc b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.icc
index ded0f6ed7cc..874d77c0fbc 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.icc
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/ExtrapolationEngine.icc
@@ -12,6 +12,7 @@
#include "TrkGeometry/TrackingVolume.h"
#include "TrkGeometry/Layer.h"
#include "TrkExInterfaces/IPropagationEngine.h"
+#include "TrkExInterfaces/INavigationEngine.h"
#include <iostream>
#include <iomanip>
@@ -30,10 +31,12 @@ template <class T> Trk::ExtrapolationCode Trk::ExtrapolationEngine::extrapolateT
// give output that you are in the master volume loop
EX_MSG_VERBOSE(eCell.navigationStep, "extrapolate", "loop", "processing volume : " << eCell.leadVolume->volumeName() );
// get the appropriate IExtrapolationEngine
- const Trk::IExtrapolationEngine* iee = m_eeAccessor[eCell.leadVolume->geometryType()];
+ Trk::GeometryType geoType = eCell.leadVolume->geometrySignature()>2 ? Trk::Dense : Trk::Static;
+ const Trk::IExtrapolationEngine* iee = m_eeAccessor[geoType];
eCode = iee ? iee->extrapolate(eCell, sf, bcheck) : Trk::ExtrapolationCode::FailureConfiguration;
// give a message about what you have
- EX_MSG_VERBOSE(eCell.navigationStep, "extrapolate", "", "returned from volume with return code : " << eCode.toString() );
+ EX_MSG_VERBOSE(eCell.navigationStep, "extrapolate", "", "returned from volume with return code : " << eCode.toString() <<
+ " and geoType:"<< geoType );
}
EX_MSG_DEBUG(eCell.navigationStep, "extrapolate", "", "extrapolation finished with return code : " << eCode.toString() );
// before you return, finalize: sets the leadParameters to endParameters and empties the garbage bin
@@ -65,26 +68,23 @@ template <class T> Trk::ExtrapolationCode Trk::ExtrapolationEngine::initNavigati
// initialize the start parameters - try association first
eCell.startLayer = eCell.startLayer ? eCell.startLayer : eCell.leadParameters->associatedSurface().associatedLayer();
eCell.startVolume = eCell.startVolume ? eCell.startVolume :
- ( eCell.startLayer ? eCell.startLayer->enclosingTrackingVolume() : m_trackingGeometry->lowestTrackingVolume(eCell.leadParameters->position()) );
+ ( eCell.startLayer ? eCell.startLayer->enclosingTrackingVolume() : 0 );
+
+ if (!eCell.startVolume || m_trackingGeometry->atVolumeBoundary(eCell.startParameters->position(),eCell.startVolume,0.001) ) {
+ Trk::ExtrapolationCode eVol = m_navigationEngine->resolvePosition(eCell,dir,true);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ eCell.startVolume = eCell.leadVolume;
+ } else {
+ eCell.leadVolume = eCell.startVolume;
+ }
// bail out of the start volume can not be resolved
if (!eCell.startVolume) return Trk::ExtrapolationCode::FailureNavigation;
// screen output
EX_MSG_VERBOSE( eCell.navigationStep, "navigation", "", "start volume termined as : " << eCell.startVolume->volumeName() );
- // check if the parameters are on a volume boundary
- // - @TODO identify this case by something like eCell.leadParameters->associatedSurface().isBoundary()
- if ( eCell.startVolume->onVolumeBoundary(*eCell.startParameters) ){
- // re-evaluate the volume by stepping out of the volume
- EX_MSG_VERBOSE( eCell.navigationStep, "navigation", "", "parameters are on volume boundary, stepping out of this volume." );
- // stepping by one unit out of the volume
- eCell.startVolume = m_trackingGeometry->lowestTrackingVolume(eCell.leadParameters->position()+dir*eCell.leadParameters->momentum().unit());
- // screen output
- EX_MSG_VERBOSE( eCell.navigationStep, "navigation", "", "start volume re-evaluated as : " << eCell.startVolume->volumeName());
- }
+ // check layer association
eCell.startLayer = eCell.startLayer ? eCell.startLayer : eCell.startVolume->associatedLayer(eCell.leadParameters->position());
if (eCell.startLayer)
EX_MSG_VERBOSE( eCell.navigationStep, "navigation", "", "start layer termined with index : " << eCell.startLayer->layerIndex());
- // now you can assign the lead volume
- eCell.leadVolume = eCell.startVolume;
eCell.setRadialDirection();
// ---------- END initialization -----------------------------------------------------------------------------------------
if (sf){
@@ -95,7 +95,7 @@ template <class T> Trk::ExtrapolationCode Trk::ExtrapolationEngine::initNavigati
if ( !eCell.checkConfigurationMode(Trk::ExtrapolationMode::FATRAS) )
eCell.setRadialDirection();
// trying association via the layer : associated layer of material layer
- eCell.endLayer = sf ? ( sf->associatedLayer() ? sf->associatedLayer() : sf->materialLayer() ) : 0;
+ eCell.endLayer = sf->associatedLayer() ? sf->associatedLayer() : sf->materialLayer() ;
eCell.endVolume = eCell.endLayer ? eCell.endLayer->enclosingTrackingVolume() : 0;
// check if you found layer and volume
if (!eCell.endVolume){
@@ -105,7 +105,7 @@ template <class T> Trk::ExtrapolationCode Trk::ExtrapolationEngine::initNavigati
Trk::ExtrapolationCell<T> navCell(*eCell.leadParameters, dir);
navCell.addConfigurationMode(Trk::ExtrapolationMode::Direct);
// screen output
- Trk::ExtrapolationCode eCode = m_propagationEngine->propagate(navCell,*eCell.endSurface,Trk::anyDirection,true,eCell.navigationCurvilinear);
+ Trk::ExtrapolationCode eCode = m_propagationEngine->propagate(navCell,*eCell.endSurface,Trk::anyDirection,false,eCell.navigationCurvilinear);
// check for sucess to the destination
CHECK_ECODE_SUCCESS_NODEST(navCell, eCode);
// screen output
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticEngine.icc b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticEngine.icc
index 81540dca147..ac07bbc396b 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticEngine.icc
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticEngine.icc
@@ -316,7 +316,7 @@ template <class T> Trk::ExtrapolationCode Trk::StaticEngine::resolveLayerT(Trk::
// now loop over the surfaces:
// the surfaces will be sorted
for (auto& csf : cSurfaces ) {
- EX_MSG_VERBOSE(eCell.navigationStep, "layer", eCell.leadLayer->layerIndex().value(), "trying candidate surfaces with straight line path length " << csf.intersection.pathLength);
+ EX_MSG_VERBOSE(eCell.navigationStep, "layer", eCell.leadLayer->layerIndex().value(), "trying candidate surfaces with straight line path length " << csf.intersection.pathLength);
// propagate to the compatible surface, return types are
// - InProgress : propagation to compatible surface worked
// - Recovered : propagation to compatible surface did not work, leadParameters stay the same
@@ -364,7 +364,7 @@ template <class T> Trk::ExtrapolationCode Trk::StaticEngine::resolveLayerT(Trk::
// Possible return types:
// - SuccessDestination : great, desintation surface hit - but post-update needs to be done
// - SuccessPathLimit : pathlimit was reached on the way to the destination surface
- eCode = m_propagationEngine->propagate(eCell,*sf,pDir,bcheck,eCell.destinationCurvilinear);
+ eCode = m_propagationEngine->propagate(eCell,*sf,pDir,false,eCell.destinationCurvilinear);
// check for success return path limit
CHECK_ECODE_SUCCESS_NODEST(eCell,eCode);
EX_MSG_VERBOSE(eCell.navigationStep, "layer", eCell.leadLayer->layerIndex().value(), "attempt to hit destination surface resulted in " << eCode.toString() );
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.h b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.h
index dd31e5f0b5d..2470d3b66b2 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.h
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.h
@@ -53,18 +53,33 @@ namespace Trk {
/** avoid method shaddowing */
using INavigationEngine::resolveBoundary;
+ using INavigationEngine::resolvePosition;
/** resolve the boundary situation - for charged particles */
virtual ExtrapolationCode resolveBoundary(Trk::ExCellCharged& eCell, PropDirection dir=alongMomentum) const;
/** resolve the boundary situation - for neutral particles */
virtual ExtrapolationCode resolveBoundary(Trk::ExCellNeutral& eCelll, PropDirection dir=alongMomentum) const;
+
+ /** resolve the boundary situation - for charged particles */
+ virtual ExtrapolationCode resolvePosition(Trk::ExCellCharged& eCell, PropDirection dir=alongMomentum, bool noLoop=false) const;
+
+ /** resolve the boundary situation - for neutral particles */
+ virtual ExtrapolationCode resolvePosition(Trk::ExCellNeutral& eCelll, PropDirection dir=alongMomentum, bool noLoop=false) const;
+
+ /** acces to tracking geometry */
+ virtual const TrackingGeometry& trackingGeometry() const throw (GaudiException);
private:
/** resolve the boundary situation */
template <class T> ExtrapolationCode resolveBoundaryT(ExtrapolationCell<T>& eCell,
PropDirection dir=alongMomentum) const;
+ /** resolve position */
+ template <class T> ExtrapolationCode resolvePositionT(ExtrapolationCell<T>& eCell,
+ PropDirection dir=alongMomentum,
+ bool noLoop=false) const;
+
/** deal with the boundary Surface - called by resolveBoundary */
template <class T> ExtrapolationCode handleBoundaryT(ExtrapolationCell<T>& eCell,
const BoundarySurface<TrackingVolume>& bSurfaceTV,
@@ -75,9 +90,6 @@ namespace Trk {
//!< retrieve TrackingGeometry
StatusCode updateTrackingGeometry() const;
- //!< return and retrieve
- const TrackingGeometry& trackingGeometry() const throw (GaudiException);
-
ToolHandle<IPropagationEngine> m_propagationEngine; //!< the used propagation engine
ToolHandle<IMaterialEffectsEngine> m_materialEffectsEngine; //!< the material effects updated
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.icc b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.icc
index ca6d15345c2..89a1f305061 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.icc
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StaticNavigationEngine.icc
@@ -180,14 +180,15 @@ template <class T> Trk::ExtrapolationCode Trk::StaticNavigationEngine::handleBou
if (stopAtThisBoundary){
EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "", "geometry signature change from " << eCell.leadVolume->geometrySignature() << " to " << nextVolume->geometrySignature());
eCell.nextGeometrySignature = nextVolume->geometrySignature();
- // return the boundary reached
+ // return the boundary reached : the navigation resolved already
+ eCell.leadVolume = nextVolume;
return Trk::ExtrapolationCode::SuccessBoundaryReached;
}
// remember the last boundary surface for loop protection
eCell.lastBoundarySurface = &bSurface;
eCell.lastBoundaryParameters = eCell.leadParameters;
// set next volume and reset lead layer
- eCell.leadVolume = nextVolume;
+ eCell.leadVolume = nextVolume;
eCell.leadLayer = 0;
// we have bParameters -> break the loop over boundaryIntersections
return Trk::ExtrapolationCode::InProgress;
@@ -196,3 +197,33 @@ template <class T> Trk::ExtrapolationCode Trk::StaticNavigationEngine::handleBou
// you need to keep on trying
return Trk::ExtrapolationCode::Unset;
}
+
+
+/** handle the failure - as configured */
+template <class T> Trk::ExtrapolationCode Trk::StaticNavigationEngine::resolvePositionT(Trk::ExtrapolationCell<T>& eCell,
+ Trk::PropDirection pDir,
+ bool /*noLoop*/ ) const
+{
+ EX_MSG_DEBUG(++eCell.navigationStep, "navigation", "", "resolve position '"<< eCell.leadParameters->position()
+ << (int(pDir) > 0 ? "' along momentum." : "' opposite momentum.") );
+
+ // noLoop= True is used when we have exit from leadVolume
+
+ if (!eCell.leadVolume) eCell.leadVolume = trackingGeometry().lowestStaticTrackingVolume(eCell.leadParameters->position());
+ if (!eCell.leadVolume) return Trk::ExtrapolationCode::FailureNavigation;
+ const Trk::TrackingVolume* nextVol=0;
+ if ( trackingGeometry().atVolumeBoundary(eCell.leadParameters->position(),
+ eCell.leadParameters->momentum(),
+ eCell.leadVolume,
+ nextVol, pDir, 0.01) ) { // set tolerance globally
+
+ //if (noLoop && nextVol==eCell.leadVolume) return Trk::ExtrapolationCode::FailureLoop;
+
+ if (nextVol) {
+ eCell.leadVolume = nextVol;
+ return Trk::ExtrapolationCode::InProgress;
+ } else return Trk::ExtrapolationCode::FailureNavigation;
+ }
+
+ return Trk::ExtrapolationCode::InProgress;
+}
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.h b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.h
new file mode 100644
index 00000000000..bd16cab0e74
--- /dev/null
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.h
@@ -0,0 +1,134 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// StepEngine.h, (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+
+#ifndef TRKEXENINGE_STEPENGINE_H
+#define TRKEXENINGE_STEPENGINE_H
+
+#ifndef TRKEXENINGE_OUTPUTHELPER
+#define TRKEXENINGE_OUTPUTHELPER
+#define OH_CHECKFOUND(object) ( object ? "found" : "not found")
+#endif
+
+// Gaudi
+#include "AthenaBaseComps/AthAlgTool.h"
+#include "GaudiKernel/ToolHandle.h"
+// Trk
+#include "TrkExInterfaces/IExtrapolationEngine.h"
+#include "TrkExInterfaces/ExtrapolationMacros.h"
+#include "TrkExInterfaces/IPropagator.h"
+#include "TrkExUtils/ExtrapolationCell.h"
+#include "TrkExUtils/TargetSurfaces.h"
+#include "TrkParameters/TrackParameters.h"
+#include "TrkNeutralParameters/NeutralParameters.h"
+
+namespace Trk {
+
+ class IMaterialEffectsEngine;
+ class INavigationEngine;
+
+ /** @class StepEngine
+
+ Extrapolation engine for arbitrary tracking geometry environment. Intended as primary choice for Calo/MS.
+
+ @author sarka.todorova -at- cern.ch
+
+ */
+ class StepEngine : public AthAlgTool, virtual public IExtrapolationEngine {
+
+ public:
+
+ /** @enum ResolveLayerType
+ - use for code readability
+ */
+ enum ResolveLayerType {
+ StartLayer = 0,
+ NavigationLayer = 1,
+ PassThroughLayer = 2,
+ SubStructureLayer = 3,
+ DestinationLayer = 4,
+ StartAndDestinationLayer = 6,
+ UndefinedLayer = 5
+ };
+
+ /** Constructor */
+ StepEngine(const std::string&,const std::string&,const IInterface*);
+
+ /** Destructor */
+ ~StepEngine();
+
+ /** AlgTool initialize method */
+ StatusCode initialize();
+
+ /** AlgTool finalize method */
+ StatusCode finalize();
+
+ using IExtrapolationEngine::extrapolate;
+
+ /** charged extrapolation - public interface */
+ virtual ExtrapolationCode extrapolate(ExCellCharged& ecCharged,
+ const Surface* sf = 0,
+ BoundaryCheck bcheck = true) const;
+
+
+ /** neutral extrapolation - public interface */
+ virtual ExtrapolationCode extrapolate(ExCellNeutral& ecNeutral,
+ const Surface* sf = 0,
+ BoundaryCheck bcheck = true) const;
+
+ /** define for which GeometrySignature this extrapolator is valid - this is GLOBAL */
+ GeometryType geometryType() const;
+
+ private:
+ /** main loop extrapolation method */
+ template <class T> Trk::ExtrapolationCode targetSurfacesT(ExtrapolationCell<T>& eCell,
+ Trk::TargetSurfaceVector& ts,
+ bool trueOrderedIntersections,
+ const Surface* sf = 0,
+ BoundaryCheck bcheck = true ) const;
+
+ template<class T> Trk::ExtrapolationCode resolveFrameBoundaryT(ExtrapolationCell<T>& eCell,
+ Amg::Vector3D position,
+ unsigned int index ) const;
+ /** distance calculations */
+ void evaluateDistance(Trk::TargetSurface& tt, Amg::Vector3D pos, Amg::Vector3D mom,
+ Trk::TargetSurfaceVector&ts, bool trueOrdered) const;
+
+ /** handle extrapolation step */
+ Trk::ExtrapolationCode handleIntersection(ExCellCharged& ecCharged,
+ Trk::TargetSurfaceVector& solutions) const;
+
+ ToolHandle<IPropagator> m_propagator; //!< the used propagation engine
+ ToolHandle<IMaterialEffectsEngine> m_materialEffectsEngine; //!< the material effects updated
+ ToolHandle<INavigationEngine> m_navigationEngine; //!< access to tracking geometry
+
+ // local variables
+ double m_tolerance;
+
+ // intersection cache
+ mutable std::vector<Amg::Vector3D> m_intersections;
+
+ // target surfaces
+ mutable TargetSurfaces m_targetSurfaces;
+
+ // debugging mode & switch
+ bool m_debugAndFix; // switch
+ mutable bool m_debugCall; // mode
+
+ };
+
+ inline GeometryType StepEngine::geometryType() const
+ { return Trk::Dense; }
+
+
+} // end of namespace
+
+//!< define the templated function
+#include "StepEngine.icc"
+
+#endif // TRKEXINTERFACES_IStepEngine_H
+
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.icc b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.icc
new file mode 100644
index 00000000000..66f0003e50f
--- /dev/null
+++ b/Tracking/TrkExtrapolation/TrkExEngine/TrkExEngine/StepEngine.icc
@@ -0,0 +1,144 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// StepEngine.icc, (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+
+#include "TrkExInterfaces/IMaterialEffectsEngine.h"
+#include "TrkExInterfaces/INavigationEngine.h"
+#include "TrkSurfaces/Surface.h"
+#include "TrkGeometry/TrackingGeometry.h"
+#include "TrkGeometry/TrackingVolume.h"
+#include "TrkGeometry/Layer.h"
+#include <iostream>
+#include <iomanip>
+
+template<class T> Trk::ExtrapolationCode Trk::StepEngine::targetSurfacesT(Trk::ExtrapolationCell<T>& eCell,
+ Trk::TargetSurfaceVector& ts,
+ bool trueOrderedIntersections,
+ const Trk::Surface* sf,
+ Trk::BoundaryCheck bcheck ) const
+{
+ // evaluate the initial distance and save for the propagator/transport step
+ Amg::Vector3D gp = eCell.leadParameters->position();
+ Amg::Vector3D mom = eCell.leadParameters->momentum().normalized();
+ Trk::PropDirection dir = eCell.propDirection;
+
+ // destination surface
+ if (sf) {
+ Trk::TargetSurface tt(sf,bcheck,Trk::SurfNavigType::Target,0,0,Trk::TVNavigType::Unknown);
+ evaluateDistance(tt,gp,mom*dir,ts,trueOrderedIntersections);
+ // abort if no valid intersection
+ if (!ts.size()) return Trk::ExtrapolationCode::FailureDestination;
+ }
+
+ // static frame boundaries
+ const std::vector< Trk::SharedObject< const Trk::BoundarySurface< Trk::TrackingVolume> > > bounds = eCell.leadVolume->boundarySurfaces();
+ // count number of non-trivial solutions
+ unsigned int iInt=0;
+ for (unsigned int ib=0; ib< bounds.size(); ib++ ){
+ const Trk::Surface& surf = (bounds[ib].getPtr())->surfaceRepresentation();
+ Trk::TargetSurface bb(&surf,true,Trk::SurfNavigType::BoundaryFrame,ib,eCell.leadVolume,Trk::TVNavigType::Frame);
+ evaluateDistance(bb,gp,mom*dir,ts,trueOrderedIntersections);
+ if (ts.size() && ts.back().surf==&surf && ts.back().distanceAlongPath>m_tolerance) iInt++;
+ }
+ if (iInt==0 ) { // navigation problem : frame volume not resolved : trigger return && catch inf loop
+ return Trk::ExtrapolationCode::FailureLoop;
+ }
+
+ // frame resolved
+
+ // order intersections
+ if ( ts.size()>1 && trueOrderedIntersections ) {
+
+ std::vector<unsigned int> sols(ts.size());
+ for (unsigned int i=0;i<ts.size(); i++) { sols[i]=i; }
+
+ unsigned int itest=1;
+ while ( itest<sols.size() ) {
+ if ( ts[sols[itest]].distanceAlongPath < ts[sols[itest-1]].distanceAlongPath ) {
+ unsigned int iex = sols[itest-1];
+ sols[itest-1]=sols[itest];
+ sols[itest]=iex;
+ itest=1;
+ } else itest++;
+ }
+
+ Trk::TargetSurfaceVector tsOrd ;
+ for (unsigned int i=0;i<ts.size(); i++) tsOrd.push_back(ts[sols[i]]);
+ for (unsigned int i=0;i<ts.size(); i++) ts[i]=tsOrd[i];
+ }
+
+ return Trk::ExtrapolationCode::InProgress;
+}
+
+template<class T> Trk::ExtrapolationCode Trk::StepEngine::resolveFrameBoundaryT(Trk::ExtrapolationCell<T>& eCell, Amg::Vector3D gp,
+ unsigned int bIndex ) const
+{
+ const std::vector< Trk::SharedObject< const Trk::BoundarySurface< Trk::TrackingVolume> > > bounds = eCell.leadVolume->boundarySurfaces();
+
+ EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "","current lead volume:"<<eCell.leadVolume->volumeName());
+
+ const Trk::TrackingVolume* nextVolume = bounds[bIndex].getPtr()->attachedVolume(gp,
+ eCell.leadParameters->momentum(),
+ eCell.propDirection);
+
+ if (nextVolume) EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "", "resolveFrameBoundary:"<< nextVolume->volumeName()
+ <<" at position:"<< eCell.leadParameters->position().perp() <<","<<eCell.leadParameters->position().z()
+ << ","<<eCell.leadParameters->momentum().normalized()<<","<<eCell.propDirection);
+
+ if (!nextVolume) return Trk::ExtrapolationCode::SuccessBoundaryReached;
+
+ // - geometrySignature change and configuration to stop then triggers a Success
+ bool stopAtThisBoundary = eCell.checkConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary)
+ && (nextVolume->geometrySignature() != eCell.leadVolume->geometrySignature());
+
+ if (stopAtThisBoundary) {
+ EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "", "geometry signature change from " <<
+ eCell.leadVolume->geometrySignature() << " to " << nextVolume->geometrySignature());
+ eCell.nextGeometrySignature = nextVolume->geometrySignature();
+ // return the boundary reached
+ return Trk::ExtrapolationCode::SuccessBoundaryReached;
+ } else {
+ // fill the boundary into the cache if successfully hit boundary surface
+ eCell.stepParameters(eCell.leadParameters, Trk::ExtrapolationMode::CollectBoundary);
+ }
+ // loop protection - relaxed for the cases where you start from the boundary
+ if (eCell.leadVolume == nextVolume ) {
+ // the start parameters where on the boundary already give a relaxed return code
+ // if (&bSurface == eCell.lastBoundarySurface) return Trk::ExtrapolationCode::Unset;
+ // give some screen output as of why this happens
+ EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "",
+ "loop detected while trying to leave TrackingVolume '" << nextVolume->volumeName() << ".");
+
+ // can be a precision problem
+ if ((eCell.leadParameters->position()-eCell.lastBoundaryParameters->position()).mag()>0
+ && nextVolume->inside(gp+0.01*eCell.leadParameters->momentum().normalized()*eCell.propDirection,0.001) )
+ EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "",
+ "recovered, continue propagation in " << nextVolume->volumeName() << ".");
+ else return Trk::ExtrapolationCode::FailureLoop;
+ }
+ // remember the last boundary surface for loop protection
+ // ecNeutral.lastBoundarySurface = &bounds[index];
+ eCell.lastBoundaryParameters = eCell.leadParameters;
+ // set next volume
+ eCell.leadVolume = nextVolume;
+
+ // retrieve frame boundary surfaces, check early exit
+ if ( !m_targetSurfaces.initFrameVolume(eCell.leadParameters->position(),
+ eCell.leadParameters->momentum().normalized()*eCell.propDirection,
+ eCell.leadVolume ) ) {
+ // problem detected in the update of tracking surfaces: probably exit from the frame volume within tolerance range
+ EX_MSG_VERBOSE(eCell.navigationStep, "navigation", "","early exit detected at boundary index:"<< m_targetSurfaces.nextSf());
+ Amg::Vector3D gp = eCell.leadParameters->position()+
+ eCell.leadParameters->momentum().normalized()*eCell.propDirection*m_targetSurfaces.distanceToNext();
+ return resolveFrameBoundaryT(eCell,gp,m_targetSurfaces.nextSf());
+
+ }
+
+ EX_MSG_VERBOSE(eCell.navigationStep, "navigation", ""," frame volume resolved :"<< eCell.leadVolume->volumeName());
+
+ return Trk::ExtrapolationCode::InProgress;
+ }
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/python/AtlasExtrapolationEngine.py b/Tracking/TrkExtrapolation/TrkExEngine/python/AtlasExtrapolationEngine.py
index d537cdf342a..c976e4207d7 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/python/AtlasExtrapolationEngine.py
+++ b/Tracking/TrkExtrapolation/TrkExEngine/python/AtlasExtrapolationEngine.py
@@ -58,6 +58,7 @@ class AtlasExtrapolationEngine( ExEngine ):
# configure output formatting
MaterialEffectsEngine.OutputPrefix = '[ME] - '
MaterialEffectsEngine.OutputPostfix = ' - '
+ #MaterialEffectsEngine.EnergyLossCorrection = False
if ToolOutputLevel :
MaterialEffectsEngine.OutputLevel = ToolOutputLevel
# add to tool service
@@ -99,6 +100,7 @@ class AtlasExtrapolationEngine( ExEngine ):
ExEngine.__init__(self, name=nameprefix+'Extrapolation',\
ExtrapolationEngines = [ StaticExtrapolator ], \
PropagationEngine = StaticPropagator, \
+ NavigationEngine = StaticNavigator, \
TrackingGeometrySvc = AtlasTrackingGeometrySvc, \
OutputPrefix = '[ME] - ', \
OutputPostfix = ' - ')
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/src/ExtrapolationEngine.cxx b/Tracking/TrkExtrapolation/TrkExEngine/src/ExtrapolationEngine.cxx
index 6fdfb5eed6f..af6415dcb15 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/src/ExtrapolationEngine.cxx
+++ b/Tracking/TrkExtrapolation/TrkExEngine/src/ExtrapolationEngine.cxx
@@ -19,6 +19,7 @@ Trk::ExtrapolationEngine::ExtrapolationEngine(const std::string& t, const std::s
m_trackingGeometrySvc("AtlasTrackingGeometrySvc", n),
m_trackingGeometryName("AtlasTrackingGeometry"),
m_propagationEngine(""),
+ m_navigationEngine(""),
m_forceSearchInit(false)
{
declareInterface<Trk::IExtrapolationEngine>(this);
@@ -28,6 +29,7 @@ Trk::ExtrapolationEngine::ExtrapolationEngine(const std::string& t, const std::s
declareProperty("ExtrapolationEngines" , m_extrapolationEngines);
// The Tools needed
declareProperty("PropagationEngine" , m_propagationEngine);
+ declareProperty("NavigationEngine" , m_navigationEngine);
// steering of the screen outoput (SOP)
declareProperty("OutputPrefix" , m_sopPrefix);
declareProperty("OutputPostfix" , m_sopPostfix);
@@ -73,6 +75,12 @@ StatusCode Trk::ExtrapolationEngine::initialize()
} else
EX_MSG_DEBUG("", "initialize", "", "successfully propagation engine '" << m_propagationEngine << "'." );
+ if (m_navigationEngine.retrieve().isFailure()){
+ EX_MSG_FATAL("", "initialize", "", "failed to retrieve navigation engine '"<< m_navigationEngine << "'. Aborting." );
+ return StatusCode::FAILURE;
+ } else
+ EX_MSG_DEBUG("", "initialize", "", "successfully retrieved '" << m_navigationEngine << "'." );
+
return StatusCode::SUCCESS;
}
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/src/MaterialEffectsEngine.cxx b/Tracking/TrkExtrapolation/TrkExEngine/src/MaterialEffectsEngine.cxx
index d6d8a56387d..95a5e2d9461 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/src/MaterialEffectsEngine.cxx
+++ b/Tracking/TrkExtrapolation/TrkExEngine/src/MaterialEffectsEngine.cxx
@@ -155,7 +155,7 @@ const Trk::TrackParameters* Trk::MaterialEffectsEngine::updateTrackParameters(co
double sigmaP = 0.;
double kazl = 0.;
/** dE/dl ionization energy loss per path unit */
- double dEdl = sign*m_interactionFormulae.dEdl_ionization(p, &material, eCell.pHypothesis, sigmaP, kazl);
+ double dEdl = sign*dir*m_interactionFormulae.dEdl_ionization(p, &material, eCell.pHypothesis, sigmaP, kazl);
double dE = thickness*pathCorrection*dEdl;
sigmaP *= thickness*pathCorrection;
// calcuate the new momentum
@@ -166,7 +166,7 @@ const Trk::TrackParameters* Trk::MaterialEffectsEngine::updateTrackParameters(co
// update the covariance if needed
if (uCovariance)
(*uCovariance)(Trk::qOverP, Trk::qOverP) += sign*sigmaQoverP*sigmaQoverP;
- }
+ }
// (B) - update the covariance if needed
if (uCovariance && m_mscCorrection){
/** multiple scattering as function of dInX0 */
@@ -189,7 +189,7 @@ const Trk::TrackParameters* Trk::MaterialEffectsEngine::updateTrackParameters(co
// now either create new ones or update - only start parameters can not be updated
if (eCell.leadParameters != eCell.startParameters ){
EX_MSG_VERBOSE(eCell.navigationStep, "layer", layer->layerIndex().value(), "material update on non-initial parameters.");
- parameters.updateParameters(uParameters,uCovariance);
+ const_cast<Trk::TrackParameters*>(¶meters)->updateParameters(uParameters,uCovariance);
} else {
EX_MSG_VERBOSE(eCell.navigationStep, "layer", layer->layerIndex().value(), "material update on initial parameters, creating new ones.");
// create new parameters
@@ -200,8 +200,8 @@ const Trk::TrackParameters* Trk::MaterialEffectsEngine::updateTrackParameters(co
uParameters[Trk::theta],
uParameters[Trk::qOverP],
uCovariance);
- // these are newly created
- return tParameters;
+ // these are newly created
+ return tParameters;
}
}
return (¶meters);
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/src/PropagationEngine.cxx b/Tracking/TrkExtrapolation/TrkExEngine/src/PropagationEngine.cxx
index 2fee90054e2..8e5f8305513 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/src/PropagationEngine.cxx
+++ b/Tracking/TrkExtrapolation/TrkExEngine/src/PropagationEngine.cxx
@@ -12,6 +12,8 @@
#include "TrkExEngine/PropagationEngine.h"
#include "TrkExInterfaces/IPropagator.h"
#include "TrkSurfaces/Surface.h"
+#include "TrkEventPrimitives/PropDirection.h"
+//https://svnweb.cern.ch/trac/atlasoff/browser/Tracking/TrkEvent/TrkEventPrimitives/trunk/TrkEventPrimitives/PropDirection.h
// constructor
Trk::PropagationEngine::PropagationEngine(const std::string& t, const std::string& n, const IInterface* p)
@@ -37,7 +39,6 @@ Trk::PropagationEngine::~PropagationEngine()
// the interface method initialize
StatusCode Trk::PropagationEngine::initialize()
{
-
if (m_propagator.retrieve().isFailure()){
EX_MSG_FATAL( "", "initialize", "", "failed to retrieve propagator '"<< m_propagator << "'. Aborting." );
return StatusCode::FAILURE;
@@ -130,7 +131,7 @@ Trk::ExtrapolationCode Trk::PropagationEngine::propagate(Trk::ExCellNeutral& eCe
// it is the final propagation if it is the endSurface
bool finalPropagation = (eCell.endSurface == (&sf));
// intersect the surface
- Trk::Intersection sfIntersection = pDir ? sf.straightLineIntersection(eCell.leadParameters->position(),
+ Trk::Intersection sfIntersection = (pDir!=Trk::anyDirection) ? sf.straightLineIntersection(eCell.leadParameters->position(),
pDir*eCell.leadParameters->momentum().unit(),
true, bcheck) :
sf.straightLineIntersection(eCell.leadParameters->position(),
@@ -138,17 +139,17 @@ Trk::ExtrapolationCode Trk::PropagationEngine::propagate(Trk::ExCellNeutral& eCe
false, bcheck);
// we have a valid intersection
if (sfIntersection.valid){
- // fill the transport information - only if the propation direction is not 0
- if (pDir){
+ // fill the transport information - only if the propation direction is not 0 ('anyDirection')
+ if (pDir!=Trk::anyDirection){
double pLength = (sfIntersection.position-eCell.leadParameters->position()).mag();
EX_MSG_VERBOSE(eCell.navigationStep,"propagate", "neut", "path length of " << pLength << " added to the extrapolation cell (limit = " << eCell.pathLimit << ")" );
eCell.stepTransport(sf,pLength);
}
// now check if it is valud it's further away than the pathLimit
if (eCell.pathLimitReached(m_pathLimitTolerance)){
- // cache the last lead parameters
- eCell.lastLeadParameters = eCell.leadParameters;
- // create new neutral curvilinear parameters at the path limit reached
+ // cache the last lead parameters
+ eCell.lastLeadParameters = eCell.leadParameters;
+ // create new neutral curvilinear parameters at the path limit reached
double pDiff = eCell.pathLimit - cPath;
eCell.leadParameters = new Trk::NeutralCurvilinearParameters(eCell.leadParameters->position()+pDiff*eCell.leadParameters->momentum().unit(),
eCell.leadParameters->momentum(),
@@ -156,8 +157,8 @@ Trk::ExtrapolationCode Trk::PropagationEngine::propagate(Trk::ExCellNeutral& eCe
EX_MSG_VERBOSE(eCell.navigationStep,"propagate", "neut", "path limit of " << eCell.pathLimit << " reached. Stopping extrapolation.");
return Trk::ExtrapolationCode::SuccessPathLimit;
}
- // cache the last lead parameters
- eCell.lastLeadParameters = eCell.leadParameters;
+ // cache the last lead parameters
+ eCell.lastLeadParameters = eCell.leadParameters;
// now exchange the lead parameters
// create the new curvilinear paramters at the surface intersection -> if so, trigger the success
eCell.leadParameters = returnCurvilinear ? new Trk::NeutralCurvilinearParameters(sfIntersection.position, eCell.leadParameters->momentum(), 0.) :
@@ -165,9 +166,9 @@ Trk::ExtrapolationCode Trk::PropagationEngine::propagate(Trk::ExCellNeutral& eCe
// check if the propagation was called with directly, then lead parameters become end parameters
if (eCell.checkConfigurationMode(Trk::ExtrapolationMode::Direct))
- eCell.endParameters = eCell.leadParameters;
+ eCell.endParameters = eCell.leadParameters;
- // return success for the final destination or in progress
+ // return success for the final destination or in progress
return (finalPropagation ? Trk::ExtrapolationCode::SuccessDestination : Trk::ExtrapolationCode::InProgress);
} else {
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/src/StaticNavigationEngine.cxx b/Tracking/TrkExtrapolation/TrkExEngine/src/StaticNavigationEngine.cxx
index a18d05c646e..85e66703fe8 100644
--- a/Tracking/TrkExtrapolation/TrkExEngine/src/StaticNavigationEngine.cxx
+++ b/Tracking/TrkExtrapolation/TrkExEngine/src/StaticNavigationEngine.cxx
@@ -69,6 +69,14 @@ Trk::ExtrapolationCode Trk::StaticNavigationEngine::resolveBoundary(Trk::ExCellC
Trk::ExtrapolationCode Trk::StaticNavigationEngine::resolveBoundary(Trk::ExCellNeutral& ecNeutral, PropDirection dir) const
{ return resolveBoundaryT<Trk::NeutralParameters>(ecNeutral,dir); }
+/** charged */
+Trk::ExtrapolationCode Trk::StaticNavigationEngine::resolvePosition(Trk::ExCellCharged& ecCharged, PropDirection dir, bool noLoop) const
+{ return resolvePositionT<Trk::TrackParameters>(ecCharged,dir, noLoop); }
+
+/** neutral */
+Trk::ExtrapolationCode Trk::StaticNavigationEngine::resolvePosition(Trk::ExCellNeutral& ecNeutral, PropDirection dir, bool noLoop) const
+{ return resolvePositionT<Trk::NeutralParameters>(ecNeutral,dir, noLoop); }
+
StatusCode Trk::StaticNavigationEngine::updateTrackingGeometry() const {
// retrieve the TrackingGeometry from the detector store
if (detStore()->retrieve(m_trackingGeometry, m_trackingGeometryName).isFailure()){
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/src/StepEngine.cxx b/Tracking/TrkExtrapolation/TrkExEngine/src/StepEngine.cxx
new file mode 100644
index 00000000000..85ad067e7ce
--- /dev/null
+++ b/Tracking/TrkExtrapolation/TrkExEngine/src/StepEngine.cxx
@@ -0,0 +1,324 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+///////////////////////////////////////////////////////////////////
+// StaticEngine.cxx, (c) ATLAS Detector software
+///////////////////////////////////////////////////////////////////
+
+// STL
+#include <sstream>
+// Trk include
+#include "TrkExEngine/StepEngine.h"
+
+// constructor
+Trk::StepEngine::StepEngine(const std::string& t, const std::string& n, const IInterface* p)
+: AthAlgTool(t,n,p),
+ m_propagator(""),
+ m_materialEffectsEngine(""),
+ m_navigationEngine(""),
+ m_tolerance(0.01),
+ m_debugAndFix(true),
+ m_debugCall(false)
+{
+ declareInterface<Trk::IExtrapolationEngine>(this);
+ // The Tools needed
+ declareProperty("Propagator" , m_propagator);
+ declareProperty("MaterialEffectsEngine" , m_materialEffectsEngine);
+ declareProperty("NavigationEngine" , m_navigationEngine);
+ // steering of the screen outoput (SOP)
+ declareProperty("OutputPrefix" , m_sopPrefix);
+ declareProperty("OutputPostfix" , m_sopPostfix);
+ declareProperty("DebugAndFixMode" , m_debugAndFix);
+}
+
+// destructor
+Trk::StepEngine::~StepEngine()
+{}
+
+
+// the interface method initialize
+StatusCode Trk::StepEngine::initialize()
+{
+ if (m_propagator.retrieve().isFailure()){
+ EX_MSG_FATAL("", "initialize", "", "failed to retrieve propagator '"<< m_propagator << "'. Aborting." );
+ return StatusCode::FAILURE;
+ } else
+ EX_MSG_DEBUG("", "initialize", "", "successfully retrieved '" << m_propagator << "'." );
+
+ if (m_materialEffectsEngine.retrieve().isFailure()){
+ EX_MSG_FATAL("", "initialize", "", "failed to retrieve material effect engine '"<< m_materialEffectsEngine << "'. Aborting." );
+ return StatusCode::FAILURE;
+ } else
+ EX_MSG_DEBUG("", "initialize", "", "successfully retrieved '" << m_materialEffectsEngine << "'." );
+
+ if (m_navigationEngine.retrieve().isFailure()){
+ EX_MSG_FATAL("", "initialize", "", "failed to retrieve navigation engine '"<< m_navigationEngine << "'. Aborting." );
+ return StatusCode::FAILURE;
+ } else
+ EX_MSG_DEBUG("", "initialize", "", "successfully retrieved '" << m_navigationEngine << "'." );
+
+ EX_MSG_DEBUG("", "initialize", "", "successful." );
+
+ m_targetSurfaces.setDebugModeOff();
+
+ return StatusCode::SUCCESS;
+}
+
+// the interface method finalize
+StatusCode Trk::StepEngine::finalize()
+{
+ EX_MSG_DEBUG("", "finalize", "", "successful." );
+ return StatusCode::SUCCESS;
+}
+
+
+/** charged extrapolation */
+Trk::ExtrapolationCode Trk::StepEngine::extrapolate(ExCellCharged& ecCharged,
+ const Surface* sf,
+ BoundaryCheck bcheck) const
+{
+ // check the input
+ if (!ecCharged.leadParameters) return Trk::ExtrapolationCode::FailureConfiguration;
+ if (!ecCharged.leadVolume) { // should not happen : try to recover anyway
+ Trk::ExtrapolationCode eVol = m_navigationEngine->resolvePosition(ecCharged,ecCharged.propDirection,false);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ }
+
+ // cache the last lead parameters, useful in case a navigation error occured
+ ecCharged.lastLeadParameters = ecCharged.leadParameters;
+ // new extrapolation call : find possible target surfaces
+ //Trk::TargetSurfaceVector ts;
+ //Trk::ExtrapolationCode eCode = targetSurfacesT(ecCharged,ts,false,sf,bcheck);
+
+ // test new class
+ Trk::ExtrapolationCode exC = m_targetSurfaces.setOnInput(ecCharged,sf,bcheck);
+
+ // catch an infinite loop in frame navigation
+ if (exC==Trk::ExtrapolationCode::FailureLoop) { // should not happen : try to recover once
+
+ EX_MSG_WARNING("", "extrapolate", "", "frame navigation recovery attempt for:"<< ecCharged.leadVolume->volumeName());
+
+ Trk::ExtrapolationCode eVol = m_navigationEngine->resolvePosition(ecCharged,ecCharged.propDirection,true);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ //ts.clear();
+ //eVol = targetSurfacesT(ecCharged,ts,false,sf,bcheck);
+ exC = m_targetSurfaces.setOnInput(ecCharged,sf,bcheck);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ }
+
+ // build PathLimit from eCell info
+ Trk::PathLimit pathLim(ecCharged.materialLimitX0-ecCharged.materialX0, ecCharged.materialProcess);
+ // build TimeLimit from eCell info
+ //Trk::TimeLimit timeLim(ecCharged.timeLimit,ecCharged.time,201); // Fatras decay code
+ Trk::TimeLimit timeLim(-1.,ecCharged.time,201); // temporary till time limit absent in eCell
+
+ // loop over intersections
+ // std::vector<unsigned int> solutions;
+ Trk::TargetSurfaceVector solutions;
+ while (ecCharged.leadParameters ) {
+ solutions.clear();
+ // TODO pass TargetSurfaces directly to the propagator
+
+ EX_MSG_VERBOSE(ecCharged.navigationStep, "propagate", "loop", "starting propagation at position : " << ecCharged.leadParameters->position()<<","<<ecCharged.leadParameters->momentum() );
+ if (m_debugCall) EX_MSG_DEBUG(ecCharged.navigationStep, "propagate", "debug loop", "starting propagation at position : " << ecCharged.leadParameters->position()<<","<<ecCharged.leadParameters->momentum() );
+ //const Trk::TrackParameters* nextPar = m_propagator->propagateT(*ecCharged.leadParameters,m_targetSurfaces,ecCharged.propDirection,
+ // ecCharged.mFieldMode, ecCharged.pHypothesis, solutions,
+ // pathLim,timeLim,ecCharged.navigationCurvilinear,hitVector);
+ // cache the last lead parameters, useful in case a navigation error occured
+ //ecCharged.lastLeadParameters = ecCharged.leadParameters;
+
+ Trk::ExtrapolationCode eCode = m_propagator->propagate(ecCharged,m_targetSurfaces,solutions);
+
+ // enforced debugging
+ if (eCode==Trk::ExtrapolationCode::FailureNavigation && m_debugAndFix) {
+ if (m_debugCall) {
+ EX_MSG_INFO(ecCharged.navigationStep, "extrapolate", "loop:debug mode:"," stopping execution for further debugging");
+ exit(0);
+ } else {
+ // repeat the call with printout && stop the execution
+ EX_MSG_INFO(ecCharged.navigationStep, "extrapolate", "loop:debug mode:","rerun last extrapolation call"<<ecCharged.lastLeadParameters->position());
+ ecCharged.leadParameters = ecCharged.lastLeadParameters;
+ ecCharged.leadVolume = 0;
+ m_debugCall = true;
+ m_targetSurfaces.setDebugModeOn();
+ eCode = extrapolate(ecCharged,sf,bcheck);
+ exit(0);
+ }
+ }
+
+ const Trk::TrackParameters* nextPar = (eCode.isSuccess() || eCode.inProgress()) ? ecCharged.leadParameters : 0 ;
+ if (!nextPar) EX_MSG_VERBOSE(ecCharged.navigationStep, "propagate", "loop", "propagation failed ");
+ else EX_MSG_VERBOSE(ecCharged.navigationStep, "propagate", "loop", "propagated to :"<< nextPar->position());
+
+ if (!nextPar) return ( sf ? Trk::ExtrapolationCode::FailureDestination : eCode ) ;
+
+ EX_MSG_VERBOSE(ecCharged.navigationStep, "propagate", "loop", "propagation arrived at position : " << ecCharged.leadParameters->position() );
+ if (m_debugCall) EX_MSG_DEBUG(ecCharged.navigationStep, "propagate", "loop", "propagation arrived at position : " << ecCharged.leadParameters->position() );
+
+ // record the parameters as sensitive or passive depending on the surface
+ //Trk::ExtrapolationMode::eMode emode = csf.object->isActive() ? Trk::ExtrapolationMode::CollectSensitive : Trk::ExtrapolationMode::CollectPassive;
+ // fill the corresponding parameters, the material effects updator can attach material to them
+ // ecCharged.stepParameters(ecCharged.leadParameters, Trk::ExtrapolationMode::CollectPassive);
+
+ // handle extrapolation step
+ Trk::ExtrapolationCode eStep = handleIntersection(ecCharged,solutions);
+
+ if (eStep != Trk::ExtrapolationCode::InProgress ) return eStep;
+ }
+
+ return Trk::ExtrapolationCode::FailureConfiguration;
+ }
+
+
+/** neutral extrapolation */
+Trk::ExtrapolationCode Trk::StepEngine::extrapolate(ExCellNeutral& ecNeutral,
+ const Surface* sf,
+ BoundaryCheck bcheck) const
+{
+ // check the input
+
+ if (!ecNeutral.leadParameters) return Trk::ExtrapolationCode::FailureConfiguration;
+ if (!ecNeutral.leadVolume) { // should not happen : try to recover anyway
+ Trk::ExtrapolationCode eVol = m_navigationEngine->resolvePosition(ecNeutral,ecNeutral.propDirection,false);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ }
+
+ // new extrapolation call : find valid intersections
+ Trk::TargetSurfaceVector ts; ts.clear();
+ Trk::ExtrapolationCode eCode = targetSurfacesT(ecNeutral,ts,true,sf,bcheck);
+
+ CHECK_ECODE_SUCCESS(ecNeutral, eCode);
+
+ // catch an infinite loop in frame navigation
+ if (eCode==Trk::ExtrapolationCode::FailureLoop) { // should not happen : try to recover once
+ Trk::ExtrapolationCode eVol = m_navigationEngine->resolvePosition(ecNeutral,ecNeutral.propDirection,true);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ ts.clear();
+ eVol = targetSurfacesT(ecNeutral,ts,true,sf,bcheck);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+ }
+
+ // catch an infinite loop in frame navigation
+ if (!ecNeutral.leadVolume) { // should not happen : try to recover once
+ Trk::ExtrapolationCode eVol = m_navigationEngine->resolvePosition(ecNeutral,ecNeutral.propDirection);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+
+ eVol = targetSurfacesT(ecNeutral,ts,true,sf,bcheck);
+ if (!eVol.isSuccessOrRecovered() && !eVol.inProgress()) return eVol;
+
+ }
+
+ // loop over intersections
+ for (unsigned int ii=0; ii<ts.size(); ii++) {
+
+
+ ////// STEP ACROSS FRAME BOUNDARY /////////////// TO DO:TEMPLATE //////////////////////////
+ if ( ts[ii].sfType == Trk::SurfNavigType::BoundaryFrame ) {
+ const std::vector< Trk::SharedObject< const Trk::BoundarySurface< Trk::TrackingVolume> > > bounds = ecNeutral.leadVolume->boundarySurfaces();
+ const Trk::TrackingVolume* nextVolume = bounds[ts[ii].index].getPtr()->attachedVolume(ts[ii].intersection,
+ ecNeutral.leadParameters->momentum(),
+ ecNeutral.propDirection);
+
+ if (!nextVolume) return Trk::ExtrapolationCode::SuccessBoundaryReached;
+
+ ecNeutral.leadParameters = new Trk::NeutralCurvilinearParameters(ts[ii].intersection,ecNeutral.leadParameters->momentum(),0.);
+
+ // - geometrySignature change and configuration to stop then triggers a Success
+ bool stopAtThisBoundary = ecNeutral.checkConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary)
+ && (nextVolume->geometrySignature() != ecNeutral.leadVolume->geometrySignature());
+ // fill the boundary into the cache if successfully hit boundary surface
+ // - only cache if those are not the final parameters caused by a StopAtBoundary
+ if (!stopAtThisBoundary)
+ ecNeutral.stepParameters(ecNeutral.leadParameters, Trk::ExtrapolationMode::CollectBoundary);
+ // loop protection - relaxed for the cases where you start from the boundary
+ if (ecNeutral.leadVolume == nextVolume ) {
+ // the start parameters where on the boundary already give a relaxed return code
+ // if (&bSurface == eCell.lastBoundarySurface) return Trk::ExtrapolationCode::Unset;
+ // give some screen output as of why this happens
+ EX_MSG_VERBOSE(ecNeutral.navigationStep, "navigation", "", "loop detected while trying to leave TrackingVolume '" << nextVolume->volumeName() << ".");
+ // return a loop failure, parameter deletion will be done by cache
+ return Trk::ExtrapolationCode::FailureLoop;
+ }
+ // break if configured to break at volume boundary and signature change
+ if (stopAtThisBoundary){
+ EX_MSG_VERBOSE(ecNeutral.navigationStep, "navigation", "", "geometry signature change from " <<
+ ecNeutral.leadVolume->geometrySignature() << " to " << nextVolume->geometrySignature());
+ ecNeutral.nextGeometrySignature = nextVolume->geometrySignature();
+ // return the boundary reached
+ return Trk::ExtrapolationCode::SuccessBoundaryReached;
+ }
+ // remember the last boundary surface for loop protection
+ // ecNeutral.lastBoundarySurface = &bounds[index];
+ ecNeutral.lastBoundaryParameters = ecNeutral.leadParameters;
+ // set next volume
+ ecNeutral.leadVolume = nextVolume;
+ return Trk::ExtrapolationCode::InProgress;
+ }
+ }
+
+ return Trk::ExtrapolationCode::FailureConfiguration;
+ }
+
+
+void Trk::StepEngine::evaluateDistance(Trk::TargetSurface& tt, Amg::Vector3D pos, Amg::Vector3D mom,
+ Trk::TargetSurfaceVector&ts, bool trueOrdered) const
+{
+ Trk::DistanceSolution distSol = tt.surf->straightLineDistanceEstimate(pos,mom);
+
+ double dist = distSol.first();
+ Amg::Vector3D posi = pos + dist*mom;
+ // skip trivial solutions
+ if (distSol.numberOfSolutions()>1 && dist<m_tolerance && distSol.second()>m_tolerance) {
+ dist = distSol.second();
+ posi = pos + dist*mom;
+ if (trueOrdered && !tt.surf->isOnSurface(posi,tt.bcheck,m_tolerance,m_tolerance) ) return;
+ double dAbs = distSol.currentDistance(true);
+ tt.setDistance(dist,fabs(dAbs),distSol.signedDistance() && dAbs!=0. ? dAbs/fabs(dAbs) : 0.);
+ tt.setPosition(posi);
+ ts.push_back(tt);
+ return;
+ }
+ // save closest solution
+ if (!trueOrdered || dist>m_tolerance ) {
+ double dAbs = distSol.currentDistance(true);
+ tt.setDistance(dist,fabs(dAbs),distSol.signedDistance() && dAbs!=0. ? dAbs/fabs(dAbs) : 0.);
+ tt.setPosition(posi);
+ ts.push_back(tt);
+ }
+
+ // save multiple intersection for neutral transport
+ if (distSol.numberOfSolutions()>1 && distSol.second()>m_tolerance && trueOrdered) {
+ dist = distSol.second();
+ posi = pos + dist*mom;
+ if ( tt.surf->isOnSurface(posi,tt.bcheck,m_tolerance,m_tolerance) ) {
+ double dAbs = distSol.currentDistance(true);
+ tt.setDistance(dist,fabs(dAbs),distSol.signedDistance() && dAbs!=0. ? dAbs/fabs(dAbs) : 0.);
+ tt.setPosition(posi);
+ ts.push_back(tt);
+ }
+ }
+
+ return;
+}
+
+// handle extrapolation step
+Trk::ExtrapolationCode Trk::StepEngine::handleIntersection(ExCellCharged& ecCharged,
+ Trk::TargetSurfaceVector& solutions) const
+{
+ std::vector<Trk::TargetSurface>::iterator is=solutions.begin();
+ while ( is!=solutions.end() ) {
+
+ // destination
+ if ( (*is).sfType == Trk::SurfNavigType::Target ) return Trk::ExtrapolationCode::SuccessDestination;
+
+ // frame boundary : update of target surfaces
+ if ( (*is).sfType == Trk::SurfNavigType::BoundaryFrame ) {
+ return resolveFrameBoundaryT(ecCharged,ecCharged.leadParameters->position(),(*is).index);
+ }
+
+ is++;
+ }
+
+ return Trk::ExtrapolationCode::InProgress;
+}
diff --git a/Tracking/TrkExtrapolation/TrkExEngine/src/components/TrkExEngine_entries.cxx b/Tracking/TrkExtrapolation/TrkExEngine/src/components/TrkExEngine_entries.cxx
index 315bf68aaa4..e9f40a9439d 100755
--- a/Tracking/TrkExtrapolation/TrkExEngine/src/components/TrkExEngine_entries.cxx
+++ b/Tracking/TrkExtrapolation/TrkExEngine/src/components/TrkExEngine_entries.cxx
@@ -4,6 +4,7 @@
#include "TrkExEngine/StaticNavigationEngine.h"
#include "TrkExEngine/MaterialEffectsEngine.h"
#include "TrkExEngine/PropagationEngine.h"
+#include "TrkExEngine/StepEngine.h"
using namespace Trk;
@@ -12,6 +13,7 @@ DECLARE_TOOL_FACTORY( MaterialEffectsEngine )
DECLARE_TOOL_FACTORY( StaticEngine )
DECLARE_TOOL_FACTORY( StaticNavigationEngine )
DECLARE_TOOL_FACTORY( PropagationEngine )
+DECLARE_TOOL_FACTORY( StepEngine )
/** factory entries need to have the name of the package */
DECLARE_FACTORY_ENTRIES( TrkExEngine )
@@ -21,6 +23,7 @@ DECLARE_FACTORY_ENTRIES( TrkExEngine )
DECLARE_TOOL( StaticEngine )
DECLARE_TOOL( StaticNavigationEngine )
DECLARE_TOOL( PropagationEngine )
+ DECLARE_TOOL( StepEngine )
}
--
GitLab