Commit 67d2e782 authored by Whitney Armstrong's avatar Whitney Armstrong
Browse files

Merge branch 'acts_v1' into 'master'

Updating for acts v1

See merge request !10
parents 3e0effb9 55e800b1
......@@ -10,7 +10,7 @@
#include "JugReco/Utilities/GroupBy.hpp"
#include "JugReco/Utilities/Range.hpp"
#include "Acts/Geometry/GeometryID.hpp"
#include "Acts/Geometry/GeometryIdentifier.hpp"
#include <algorithm>
#include <cstddef>
......@@ -24,23 +24,23 @@ namespace detail {
// extract the geometry identifier from a variety of types
struct GeometryIdGetter {
// explicit geometry identifier are just forwarded
constexpr Acts::GeometryID operator()(Acts::GeometryID geometryId) const {
constexpr Acts::GeometryIdentifier operator()(Acts::GeometryIdentifier geometryId) const {
return geometryId;
}
// encoded geometry ids are converted back to geometry identifiers.
constexpr Acts::GeometryID operator()(Acts::GeometryID::Value encoded) const {
return Acts::GeometryID(encoded);
constexpr Acts::GeometryIdentifier operator()(Acts::GeometryIdentifier::Value encoded) const {
return Acts::GeometryIdentifier(encoded);
}
// support elements in map-like structures.
template <typename T>
constexpr Acts::GeometryID operator()(
const std::pair<Acts::GeometryID, T>& mapItem) const {
constexpr Acts::GeometryIdentifier operator()(
const std::pair<Acts::GeometryIdentifier, T>& mapItem) const {
return mapItem.first;
}
// support elements that implement `.geometryId()`.
template <typename T>
inline auto operator()(const T& thing) const
-> decltype(thing.geometryId(), Acts::GeometryID()) {
-> decltype(thing.geometryId(), Acts::GeometryIdentifier()) {
return thing.geometryId();
}
};
......@@ -84,17 +84,17 @@ using GeometryIdMultiset =
/// }
///
template <typename T>
using GeometryIdMultimap = GeometryIdMultiset<std::pair<Acts::GeometryID, T>>;
using GeometryIdMultimap = GeometryIdMultiset<std::pair<Acts::GeometryIdentifier, T>>;
/// Select all elements within the given volume.
template <typename T>
inline Range<typename GeometryIdMultiset<T>::const_iterator> selectVolume(
const GeometryIdMultiset<T>& container, Acts::GeometryID::Value volume) {
auto cmp = Acts::GeometryID().setVolume(volume);
const GeometryIdMultiset<T>& container, Acts::GeometryIdentifier::Value volume) {
auto cmp = Acts::GeometryIdentifier().setVolume(volume);
auto beg = std::lower_bound(container.begin(), container.end(), cmp,
detail::CompareGeometryId{});
// WARNING overflows to volume==0 if the input volume is the last one
cmp = Acts::GeometryID().setVolume(volume + 1u);
cmp = Acts::GeometryIdentifier().setVolume(volume + 1u);
// optimize search by using the lower bound as start point. also handles
// volume overflows since the geo id would be located before the start of
// the upper edge search window.
......@@ -104,20 +104,20 @@ inline Range<typename GeometryIdMultiset<T>::const_iterator> selectVolume(
}
template <typename T>
inline auto selectVolume(const GeometryIdMultiset<T>& container,
Acts::GeometryID id) {
Acts::GeometryIdentifier id) {
return selectVolume(container, id.volume());
}
/// Select all elements within the given layer.
template <typename T>
inline Range<typename GeometryIdMultiset<T>::const_iterator> selectLayer(
const GeometryIdMultiset<T>& container, Acts::GeometryID::Value volume,
Acts::GeometryID::Value layer) {
auto cmp = Acts::GeometryID().setVolume(volume).setLayer(layer);
const GeometryIdMultiset<T>& container, Acts::GeometryIdentifier::Value volume,
Acts::GeometryIdentifier::Value layer) {
auto cmp = Acts::GeometryIdentifier().setVolume(volume).setLayer(layer);
auto beg = std::lower_bound(container.begin(), container.end(), cmp,
detail::CompareGeometryId{});
// WARNING resets to layer==0 if the input layer is the last one
cmp = Acts::GeometryID().setVolume(volume).setLayer(layer + 1u);
cmp = Acts::GeometryIdentifier().setVolume(volume).setLayer(layer + 1u);
// optimize search by using the lower bound as start point. also handles
// volume overflows since the geo id would be located before the start of
// the upper edge search window.
......@@ -127,25 +127,25 @@ inline Range<typename GeometryIdMultiset<T>::const_iterator> selectLayer(
}
template <typename T>
inline auto selectLayer(const GeometryIdMultiset<T>& container,
Acts::GeometryID id) {
Acts::GeometryIdentifier id) {
return selectLayer(container, id.volume(), id.layer());
}
/// Select all elements for the given module / sensitive surface.
template <typename T>
inline Range<typename GeometryIdMultiset<T>::const_iterator> selectModule(
const GeometryIdMultiset<T>& container, Acts::GeometryID geoId) {
const GeometryIdMultiset<T>& container, Acts::GeometryIdentifier geoId) {
// module is the lowest level and defines a single geometry id value
return makeRange(container.equal_range(geoId));
}
template <typename T>
inline auto selectModule(const GeometryIdMultiset<T>& container,
Acts::GeometryID::Value volume,
Acts::GeometryID::Value layer,
Acts::GeometryID::Value module) {
Acts::GeometryIdentifier::Value volume,
Acts::GeometryIdentifier::Value layer,
Acts::GeometryIdentifier::Value module) {
return selectModule(
container,
Acts::GeometryID().setVolume(volume).setLayer(layer).setSensitive(
Acts::GeometryIdentifier().setVolume(volume).setLayer(layer).setSensitive(
module));
}
......
......@@ -24,7 +24,7 @@ namespace Jug {
std::size_t hitCount;
};
using IndexedParams = std::unordered_map<size_t, Acts::BoundParameters>;
using IndexedParams = std::unordered_map<size_t, Acts::BoundTrackParameters>;
/// @brief Struct for truth track fitting/finding result with
/// Acts::KalmanFitter/Acts::CombinatorialKalmanFilter
......@@ -128,7 +128,7 @@ struct SimMultiTrajectory {
/// @param entryIndex The trajectory entry index
///
/// @return The fitted track parameters of the trajectory
const Acts::BoundParameters& trackParameters(const size_t& entryIndex) const {
const Acts::BoundTrackParameters& trackParameters(const size_t& entryIndex) const {
auto it = m_trackParameters.find(entryIndex);
if (it != m_trackParameters.end()) {
return it->second;
......
......@@ -40,7 +40,7 @@
// SimSourceLink& operator=(SimSourceLink&&) = default;
// SimSourceLink& operator=(const SimSourceLink&) = default;
//
// constexpr Acts::GeometryID geometryId() const { return m_geometryId; }
// constexpr Acts::GeometryIdentifier geometryId() const { return m_geometryId; }
// constexpr const Acts::Surface& referenceSurface() const { return *m_surface; }
// constexpr const ActsFatras::Hit& truthHit() const { return *m_truthHit; }
//
......@@ -68,7 +68,7 @@
// Acts::BoundMatrix m_cov;
// size_t m_dim = 0u;
// // store geo id copy to avoid indirection via truth hit
// Acts::GeometryID m_geometryId;
// Acts::GeometryIdentifier m_geometryId;
// // need to store pointers to make the object copyable
// const Acts::Surface* m_surface;
// const ActsFatras::Hit* m_truthHit;
......
......@@ -26,7 +26,7 @@ class SourceLink {
Acts::BoundMatrix m_cov;
size_t m_dim = 0u;
// store geo id copy to avoid indirection via truth hit
Acts::GeometryID m_geometryId;
Acts::GeometryIdentifier m_geometryId;
// need to store pointers to make the object copyable
const Acts::Surface* m_surface;
//const ActsFatras::Hit* m_truthHit;
......@@ -47,7 +47,7 @@ class SourceLink {
SourceLink& operator=(SourceLink&&) = default;
SourceLink& operator=(const SourceLink&) = default;
constexpr Acts::GeometryID geometryId() const { return m_geometryId; }
constexpr Acts::GeometryIdentifier geometryId() const { return m_geometryId; }
constexpr const Acts::Surface& referenceSurface() const { return *m_surface; }
//constexpr const ActsFatras::Hit& truthHit() const { return *m_truthHit; }
......@@ -55,12 +55,12 @@ class SourceLink {
if (m_dim == 0) {
throw std::runtime_error("Cannot create dim 0 measurement");
} else if (m_dim == 1) {
return Acts::Measurement<SourceLink, Acts::BoundParametersIndices,
return Acts::Measurement<SourceLink, Acts::BoundIndices,
Acts::eBoundLoc0>{
m_surface->getSharedPtr(), *this, m_cov.topLeftCorner<1, 1>(),
m_values[0]};
} else if (m_dim == 2) {
return Acts::Measurement<SourceLink, Acts::BoundParametersIndices,
return Acts::Measurement<SourceLink, Acts::BoundIndices,
Acts::eBoundLoc0, Acts::eBoundLoc1>{
m_surface->getSharedPtr(), *this, m_cov.topLeftCorner<2, 2>(),
m_values[0], m_values[1]};
......
......@@ -23,7 +23,7 @@ namespace Jug {
/// (Reconstructed) track parameters e.g. close to the vertex.
using TrackParameters = Acts::CurvilinearParameters;
using TrackParameters = Acts::CurvilinearTrackParameters;
/// Container of reconstructed track states for multiple tracks.
using TrackParametersContainer = std::vector<TrackParameters>;
......
......@@ -63,15 +63,15 @@ namespace Jug {
//debug() << trackstate.hasPredicted() << endmsg;
//debug() << trackstate.predicted() << endmsg;
auto params = trackstate.predicted() ;//<< endmsg;
debug() << 1.0/params[Acts::eQOP] << " GeV" << endmsg;
if ( std::abs(1.0 / params[Acts::eQOP]) > 10) {
debug() << 1.0/params[Acts::eBoundQOverP] << " GeV" << endmsg;
if ( std::abs(1.0 / params[Acts::eBoundQOverP]) > 10) {
debug() << "skipping" << endmsg;
return;
}
eic::Particle p({params[Acts::ePHI], params[Acts::eTHETA], 1.0 / params[Acts::eQOP], 0.105},
{0.0, 0.0, 0.0, params[Acts::eT]},
(long long)13 * params[Acts::eQOP] / std::abs(params[Acts::eQOP]), 0);
eic::Particle p({params[Acts::eBoundPhi], params[Acts::eBoundTheta], 1.0 / params[Acts::eBoundQOverP], 0.105},
{0.0, 0.0, 0.0, params[Acts::eBoundTime]},
(long long)13 * params[Acts::eBoundQOverP] / std::abs(params[Acts::eBoundQOverP]), 0);
debug() << p << endmsg;
rec_parts->push_back(p);
});
......
......@@ -37,11 +37,11 @@
#include <vector>
#include "Acts/Geometry/TrackingGeometry.hpp"
#include "Acts/TrackFinder/CKFSourceLinkSelector.hpp"
#include "Acts/TrackFinder/CombinatorialKalmanFilter.hpp"
#include "Acts/TrackFinding/CKFSourceLinkSelector.hpp"
#include "Acts/TrackFinding/CombinatorialKalmanFilter.hpp"
#include "Acts/Fitter/GainMatrixSmoother.hpp"
#include "Acts/Fitter/GainMatrixUpdater.hpp"
#include "Acts/TrackFitting/GainMatrixSmoother.hpp"
#include "Acts/TrackFitting/GainMatrixUpdater.hpp"
#include "Acts/MagneticField/ConstantBField.hpp"
#include "Acts/MagneticField/InterpolatedBFieldMap.hpp"
#include "Acts/MagneticField/SharedBField.hpp"
......@@ -75,7 +75,7 @@ namespace Jug::Reco {
m_BField = std::make_shared<Acts::ConstantBField>(Acts::Vector3D{0.0, 0.0, 1.0_T});
m_fieldctx = BFieldVariant(m_BField);
m_sourcelinkSelectorCfg = {
{Acts::GeometryID(), {15, 10}},
{Acts::GeometryIdentifier(), {15, 10}},
};
findTracks = TrackFindingAlgorithm::makeTrackFinderFunction(m_geoSvc->trackingGeometry(), m_BField);
......@@ -109,9 +109,12 @@ namespace Jug::Reco {
for (std::size_t iseed = 0; iseed < init_trk_params->size(); ++iseed) {
const auto& initialParams = (*init_trk_params)[iseed];
Acts::PropagatorPlainOptions pOptions;
pOptions.maxSteps = 10000;
// Set the CombinatorialKalmanFilter options
TrackFindingAlgorithm::CKFOptions ckfOptions(m_geoctx, m_fieldctx, m_calibctx, m_sourcelinkSelectorCfg,
Acts::LoggerWrapper{logger()}, &(*pSurface));
CKFOptions ckfOptions(m_geoctx, m_fieldctx, m_calibctx, m_sourcelinkSelectorCfg,
Acts::LoggerWrapper{logger()}, pOptions,
&(*pSurface));
// TrackFindingAlgorithm::CKFOptions ckfOptions(ctx.geoContext, ctx.magFieldContext, ctx.calibContext,
// m_cfg.sourcelinkSelectorCfg, Acts::LoggerWrapper{logger()},
// &(*pSurface));
......
......@@ -32,8 +32,8 @@
//#include "Acts/Surfaces/PerigeeSurface.hpp"
#include "Acts/Geometry/TrackingGeometry.hpp"
#include "Acts/TrackFinder/CKFSourceLinkSelector.hpp"
#include "Acts/TrackFinder/CombinatorialKalmanFilter.hpp"
#include "Acts/TrackFinding/CKFSourceLinkSelector.hpp"
#include "Acts/TrackFinding/CombinatorialKalmanFilter.hpp"
//#include "Acts/Fitter/GainMatrixSmoother.hpp"
//#include "Acts/Fitter/GainMatrixUpdater.hpp"
......@@ -57,6 +57,7 @@ namespace Jug::Reco {
/// Track finding function that takes input measurements, initial trackstate
/// and track finder options and returns some track-finding-specific result.
using CKFOptions = Acts::CombinatorialKalmanFilterOptions<Acts::CKFSourceLinkSelector>;
using TrackFinderFunction =
std::function<TrackFinderResult(const SourceLinkContainer&, const TrackParameters&, const CKFOptions&)>;
......
......@@ -57,17 +57,19 @@ namespace Jug::Reco {
// build some track cov matrix
Acts::BoundSymMatrix cov = Acts::BoundSymMatrix::Zero();
cov(Acts::eLOC_0, Acts::eLOC_0) = 0.1 * mm*0.1 * mm;
cov(Acts::eLOC_1, Acts::eLOC_1) = 0.1 * mm*0.1 * mm;
cov(Acts::ePHI, Acts::ePHI) = M_PI / 180.0;
cov(Acts::eTHETA, Acts::eTHETA) = M_PI / 180.0;
cov(Acts::eQOP, Acts::eQOP) = 1.0 / (0.3 * GeV* 0.3 * GeV);
cov(Acts::eT, Acts::eT) = Acts::UnitConstants::ns;
init_trk_params->emplace_back(std::make_optional(std::move(cov)),
Acts::Vector3D(part.vsx() * mm, part.vsy() * mm, part.vsz() * mm),
cov(Acts::eBoundLoc0, Acts::eBoundLoc0) = 0.1 * mm*0.1 * mm;
cov(Acts::eBoundLoc1, Acts::eBoundLoc1) = 0.1 * mm*0.1 * mm;
cov(Acts::eBoundPhi, Acts::eBoundPhi) = M_PI / 180.0;
cov(Acts::eBoundTheta, Acts::eBoundTheta) = M_PI / 180.0;
cov(Acts::eBoundQOverP, Acts::eBoundQOverP) = 1.0 / (0.3 * GeV* 0.3 * GeV);
cov(Acts::eBoundTime, Acts::eBoundTime) = Acts::UnitConstants::ns;
init_trk_params->emplace_back(Acts::Vector4D(part.vsx() * mm, part.vsy() * mm, part.vsz() * mm, part.time() * Acts::UnitConstants::ns),
Acts::Vector3D(part.psx() * GeV, part.psy() * GeV, part.psz() * GeV),
((part.pdgID() > 0) ? -1 : 1), part.time() * Acts::UnitConstants::ns);
std::sqrt(part.psx() *part.psx() + part.psy() * part.psy() + part.psz() * part.psz())*GeV,
((part.pdgID() > 0) ? -1 : 1),
std::make_optional(std::move(cov))
);
//part .charge()
debug() << "Invoke track finding seeded by truth particle " << part << endmsg;
......
......@@ -98,8 +98,8 @@ namespace Jug::Reco {
for(const auto& ahit : *hits) {
Acts::BoundMatrix cov = Acts::BoundMatrix::Zero();
cov(Acts::eLOC_0, Acts::eLOC_0) = ahit.covMatrix(0)*Acts::UnitConstants::mm;//*ahit.covMatrix(0);
cov(Acts::eLOC_1, Acts::eLOC_1) = ahit.covMatrix(1)*Acts::UnitConstants::mm;//*ahit.covMatrix(1);
cov(Acts::eBoundLoc0, Acts::eBoundLoc0) = ahit.covMatrix(0)*Acts::UnitConstants::mm;//*ahit.covMatrix(0);
cov(Acts::eBoundLoc1, Acts::eBoundLoc1) = ahit.covMatrix(1)*Acts::UnitConstants::mm;//*ahit.covMatrix(1);
debug() << "cell ID : " << ahit.cellID0() << endmsg;
auto vol_ctx = m_geoSvc->cellIDPositionConverter()->findContext(ahit.cellID0());
......@@ -117,7 +117,7 @@ namespace Jug::Reco {
Acts::Vector2D pos(0, 0);
// geometry context contains nothing here
surface->globalToLocal(Acts::GeometryContext(), {ahit.position(0), ahit.position(1), ahit.position(2)},
{0, 0, 0}, pos);
{0, 0, 0});//, pos);
//// smear truth to create local measurement
Acts::BoundVector loc = Acts::BoundVector::Zero();
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment